Thursday, June 23, 2011

Hybrid CMOS-MQCA Architecture using Multi-layer Spintronic Devices (WIP/DAC)

J. Das, S. M. Alam, S. Rajaram and S. Bhanja, “Hybrid CMOS-MQCA Architecture using
Multi-layer Spintronic Devices, Accepted as WIP for IEEE/ACM Design Automation
Conference (DAC), 2011

Wednesday, June 15, 2011

Dr. Bhanja is a member of the Editorial Board of ACM JETC

Dr. Bhanja is a member of the Editorial Board IEEE Transactions on VLSI Systems (2011-2012)

Dr. Bhanja receives USF Outstanding Undergraduate Teacher award 2010

Dr. Bhanja receives William R. Jones Outstanding Mentor Awards, Florida Education Foundation, 2010

Low Power CMOS-Magnetic Nano-Logic With Increased Bit Controllability (IEEE Nano 2011)

J. Das, S. M. Alam, and S. Bhanja, “Low Power CMOS-Magnetic Nano-Logic With Increased Bit Controllability”, to IEEE Conference on Nanotechnology, 2011

Experimental Demonstration of Viability of Energy Minimizing Computing using Nano-magnets (IEEE Nano 2011)

J. Pulecio, S. Sarkar and S. Bhanja, “Experimental Demonstration of Viability of Energy Minimizing Computing using Nano-magnets”, to IEEE Conference on Nanotechnology, 2011

Novel knowledge module on fusion of logic and memory to undergraduate students (IEEE MSE)

D. Karunaratne, S. Rajaram, P. De, K. Kusmierek and S. Bhanja, “Novel knowledge module on fusion of logic and memory to undergraduate students”, to IEEE Microelectronic System Education Conference, 2011.

Tool for Analysis and Quantification of Fabrication Layouts in Nanomagnet-based Computing (IEEE NANO 2011)

R. Panchumarthy, D. Karunaratne, S. Sarkar and S. Bhanja, “Tool for Analysis and Quantification of Fabrication Layouts in Nanomagnet-based Computing”, to IEEE Conference on Nanotechnology, 2011.

A Review of Magnetic Cellular Automata Systems (Invited to IEEE ISCAS), 2011

S. Bhanja and J. Pulecio, “A Review of Magnetic Cellular Automata Systems”, (Invited paper) to IEEE International Symposium on Circuits and Systems (ISCAS), 2011.

QCAPro - An Error-Power Estimation Tool for QCA Circuit Design (Invited paper to IEEE ISCAS), 2011

S. Srivastava, A. Asthana, S. Bhanja and S. Sarkar “QCAPro - An Error-Power Estimation Tool for QCA Circuit Design”, (Invited paper) to IEEE International Symposium on Circuits and Systems (ISCAS), 2011. Download the tool

Javier Pulecio recives USF Hispanic Pathway Award, 2010

Low Power Magnetic Quantum Cellular Automata Realization Using Magnetic Multi-Layer Structures (Invited Paper in IEEE JETCAS) 2011

J. Das, S.M. Alam and S. Bhanja, "Low Power Magnetic Quantum Cellular Automata Realization Using Magnetic Multi-Layer Structures", in IEEE Journal of Emerging Technologies CAS, June 2011 (Invited Paper).

Abstract: In this paper, we report Magnetic Quantum Cellular Automata (MQCA) realization using multi-layer cells with tilted polarizer reference layer with a particular focus on the critical
need to shift towards the multi-layer cells as elemental entities from the conventional single-domain nanomagnets. We have reported a novel spin-transfer torque current induced clocking
scheme, theoretically derived the clocking current, and shown the reduction in power consumption achieved against the traditional mechanism of clocking using magnetic fields typically generated from overhead or underneath wires. We have modeled the multilayer
cell behavior in Verilog-A along with the underlying algorithm used in implementing the neighbor interaction between the cells. The paper reports the switching and clocking current magnitudes, their direction and the power consumption associated with switching and clocking operation. Finally, we present the simulation results from Verilog-A model of switching, clocking
and neighbor interaction. Low power consumption due to spin transfer torque current induced switching and clocking along with the reasonable Magneto-Resistance (MR) distinguishing the
two energy minimum states of the device, make these devices a promising candidate in MQCA realization.

Wednesday, June 8, 2011

Verilog-A Model File for Magnetic Tunnel Junction (MTJ) Cells in a Vertical Array

//This VerilogA program models the behavior of a multi-layer cell (MTJ) in a horizontal array.


`include "constants.vams"

nature Current
units = "A";
access = I;
abstol = 1e-12;
blowup = 1e32;
endnature

nature Voltage
units = "V";
access = V;
abstol = 1e-12;
blowup = 1e32;
endnature

discipline electrical
potential Voltage;
flow Current;
enddiscipline

discipline voltage
potential Voltage;
enddiscipline

discipline current
flow Current;
enddiscipline

nature Magnetic_induction
abstol = 1e-12;
access = H;
units = "A/m";
endnature

discipline magnetic_ports
potential Magnetic_induction;
enddiscipline

//rf and rp are the terminals of the free and the fixed layers.
//mag_a, mag_b, mag_c and mag_d are the 2-bit magnetic ports that are used to communicate with neighboring cells.
//The magnetic ports are used to model the neighbor interaction among the free layers of the multi-layer cells.
module mtj_vert_cell(rf,rp,mag_porta,mag_portb,mag_portc,mag_portd);

inout rf,rp;
inout [0:1]mag_porta,mag_portb,mag_portc,mag_portd;

electrical rf,rp;
magnetic_ports [0:1] mag_porta,mag_portb,mag_portc,mag_portd;


//The gyromagnetic ratio
`define gamma -1.76e11

// Paramter constants:
//Length, Width, thickness --- dimensions of the free layer of the MTJ.
//alpha --- Gilbert damping factor.
//R0 and R1 --- resistances of the logic 0 and logic 1 states.
//P --- polarization factor.
//gamma_p, alpha_p --- direction cosines of the magnetization of the reference layer along the x and z directioni.
//Hd --- Coupling field from the reference layer onto the free layer.
//Hk --- Anisotropy field.
//Gp,Gap --- conductances of MTJ when the magnetization of the free layer is parallel and anti-parallel to the magnetization of the reference layer.
//theta0 --- angle subtended by the magnetization of the reference layer with the z-axis.
//I_delta --- deviation in the clocking current magnitude about the calculated value.

//low, high --- low and high signal level for the magnetic ports.
//V_delta --- differential voltage to determine the voltage across the device crossing the zero level.
//tdelay, trise, tfall --- transition delay, rise time and fall time of the signal at the magnetic ports.
//t_delta --- model the delay in the device to respond to the switching and clocking voltage pulse.


parameter real Length = 100e-9,
Width = 50e-9,
thickness = 2e-9,
alpha = 0.01,
Ms = 8e5, //in A/m
R0 = 412.72, //in ohms
R1 = 423.13, //in ohms
P=0.35,
gamma_p = 0.707,
alpha_p = 0.707,
Hd = 5000, //in A/m
Hk = 795.77, //in A/m
I_upper=200e-6, //in A
Gp = 2.82e-3, //in mho
Gap = 1.967e-3, //in mho
theta0 = `M_PI/4,
I_delta = 1e-6; //in A

parameter real low=0,
high=1,
V_delta=10e-3;

parameter real tdelay=2p from [0:inf),
trise=1p from [0:inf),
tfall = 1p from [0:inf),
t_delta = 1p from [0:inf);


electrical int_nodea,int_nodeb;
real Resistance,I_max,I_prev,tpulse,t_prev_pulse;
real theta,Volume,t_precession,I_switch,fsm_out;
real trise_pos,tfall_neg,tfall_pos,trise_neg,tswitch;
real g,geff,aj_piby2,Iclk;

integer out[0:1],out_prev[0:1];
integer value_porta[0:1],value_portb[0:1],value_portc[0:1];
integer switch_out,switch_prev_out,clocked;



// This function computes the critical switching current for a multi-layer stack with the given device specifications.


analog function real critical_switching_current;

input Volume,I_rfrp;
real Volume,I_rfrp,eta;
begin
if(I_rfrp < 0) begin
eta = P/(1+pow(P,2));
critical_switching_current = (2*`P_Q/(`P_H/(2*`M_PI)))*(alpha*Ms*Volume/eta)*`P_U0*(Ms/(2*sqrt(2)));
end
else begin
eta = P/(1-pow(P,2));
critical_switching_current = (2*`P_Q/(`P_H/(2*`M_PI)))*(alpha*Ms*Volume/eta)*`P_U0*(Ms/(2*sqrt(2)));
end
end
endfunction



// This function computes the angle of magnetization of the free layer of the multi-layer device as function of the current through it.


analog function real angle;

input aj_piby2;
real aj_piby2,cos_angle;
begin
cos_angle = (Hd*gamma_p*(4*`M_PI/1e7) - aj_piby2/alpha)/((4*`M_PI*Ms + 0.5*Hk)*(4*`M_PI/1e7));
angle = acos(cos_angle);
end
endfunction



// This function computes the resistance of the multi-layer stack as a function of the angle of magnetization of the free layer.


analog function real resistance;

input theta;
real theta,G_theta;
begin
G_theta = 0.5*((1+cos(theta))*Gp + (1-cos(theta))*Gap);
resistance = 1/G_theta;
end
endfunction


// This function computes the magnitude of the clocking current that is required to align the magnetization of the free layer along the y-axis into a
// stationary magnetization state.


analog function real clocking_current;

input Volume;
real Volume;
real Hdz,G,Jp;
begin
Hdz = Hd*cos(theta0);
G = 1/(-4 + (pow((1+P),3))*3/(4*pow(P,1.5)));
Jp = `P_U0*(pow(Ms,2))*(`P_Q)*Volume/((`P_H)/(2*`M_PI));

clocking_current = Hdz*Jp/(G*(sin(theta0))*Ms);
end
endfunction


// This function is used to compute the post-clocking neighbor interaction among the free layers of the neighboring multi-layer cells.
// The neighbor interaction is implemented as a Finite State Machine.


analog function real mtj_fsm_vert_cell;

input [0:1] value_porta,value_portb,value_portc;
integer value_porta[0:1],value_portb[0:1],value_portc[0:1];
integer portb_and,porta_and,portb_xor_porta;
integer bitwise_xor[0:2];

integer out1[0:1],out2[0:1],out3[0:1],out[0:1];
begin
bitwise_xor[0] = (value_porta[0] ^ value_porta[1]);
bitwise_xor[1] = (value_portb[0] ^ value_portb[1]);
bitwise_xor[2] = (value_portc[0] ^ value_portc[1]);

portb_and = (value_portb[0] && value_portb[1]);
porta_and = (value_porta[0] && value_porta[1]);
portb_xor_porta = portb_and ^ porta_and;

out1[0] = value_portb[0] && !value_portc[0];
out1[1] = value_portb[1] && !value_portc[1];

out2[0] = (value_portb[0] || !value_portc[0]);
out2[1] = (value_portb[1] || !value_portc[1]);
out3[0] = !value_porta[0];
out3[1] = !value_porta[1];

out2[0] = out2[0] && out3[0];
out2[1] = out2[1] && out3[1];
out1[0] = out1[0] || out2[0];
out1[1] = out1[1] || out2[1];

if(bitwise_xor[0]) begin
if(bitwise_xor[1] && bitwise_xor[2]) begin
out[0] = low;
out[1] = high;
end
else if(bitwise_xor[1]) begin
out[0] = value_portc[0];
out[1] = value_portc[1];
end
else if(bitwise_xor[2]) begin
out[0] = !value_portb[0];
out[1] = !value_portb[1];
end
else begin
out[0] = value_portc[0];
out[1] = value_portc[1];
end
end
else begin
if(bitwise_xor[1] && bitwise_xor[2]) begin
out[0] = !(value_porta[0]);
out[1] = !(value_porta[1]);
end
else if(bitwise_xor[1])
begin
out[0] = value_portc[0];
out[1] = value_portc[1];
end
else if(bitwise_xor[2]) begin
if(portb_xor_porta) begin
out[0] = low;
out[1] = high;
end
else begin
out[0] = !(value_porta[0]);
out[1] = !(value_porta[1]);
end
end
else begin
out[0] = out1[0];
out[1] = out1[1];
end
end

if(out[0] ^ out[1])
mtj_fsm_vert_cell = 0.5;
else begin
if(out[0] == 0)
mtj_fsm_vert_cell = 0.0;
else
mtj_fsm_vert_cell = 1.0;
end
end
endfunction


// This function is used to decide the specific writing or clocking or no operation that needs to be performed on the cell.
// The specific operation that needs to be performed depends on the magnitude of the current through the cell and the duration of the current pulse.


analog function integer switch_output;

input [0:1]out_prev,switch_prev_out;
input I_switch, Iclk, t_precession, tpulse, I_max;
integer out_prev[0:1],switch_prev_out;
real I_switch,Iclk,t_precession,tpulse,I_max;
integer out[0:1];
real Imax_abs,Iswitch_abs;

begin
Imax_abs = abs(I_max);
Iswitch_abs = abs(I_switch);

if(abs(Imax_abs - Iclk) < I_delta) begin
if(I_max > 0) begin
out[0] = high;
out[1] = low;
switch_output = 01;
end
else begin
out[0] = out_prev[0];
out[1] = out_prev[1];
switch_output = 10;
end
end
else if((Imax_abs > Iswitch_abs) && (tpulse > t_precession/2)) begin
if(I_max < 0) begin
out[0] = high;
out[1] = high;
switch_output = 11;
end
else begin
out[0] = low;
out[1] = low;
switch_output = 00;
end
end
else begin
out[0] = out_prev[0];
out[1] = out_prev[1];
switch_output = switch_prev_out;
end
end
endfunction


// Main body of the program.


analog
begin
Volume = Length*Width*thickness;
t_precession = abs(`M_PI/`gamma);
geff = sqrt(P) + 1/sqrt(P);
g = 1/(-4 + (pow(geff,3))*0.75);
aj_piby2 = ((`P_H/(2*`M_PI))/(2*`P_Q))*g*I(rf,rp)/(Ms*Volume);
theta = angle(aj_piby2);
Resistance = resistance(theta) ;
I_switch = critical_switching_current(Volume,I(rf,rp));
Iclk = clocking_current(Volume);

@(initial_step) begin
trise_pos = 0;
tfall_pos = 0;
tfall_neg = 0;
trise_neg = 0;
tswitch = 0;
out[0] = low;
out[1] = low;
out_prev[0] = out[0];
out_prev[1] = out[1];
clocked = 0;
Resistance = R0;
I_prev = 0;
t_prev_pulse = 0;
switch_prev_out = 10;
end

@(cross(V(rf,rp) - V_delta,+1)) begin
trise_pos = $abstime;
end

@(cross(V(rf,rp) + V_delta,-1)) begin
tfall_neg = $abstime;
end

@(cross(V(rf,rp) - V_delta,-1)) begin
tfall_pos = $abstime;
t_prev_pulse = tpulse;
tpulse = tfall_pos - trise_pos;
tswitch = tfall_pos + t_delta;
end

@(cross(V(rf,rp) + V_delta,+1)) begin
trise_neg = $abstime;
t_prev_pulse = tpulse;
tpulse = trise_neg - tfall_neg;
tswitch = trise_neg + t_delta;
end

if (tpulse < 1.5e-12) begin
tpulse = t_prev_pulse;
end


@(timer(tswitch)) begin

switch_out = switch_output(out_prev,switch_prev_out,I_switch,Iclk,t_precession,tpulse,I_max);

switch_prev_out = switch_out;
case (switch_out)
00 : begin out[0] = low; out[1] = low; clocked=0; end
11 : begin out[0] = high; out[1] = high; clocked=0; end
01 : begin out[0] = high; out[1] = low; clocked=1; end
default : begin out[0] = out_prev[0]; out[1] = out_prev[1]; clocked=0; end
endcase
end


V(int_nodea,int_nodeb) <+ ddt(V(rf,rp));
if((abs(V(int_nodea,int_nodeb)) < 2e10) && ((abs(I(rf,rp))>I_upper) || (abs(I(rf,rp)-Iclk)
I_max = I(rf,rp);
I_prev = I_max;
end
else begin
I_max = I_prev;
end


value_porta[0] = (H(mag_porta[0]) > 0);
value_porta[1] = (H(mag_porta[1]) > 0);
value_portb[0] = (H(mag_portb[0]) > 0);
value_portb[1] = (H(mag_portb[1]) > 0);
value_portc[0] = (H(mag_portc[0]) > 0);
value_portc[1] = (H(mag_portc[1]) > 0);


H(mag_portd[0]) <+ transition(out[0],tdelay,trise,tfall);
H(mag_portd[1]) <+ transition(out[1],tdelay,trise,tfall);



if(clocked == 1) begin

fsm_out = mtj_fsm_vert_cell(value_porta,value_portb,value_portc);
if((fsm_out - 0.5) == 0) begin
out[0] = low;
out[1] = high;
end
else if(fsm_out == 1) begin
out[0] = high;
out[1] = high;
end
else begin
out[0] = low;
out[1] = low;
end
clocked = 0;
end
out_prev[0] = out[0];
out_prev[1] = out[1];


if(out[0] ^ out[1]) begin
Resistance = resistance(theta);
end
else begin
if(out[0] == 0) begin
Resistance = R0;
end
else begin
Resistance = R1;
end
end
I(rf,rp) <+ V(rf,rp)/Resistance;

end
endmodule

Wednesday, June 1, 2011

Verilog-A Model File for Magnetic Tunnel Junction (MTJ) Cells in a Horizontal Array

//This VerilogA program models the behavior of a multi-layer cell (MTJ) in a horizontal array.

`include "constants.vams"

nature Current
units = "A";
access = I;
abstol = 1e-12;
blowup = 1e32;
endnature

nature Voltage
units = "V";
access = V;
abstol = 1e-12;
blowup = 1e32;
endnature

discipline electrical
potential Voltage;
flow Current;
enddiscipline

discipline voltage
potential Voltage;
enddiscipline

discipline current
flow Current;
enddiscipline

nature Magnetic_induction
abstol = 1e-12;
access = H;
units = "A/m";
endnature

discipline magnetic_ports
potential Magnetic_induction;
enddiscipline

//rf and rp are the terminals of the free and the fixed layers.
//a, mag_b, mag_c and mag_d are the 2-bit magnetic ports that are used to communicate with neighboring cells.
//The magnetic ports are used to model the neighbor interaction among the free layers of the multi-layer cells.
module mtj_horz_cell(rf,rp,mag_porta,mag_portb,mag_portc,mag_portd);

inout rf,rp;
inout [0:1]mag_porta,mag_portb,mag_portc,mag_portd;

electrical rf,rp;
magnetic_ports [0:1] mag_porta,mag_portb,mag_portc,mag_portd;


//The gyromagnetic ratio
`define gamma -1.76e11

// Paramter constants:
//Length, Width, thickness --- dimensions of the free layer of the MTJ.
//alpha --- Gilbert damping factor.
//R0 and R1 --- resistances of the logic 0 and logic 1 states.
//P --- polarization factor.
//gamma_p, alpha_p --- direction cosines of the magnetization of the reference layer along the x and z directioni.
//Hd --- Coupling field from the reference layer onto the free layer.
//Hk --- Anisotropy field.
//Gp,Gap --- conductances of MTJ when the magnetization of the free layer is parallel and anti-parallel to the magnetization of the reference layer.
//theta0 --- angle subtended by the magnetization of the reference layer with the z-axis.
//I_delta --- deviation in the clocking current magnitude about the calculated value.

//low, high --- low and high signal level for the magnetic ports.
//V_delta --- differential voltage to determine the voltage across the device crossing the zero level.
//tdelay, trise, tfall --- transition delay, rise time and fall time of the signal at the magnetic ports.
//t_delta --- model the delay in the device to respond to the switching and clocking voltage pulse.


parameter real Length = 100e-9,
Width = 50e-9,
thickness = 2e-9,
alpha = 0.01,
Ms = 8e5, //in A/m
R0 = 412.72, //in ohms
R1 = 423.13, //in ohms
P=0.35,
gamma_p = 0.707,
alpha_p = 0.707,
Hd = 5000, //in A/m
Hk = 795.77, //in A/m
I_upper=200e-6, //in A
Gp = 2.82e-3, //in mho
Gap = 1.967e-3, //in mho
theta0 = `M_PI/4,
I_delta = 1e-6; //in A

parameter real low=0,
high=1,
V_delta=10e-3;

parameter real tdelay=2p from [0:inf),
trise=1p from [0:inf),
tfall = 1p from [0:inf),
t_delta = 1p from [0:inf);


electrical int_nodea,int_nodeb;
real Resistance,I_max,I_prev,tpulse,t_prev_pulse;
real theta,Volume,t_precession,I_switch,fsm_out;
real trise_pos,tfall_neg,tfall_pos,trise_neg,tswitch;
real g,geff,aj_piby2,Iclk;

integer out[0:1],out_prev[0:1];
integer value_porta[0:1],value_portb[0:1],value_portc[0:1];
integer switch_out,switch_prev_out,clocked;



// This function computes the critical switching current for a multi-layer stack with the given device specifications.


analog function real critical_switching_current;

input Volume,I_rfrp;
real Volume,I_rfrp,eta;
begin
if(I_rfrp < 0) begin
eta = P/(1+pow(P,2));
critical_switching_current = (2*`P_Q/(`P_H/(2*`M_PI)))*(alpha*Ms*Volume/eta)*`P_U0*(Ms/(2*sqrt(2)));
end
else begin
eta = P/(1-pow(P,2));
critical_switching_current = (2*`P_Q/(`P_H/(2*`M_PI)))*(alpha*Ms*Volume/eta)*`P_U0*(Ms/(2*sqrt(2)));
end
end
endfunction



// This function computes the angle of magnetization of the free layer of the multi-layer device as function of the current through it.


analog function real angle;

input aj_piby2;
real aj_piby2,cos_angle;
begin
cos_angle = (Hd*gamma_p*(4*`M_PI/1e7) - aj_piby2/alpha)/((4*`M_PI*Ms + 0.5*Hk)*(4*`M_PI/1e7));
angle = acos(cos_angle);
end
endfunction



// This function computes the resistance of the multi-layer stack as a function of the angle of magnetization of the free layer.


analog function real resistance;

input theta;
real theta,G_theta;
begin
G_theta = 0.5*((1+cos(theta))*Gp + (1-cos(theta))*Gap);
resistance = 1/G_theta;
end
endfunction


// This function computes the magnitude of the clocking current that is required to align the magnetization of the free layer along the y-axis into a
// stationary magnetization state.


analog function real clocking_current;

input Volume;
real Volume;
real Hdz,G,Jp;
begin
Hdz = Hd*cos(theta0);
G = 1/(-4 + (pow((1+P),3))*3/(4*pow(P,1.5)));
Jp = `P_U0*(pow(Ms,2))*(`P_Q)*Volume/((`P_H)/(2*`M_PI));

clocking_current = Hdz*Jp/(G*(sin(theta0))*Ms);
end
endfunction


// This function is used to compute the post-clocking neighbor interaction among the free layers of the neighboring multi-layer cells.
// The neighbor interaction is implemented as a Finite State Machine.


analog function real mtj_fsm_horz_cell;

input [0:1] value_porta,value_portb,value_portc;
integer value_porta[0:1],value_portb[0:1],value_portc[0:1];
integer portb_and,porta_and,portb_xor_porta;
integer bitwise_xor[0:2];

integer out1[0:1],out2[0:1],out3[0:1],out[0:1];
begin
bitwise_xor[0] = (value_porta[0] ^ value_porta[1]);
bitwise_xor[1] = (value_portb[0] ^ value_portb[1]);
bitwise_xor[2] = (value_portc[0] ^ value_portc[1]);

portb_and = (value_portb[0] && value_portb[1]);
porta_and = (value_porta[0] && value_porta[1]);
portb_xor_porta = portb_and ^ porta_and;

out1[0] = value_portc[0] && value_portb[0];
out1[1] = value_portc[1] && value_portb[1];

out2[0] = (value_portb[0] ^ value_portc[0]);
out2[1] = (value_portb[1] ^ value_portc[1]);
out3[0] = !value_porta[0];
out3[1] = !value_porta[1];

out2[0] = out2[0] && out3[0];
out2[1] = out2[1] && out3[1];
out1[0] = out1[0] || out2[0];
out1[1] = out1[1] || out2[1];

if(bitwise_xor[0]) begin
if(bitwise_xor[1] && bitwise_xor[2]) begin
out[0] = low;
out[1] = high;
end
else if(bitwise_xor[1]) begin
out[0] = value_portc[0];
out[1] = value_portc[1];
end
else if(bitwise_xor[2]) begin
out[0] = !value_portb[0];
out[1] = !value_portb[1];
end
else begin
out[0] = value_portc[0];
out[1] = value_portc[1];
end
end
else begin
if(bitwise_xor[1] && bitwise_xor[2]) begin
out[0] = !(value_porta[0]);
out[1] = !(value_porta[1]);
end
else if(bitwise_xor[1])
begin
out[0] = value_portc[0];
out[1] = value_portc[1];
end
else if(bitwise_xor[2]) begin
if(portb_xor_porta) begin
out[0] = low;
out[1] = high;
end
else begin
out[0] = !(value_porta[0]);
out[1] = !(value_porta[1]);
end
end
else begin
out[0] = out1[0];
out[1] = out1[1];
end
end

if(out[0] ^ out[1])
mtj_fsm_horz_cell = 0.5;
else begin
if(out[0] == 0)
mtj_fsm_horz_cell = 0.0;
else
mtj_fsm_horz_cell = 1.0;
end
end
endfunction


// This function is used to decide the specific writing or clocking or no operation that needs to be performed on the cell.
// The specific operation that needs to be performed depends on the magnitude of the current through the cell and the duration of the current pulse.


analog function integer switch_output;

input [0:1]out_prev,switch_prev_out;
input I_switch, Iclk, t_precession, tpulse, I_max;
integer out_prev[0:1],switch_prev_out;
real I_switch,Iclk,t_precession,tpulse,I_max;
integer out[0:1];
real Imax_abs,Iswitch_abs;

begin
Imax_abs = abs(I_max);
Iswitch_abs = abs(I_switch);

if(abs(Imax_abs - Iclk) < I_delta) begin
if(I_max > 0) begin
out[0] = high;
out[1] = low;
switch_output = 01;
end
else begin
out[0] = out_prev[0];
out[1] = out_prev[1];
switch_output = 10;
end
end
else if((Imax_abs > Iswitch_abs) && (tpulse > t_precession/2)) begin
if(I_max < 0) begin
out[0] = high;
out[1] = high;
switch_output = 11;
end
else begin
out[0] = low;
out[1] = low;
switch_output = 00;
end
end
else begin
out[0] = out_prev[0];
out[1] = out_prev[1];
switch_output = switch_prev_out;
end
end
endfunction


// Main body of the program.


analog
begin
Volume = Length*Width*thickness;
t_precession = abs(`M_PI/`gamma);
geff = sqrt(P) + 1/sqrt(P);
g = 1/(-4 + (pow(geff,3))*0.75);
aj_piby2 = ((`P_H/(2*`M_PI))/(2*`P_Q))*g*I(rf,rp)/(Ms*Volume);
theta = angle(aj_piby2);
Resistance = resistance(theta) ;
I_switch = critical_switching_current(Volume,I(rf,rp));
Iclk = clocking_current(Volume);

@(initial_step) begin
trise_pos = 0;
tfall_pos = 0;
tfall_neg = 0;
trise_neg = 0;
tswitch = 0;
out[0] = high;
out[1] = high;
out_prev[0] = out[0];
out_prev[1] = out[1];
clocked = 0;
Resistance = R0;
I_prev = 0;
t_prev_pulse = 0;
switch_prev_out = 10;
end

@(cross(V(rf,rp) - V_delta,+1)) begin
trise_pos = $abstime;
end

@(cross(V(rf,rp) + V_delta,-1)) begin
tfall_neg = $abstime;
end

@(cross(V(rf,rp) - V_delta,-1)) begin
tfall_pos = $abstime;
t_prev_pulse = tpulse;
tpulse = tfall_pos - trise_pos;
tswitch = tfall_pos + t_delta;
end

@(cross(V(rf,rp) + V_delta,+1)) begin
trise_neg = $abstime;
t_prev_pulse = tpulse;
tpulse = trise_neg - tfall_neg;
tswitch = trise_neg + t_delta;
end

if (tpulse < 1.5e-12) begin
tpulse = t_prev_pulse;
end


@(timer(tswitch)) begin

switch_out = switch_output(out_prev,switch_prev_out,I_switch,Iclk,t_precession,tpulse,I_max);

switch_prev_out = switch_out;
case (switch_out)
00 : begin out[0] = low; out[1] = low; clocked=0; end
11 : begin out[0] = high; out[1] = high; clocked=0; end
01 : begin out[0] = high; out[1] = low; clocked=1; end
default : begin out[0] = out_prev[0]; out[1] = out_prev[1]; clocked=0; end
endcase

end


V(int_nodea,int_nodeb) <+ ddt(V(rf,rp));
if((abs(V(int_nodea,int_nodeb)) < 2e10) && ((abs(I(rf,rp))>I_upper) || (abs(I(rf,rp)-Iclk)
I_max = I(rf,rp);
I_prev = I_max;
end
else begin
I_max = I_prev;
end


value_porta[0] = (H(mag_porta[0]) > 0);
value_porta[1] = (H(mag_porta[1]) > 0);
value_portb[0] = (H(mag_portb[0]) > 0);
value_portb[1] = (H(mag_portb[1]) > 0);
value_portc[0] = (H(mag_portc[0]) > 0);
value_portc[1] = (H(mag_portc[1]) > 0);


H(mag_portd[0]) <+ transition(out[0],tdelay,trise,tfall);
H(mag_portd[1]) <+ transition(out[1],tdelay,trise,tfall);




if(clocked == 1) begin

fsm_out = mtj_fsm_horz_cell(value_porta,value_portb,value_portc);
if((fsm_out - 0.5) == 0) begin
out[0] = low;
out[1] = high;
end
else if(fsm_out == 1) begin
out[0] = high;
out[1] = high;
end
else begin
out[0] = low;
out[1] = low;
end
clocked = 0;
end
out_prev[0] = out[0];
out_prev[1] = out[1];


if(out[0] ^ out[1]) begin
Resistance = resistance(theta);
end
else begin
if(out[0] == 0) begin
Resistance = R0;
end
else begin
Resistance = R1;
end
end
I(rf,rp) <+ V(rf,rp)/Resistance;

end
endmodule