A more serious implementation

This commit is contained in:
dwarning 2009-01-02 10:35:57 +00:00
parent 11679e737b
commit 760b55367a
2 changed files with 636 additions and 86 deletions

View File

@ -1,3 +1,8 @@
2009-01-02 Dietmar Warning
* adms/ekv/amsva/ekv.va: EPFL-EKV version 2.63, replacement of the long channel
version with a code according to the official manual (revision II) available
at http://legwww.epfl.ch/ekv, contribution of Ivan Riis Nielsen 11/2006.
2009-01-01 Dietmar Warning
* configure.in, include/missing_math.h, src/math/misc/isnan.c: POSIX conform
configure isnan, isinf macros and finite function

View File

@ -1,110 +1,655 @@
// EPFL-EKV version 2.6: A Verilog-A description.
// The intrinsic device is coded according to the official manual
// (revision II) available at http://legwww.epfl.ch/ekv.
// contribution of Ivan Riis Nielsen 11/2006, modified by Dietmar Warning 01/2009
//`include "std.va"
//`include "const.va"
//Spice
//Default simulator: Spectre
`ifdef insideADMS
`define P(txt) (*txt*)
`define PGIVEN(p) $given(p)
`define INITIAL_MODEL @(initial_model)
`define INSTANCE @(initial_instance)
`define NOISE @(noise)
`else
`define P(txt) (txt)
`define PGIVEN(p) p
`define INITIAL_MODEL
`define INSTANCE
`define NOISE
`endif
//ADS
//`include "constants.vams"
//`include "disciplines.vams"
//`include "compact.vams"
//Spectre
`include "constants.h"
`include "discipline.h"
`ifdef insideADMS
`define P(p) (*p*)
`define PGIVEN(p) $given(p)
`define INITIAL_MODEL @(initial_model)
`else
`define P(p)
`define PGIVEN(p) p
`define INITIAL_MODEL @(initial_step)
`endif
`define NMOS 1
`define PMOS -1
//dw
`define TMAX 326.85
`define TMIN -100.0
`define EPSSI `P_EPS0*11.7
`define EPSOX `P_EPS0*3.9
`define TREF 300.15
`define SQR(x) ((x)*(x))
`define VT(temp) (`P_K*temp/`P_Q)
`define EG(temp) (1.16-0.000702*`SQR(temp)/(temp+1108))
`define NI(temp) (1.45e16*(temp/`TREF)*exp(`EG(`TREF)/(2*`VT(`TREF))-`EG(temp)/(2*`VT(temp))))
// ****************************************************************
// * EKV MOS model (long channel) based on version 2.6 rev.15
// * Function name : ekv26_dc_long for verilog-a implementation
// * The model documetation: http://legwww.epfl.ch/ekv
// ****************************************************************
`define oneThird 3.3333333333333333e-01
module ekv(d,g,s,b);
//
// Node definitions
//
inout d,g,s,b ; // external nodes
electrical d,g,s,b ; // external nodes
// Constants needed in safe exponential function (called "expl")
`define se05 2.3025850929940458e+02
`define ke05 1.0e-100
`define ke05inv 1.0e100
// P3 3rd order polynomial expansion of exp()
`define P3(u) (1.0 + (u) * (1.0 + 0.5 * ((u) * (1.0 + (u) * `oneThird))))
// expl exp() with 3rd order polynomial extrapolation
// to avoid overflows and underflows and retain C-3 continuity
`define expl(x, res) \
if (abs(x) < `se05) begin\
res = exp(x); \
end else begin \
if ((x) < -`se05) begin\
res = `ke05 / `P3(-`se05 - (x)); \
end else begin\
res = `ke05inv * `P3((x) - `se05); \
end \
end
//dw
real VT,Tamb,Tdev,Tnom,dT;
parameter real tnom = 27 `P(spice:name="tnom" info="Temperature for which parameters are valid" unit="C");
parameter real dt = 0.0 `P(spice:name="dt" type="instance" info="Temperature change for particular transistor" unit="K");
module ekv (d,g,s,b);
//
//*** Local variables
//
real x, VG, VS, VD, VGprime, VP;
real beta, n, iff, ir, Ispec, Id;
//
//*** model parameter definitions
//
parameter real L = 10E-6 from[0.0:inf];
parameter real W = 10E-6 from[0.0:inf];
// Node definitions
//*** Threshold voltage
// substrate effect parameters (long-channel)
parameter real VTO = 0.5 from[0.0:inf];
parameter real GAMMA = 0.7 from[0.0:inf];
parameter real PHI = 0.5 from[0.2:inf];
inout d,g,s,b;
electrical d,g,s,b,di,si;
//*** Mobility parameters (long-channel)
parameter real KP = 20E-6 from[0.0:inf];
parameter real THETA = 50.0E-3 from[0.0:inf];
analog begin // EKV v2.6 long-channel
//dw
Tnom = tnom+273.15;
Tamb = $temperature;
Tdev = Tamb+dt; // selfheating instead dT later possible
// Limit temperature to avoid FPE's in equations
if(Tdev < `TMIN + 273.15)
Tdev = `TMIN + 273.15;
else
if (Tdev > `TMAX + 273.15)
Tdev = `TMAX + 273.15;
VT = `P_K*Tdev /`P_Q;
// Model parameters
parameter integer nmos=1 from [0:1] `P(info="MOS type : nmos:0");
parameter integer pmos=1 from [0:1] `P(info="MOS type : pmos:0");
parameter integer MTYPE=(nmos==0 ? (pmos==0 ? 0 : 1) : (pmos==0 ? -1 : 1));
parameter real TNOM=27 from (-273.15:inf)
`P(info="Nominal temperature [degC]");
parameter real IMAX=1 from (0:inf)
`P(info="Maximum forward junction current before linearization [A]");
// - intrinsic model (optional, section 4.2.1)
parameter real TOX=0 from [0:inf)
`P(info="Oxide thickness [m]");
parameter real NSUB=0 from [0:inf)
`P(info="Channel doping [cm^-3]");
parameter real VFB=1001.0 from (-inf:inf) // use 1001V as "not specified"
`P(info="Flat-band voltage [V]");
parameter real UO=0 from [0:inf)
`P(info="Low-field mobility [cm^2/Vs]");
parameter real VMAX=0 from [0:inf)
`P(info="Saturation velocity [m/s]");
parameter real THETA=0 from [0:inf)
`P(info="Mobility reduction coefficient [V^-1]");
// - intrinsic model (process related, section 4.1)
parameter real COX=((TOX>0) ? (`EPSOX/TOX) : 0.7m) from [0:inf)
`P(info="Oxide capacitance [F/m^2]");
parameter real XJ=0.1u from [1n:inf)
`P(info="Junction depth [m]");
parameter real DL=0 from (-inf:inf)
`P(info="Length correction [m]");
parameter real DW=0 from (-inf:inf)
`P(info="Width correction [m]");
// - intrinsic model (basic, section 4.2)
parameter real GAMMA=((NSUB>0) ? (sqrt(2*`P_Q*`EPSSI*NSUB*1e6)/COX) : 1) from [0:inf)
`P(info="Body effect parameter [V^0.5]");
parameter real PHI=((NSUB>0) ? (2*`VT((TNOM+273.15))*ln(max(NSUB,1)*1e6/`NI((TNOM+273.15)))) : 0.7) from [0.1:inf)
`P(info="Bulk Fermi potential (*2) [V]");
parameter real VTO=((VFB<1000.0) ? (VFB+MTYPE*(PHI+GAMMA*sqrt(PHI))) : 0.5) from (-inf:inf)
`P(info="Long-channel threshold voltage [V]");
parameter real KP=((UO>0) ? (UO*1e-4*COX) : 50u) from (0:inf)
`P(info="Transconductance parameter [A/V^2]");
parameter real UCRIT=(((VMAX>0) && (UO>0)) ? (VMAX/(UO*1e-4)) : 2e6 ) from [100k:inf)
`P(info="Longitudinal critical field [V/m]");
parameter real E0=((THETA>0) ? 0 : 1e12) from [100k:inf)
`P(info="Mobility reduction coefficient [V/m]");
// - intrinsic model (channel length modulation and charge sharing, section 4.3)
parameter real LAMBDA=0.5 from [0:inf)
`P(info="Depletion length coefficient (CLM)");
parameter real WETA=0.25 from (-inf:inf)
`P(info="Narrow-channel effect coefficient");
parameter real LETA=0.1 from (-inf:inf)
`P(info="Short-channel effect coefficient");
// - intrinsic model (reverse short channel effect, section 4.4)
parameter real Q0=0 from (-inf:inf)
`P(info="RSCE peak charge density [C/m^2]");
parameter real LK=0.29u from [10n:inf)
`P(info="RSCE characteristic length [m]");
// - intrinsic model (impact ionization, section 4.5)
parameter real IBA=0 from (-inf:inf)
`P(info="First impact ionization coefficient [m^-1]");
parameter real IBB=3e8 from [1e8:inf)
`P(info="Second impact ionization coefficient [V/m]");
parameter real IBN=1 from [0.1:inf)
`P(info="Saturation voltage factor for impact ionization");
// - intrinsic model (temperature, section 4.6)
parameter real TCV=1m from (-inf:inf)
`P(info="Threshold voltage TC [V/K]");
parameter real BEX=-1.5 from (-inf:inf)
`P(info="Mobility temperature exponent");
parameter real UCEX=0.8 from (-inf:inf)
`P(info="Longitudinal critical field temperature exponent");
parameter real IBBT=9e-4 from (-inf:inf)
`P(info="Temperature coefficient for IBB [K^-1]");
// - intrinsic model (matching, section 4.7)
parameter real AVTO=0 from (-inf:inf)
`P(info="Area related VTO mismatch parameter [Vm]");
parameter real AKP=0 from (-inf:inf)
`P(info="Area related KP mismatch parameter [m]");
parameter real AGAMMA=0 from (-inf:inf)
`P(info="Area related GAMMA mismatch parameter [V^0.5*m]");
// - intrinsic model (flicker noise, section 4.8)
parameter real KF=0 from [0:inf)
`P(info="Flicker noise coefficient");
parameter real AF=1 from (-inf:inf)
`P(info="Flicker noise exponent");
// - intrinsic model (setup, section 4.9)
parameter real NQS=0 from [0:1]
`P(info="Non-quasi-static operation switch");
parameter real SATLIM=exp(4) from (0:inf)
`P(info="Saturation limit (if/ir)");
parameter real XQC=0.4 from [0:1]
`P(info="Charge/capacitance model selector");
// - external parasitic parameters
parameter real HDIF=0 from [0:inf)
`P(info="S/D diffusion length (/2) [m]");
parameter real RSH=0 from [0:inf)
`P(info="S/D sheet resistance [ohm]");
parameter real JS=0 from [0:inf)
`P(info="S/D junction saturation current density [A/m^2]");
parameter real JSW=0 from [0:inf)
`P(info="S/D junction sidewall saturation current density [A/m]");
parameter real XTI=3 from [0:inf)
`P(info="S/D diode saturation current temperature exponent");
parameter real N=1 from [0.5:10]
`P(info="S/D diode emission coefficient");
parameter real CJ=0 from [0:inf)
`P(info="S/D zero-bias junction capacitance per area [F/m^2]");
parameter real CJSW=0 from [0:inf)
`P(info="S/D zero-bias junction capacitance per perimeter [F/m]");
parameter real PB=0.8 from (0:inf)
`P(info="S/D bottom junction builtin potential [V]");
parameter real PBSW=PB from (0:inf)
`P(info="S/D sidewall junction builtin potential [V]");
parameter real MJ=0.5 from (0:inf)
`P(info="S/D bottom junction grading coefficient");
parameter real MJSW=0.333 from (0:inf)
`P(info="S/D sidewall junction grading coefficient");
parameter real FC=0.5 from (0:inf)
`P(info="S/D bottom junction forward-bias threshold");
parameter real FCSW=FC from (0:inf)
`P(info="S/D sidewall junction forward-bias threshold");
parameter real CGSO=0 from [0:inf)
`P(info="Gate-source overlap capacitance per width [F/m]");
parameter real CGDO=0 from [0:inf)
`P(info="Gate-drain overlap capacitance per width [F/m]");
parameter real CGBO=0 from [0:inf)
`P(info="Gate-bulk overlap capacitance per length [F/m]");
VG = V(g); VS = V(s); VD = V(d);
// Instance parameters
// Effective gate voltage (33)
VGprime = VG - VTO + PHI + GAMMA * sqrt(PHI);
// - intrinsic model
parameter real L=10u from [0:inf]
`P(type="instance" info="Drawn length [m]" unit="m");
parameter real W=10u from [0:inf]
`P(type="instance" info="Drawn width [m]" unit="m");
parameter real M=1 from [0:inf]
`P(type="instance" info="Parallel multiplier" unit="m");
// parameter real N=1 from [0:inf]
// `P(type="instance" info="Series multiplier" unit="m");
// Pinch-off voltage (34)
VP = VGprime - PHI - GAMMA
* (sqrt(VGprime+(GAMMA/2.0)*(GAMMA/2.0))-(GAMMA/2.0));
// - external parasitics
parameter real AD=((HDIF>0) ? (2*HDIF*W) : 0) from [0:inf)
`P(info="Drain area [m^2]" type="instance");
parameter real AS=((HDIF>0) ? (2*HDIF*W) : 0) from [0:inf)
`P(info="Source area [m^2]" type="instance");
parameter real PD=((HDIF>0) ? (4*HDIF+2*W) : 0) from [0:inf)
`P(info="Drain perimeter [m]" type="instance");
parameter real PS=((HDIF>0) ? (4*HDIF+2*W) : 0) from [0:inf)
`P(info="Source perimeter [m]" type="instance");
parameter real NRD=((HDIF>0) ? (HDIF/W) : 0) from [0:inf)
`P(info="Drain no. squares" type="instance");
parameter real NRS=((HDIF>0) ? (HDIF/W) : 0) from [0:inf)
`P(info="Source no. squares" type="instance");
parameter real RS=((RSH>0) ? (RSH*NRS) : 0) from [0:inf)
`P(info="Source resistance [ohms]" type="instance");
parameter real RD=((RSH>0) ? (RSH*NRD) : 0) from [0:inf)
`P(info="Drain resistance [ohms]" type="instance");
// Slope factor (39)
n = 1.0 + GAMMA / (2.0*sqrt(PHI + VP + 4.0*VT));
// Mobility equation (58), (64)
beta = KP * (W/L) * (1.0/(1.0 + THETA * VP));
// Declaration of variables
integer mode;
real lc,isat_s,vexp_s,gexp_s,isat_d,vexp_d,gexp_d,fact,
weff,leff,np,ns,lmin,rd,rs,ceps,ca,xsi,dvrsce,
tempk,vt,sqrt_A,vto_a,kp_a,gamma_a,ucrit,phi,ibb,vc,qb0,
vg,vd,vs,tmp,vgprime,vp0,vsprime,vdprime,gamma0,gammaprime,vp,n,ifwd,
vdss,vdssprime,dv,vds,vip,dl,lprime,leq,irprime,irev,beta0,nau,
nq,xf,xr,qd,qs,qi,qb,qg,beta0prime,beta,vpprime,is,ids,vib,
idb,ibdj,ibsj,coxt,qdt,qst,qgt,qbt,
cbs0,cbs0sw,cbs,cbd0,cbd0sw,cbd,
fv,z0,z1,y;
// forward (44) and reverse (56) currents
x=(VP-VS)/VT; iff = (ln(1.0+exp( x /2.0)))*(ln(1.0+exp( x /2.0)));
x=(VP-VD)/VT; ir = (ln(1.0+exp( x /2.0)))*(ln(1.0+exp( x /2.0)));
real cgso,cgdo,cgbo;
// Specific current (65)
Ispec = 2 * n * beta * VT * VT;
analog begin
// Drain current (66)
Id = Ispec * (iff - ir);
`INITIAL_MODEL begin // Model Initialization
//
// Branch contributions to EKV v2.6 model (long-channel)
//
I(d,s) <+ Id;
lc = sqrt(`EPSSI/COX*XJ);
end // INITIAL_MODEL
end // analog
endmodule
`INSTANCE begin // temperature independent device initialization
weff = W+DW;
leff = L+DL;
np = M;
ns = 1;
// eq. 54
lmin = 0.1*ns*leff;
rs = ns/np*RS;
rd = ns/np*RD;
ceps = 4*22e-3*22e-3;
ca = 0.028;
xsi = ca*(10*leff/LK-1);
dvrsce = 2*Q0/COX/`SQR(1+0.5*(xsi+sqrt(xsi*xsi+ceps)));
coxt = np*ns*COX*weff*leff;
end // temperature independent
`INSTANCE begin // temperature dependent device initialization
tempk = $temperature;
vt = `VT(tempk);
sqrt_A = sqrt(np*weff*ns*leff);
vto_a = MTYPE*(VTO+TCV*(tempk-(TNOM+273.15)))+AVTO/sqrt_A;
kp_a = KP*pow(tempk/(TNOM+273.15),BEX)*(1+AKP/sqrt_A);
gamma_a = GAMMA+AGAMMA/sqrt_A;
ucrit = UCRIT*pow(tempk/(TNOM+273.15),UCEX);
phi = PHI*tempk/(TNOM+273.15)-3*vt*ln(tempk/(TNOM+273.15))-`EG(TNOM+273.15)*tempk/(TNOM+273.15)+`EG(tempk);
ibb = IBB*(1+IBBT*(tempk-(TNOM+273.15)));
vc = ucrit*ns*leff;
// eq. 60
qb0 = gamma_a*sqrt(phi);
fact = (`EG(TNOM+273.15)/`VT(TNOM+273.15)-`EG(tempk)/vt) * pow(tempk/(TNOM+273.15),XTI);
`expl(fact,tmp);
isat_s = np*ns*(JS*AS+JSW*PS)*tmp;
isat_d = np*ns*(JS*AD+JSW*PD)*tmp;
if (isat_s>0) begin
vexp_s = vt*ln(IMAX/isat_s+1);
gexp_s = (IMAX+isat_s)/vt;
end else begin
vexp_s = -1e9;
gexp_s = 0;
end
if (isat_d>0) begin
vexp_d = vt*ln(IMAX/isat_d+1);
gexp_d = (IMAX+isat_d)/vt;
end else begin
vexp_d = -1e9;
gexp_d = 0;
end
cbs0 = np*ns*CJ*AS;
cbd0 = np*ns*CJ*AD;
cbs0sw = np*ns*CJSW*PS;
cbd0sw = np*ns*CJSW*PD;
cgso = np*ns*CGSO*weff;
cgdo = np*ns*CGDO*weff;
cgbo = np*ns*CGBO*leff;
end // temperature dependent
begin //Bias-dependent model evaluation
vg = MTYPE*V(g,b);
vd = MTYPE*V(di,b);
vs = MTYPE*V(si,b);
// $strobe("vg=%e vd=%e vs=%e",vg,vd,vs);
if (vd>=vs)
mode = 1;
else begin
mode = -1;
tmp = vs;
vs = vd;
vd = tmp;
end
// eq. 33
vgprime = vg-vto_a-dvrsce+phi+gamma_a*sqrt(phi);
// eq. 35
vsprime = 0.5*(vs+phi+sqrt(`SQR(vs+phi)+16*`SQR(vt)));
vdprime = 0.5*(vd+phi+sqrt(`SQR(vd+phi)+16*`SQR(vt)));
// $strobe("vgprime=%e vdprime=%e vsprime=%e",vgprime,vdprime,vsprime);
// eq. 34
if (vgprime>=0) begin
vp0 = vgprime-phi-gamma_a*(sqrt(vgprime+0.25*`SQR(gamma_a))-0.5*gamma_a);
// eq. 36
gamma0 = gamma_a-`EPSSI/COX*(LETA/leff*(sqrt(vsprime)+sqrt(vdprime))-3*WETA/weff*sqrt(vp0+phi));
end else begin
vp0 = -phi;
// eq. 36 - skipped sqrt(vp0+phi) here, it produces inf on derivative
gamma0 = gamma_a-`EPSSI/COX*(LETA/leff*(sqrt(vsprime)+sqrt(vdprime)) );
end
// eq. 37
gammaprime = 0.5*(gamma0+sqrt(`SQR(gamma0)+0.1*vt));
// eq. 38
if (vgprime>=0)
vp = vgprime-phi-gammaprime*(sqrt(vgprime+0.25*`SQR(gammaprime))-0.5*gammaprime);
else
vp = -phi;
// $strobe("vp0=%e vp=%e gamma0=%e gammaprime=%e",vp0,vp,gamma0,gammaprime);
// eq. 39
n = 1+gamma_a*0.5/sqrt(vp+phi+4*vt);
// Forward current (43-44)
fv=(vp-vs)/vt;
if (fv >= -0.35)
z0=2.0/(1.3 + fv - ln(fv+1.6));
if (fv>=-15 && fv<-0.35) begin
`expl(-fv,tmp);
z0= 1.55 + tmp;
end else
z0=1;
z1=(2.0 + z0) / (1.0 + fv + ln(z0));
if (fv > -15.0)
y=(1.0 + fv + ln(z1)) / (2.0 + z1);
else begin
`expl(-fv,tmp);
y= 1.0 / (2.0 + tmp);
end
ifwd = y*(1.0 + y);
z0 = 1;
z1 = 1;
// eq. 46
vdss = vc*(sqrt(0.25+vt/vc*sqrt(ifwd))-0.5);
// eq. 47
vdssprime = vc*(sqrt(0.25+vt/vc*(sqrt(ifwd)-0.75*ln(ifwd)))-0.5)+vt*(ln(0.5*vc/vt)-0.6);
// $strobe("ifwd=%e vdss=%e vdssprime=%e",ifwd,vdss,vdssprime);
// eq. 48
dv = 4*vt*sqrt(LAMBDA*(sqrt(ifwd)-vdss/vt)+1.0/64);
// eq. 49
vds = 0.5*(vd-vs);
// eq. 50
vip = sqrt(`SQR(vdss)+`SQR(dv))-sqrt(`SQR(vds-vdss)+`SQR(dv));
// eq. 52
dl = LAMBDA*lc*ln(1+(vds-vip)/(lc*ucrit));
// eq. 53
lprime = ns*leff-dl+(vds+vip)/ucrit;
// eq. 55
leq = 0.5*(lprime+sqrt(`SQR(lprime)+`SQR(lmin)));
// eq. 56
fv=(vp-vds-vs-sqrt(`SQR(vdssprime)+`SQR(dv))+sqrt(`SQR(vds-vdssprime)+`SQR(dv)))/vt;
if (fv >= -0.35)
z0=2.0/(1.3 + fv - ln(fv+1.6));
if (fv>=-15 && fv<-0.35) begin
`expl(-fv,tmp);
z0= 1.55 + tmp;
end else
z0=1;
z1=(2.0 + z0) / (1.0 + fv + ln(z0));
if (fv > -15.0)
y=(1.0 + fv + ln(z1)) / (2.0 + z1);
else begin
`expl(-fv,tmp);
y= 1.0 / (2.0 + tmp);
end
irprime = y*(1.0 + y);
z0 = 1;
z1 = 1;
// eq. 57
fv=(vp-vd)/vt;
if (fv >= -0.35)
z0=2.0/(1.3 + fv - ln(fv+1.6));
if (fv>=-15 && fv<-0.35) begin
`expl(-fv,tmp);
z0= 1.55 + tmp;
end else
z0=1;
z1=(2.0 + z0) / (1.0 + fv + ln(z0));
if (fv > -15.0)
y=(1.0 + fv + ln(z1)) / (2.0 + z1);
else begin
`expl(-fv,tmp);
y= 1.0 / (2.0 + tmp);
end
irev = y*(1.0 + y);
// eq. 58
beta0 = kp_a*np*weff/leq;
// eq. 59
nau = (5+MTYPE)/12.0;
// eq. 69
nq = 1+0.5*gamma_a/sqrt(vp+phi+1e-6);
// eq. 70
xf = sqrt(0.25+ifwd);
// eq. 71
xr = sqrt(0.25+irev);
// eq. 72
qd = -nq*(4.0/15*(3*`SQR(xr)*(xr+2*xf)+2*`SQR(xf)*(xf+2*xr))/`SQR(xf+xr)-0.5);
// eq. 73
qs = -nq*(4.0/15*(3*`SQR(xf)*(xf+2*xr)+2*`SQR(xr)*(xr+2*xf))/`SQR(xf+xr)-0.5);
// eq. 74
qi = qs+qd;
// eq. 75
if (vgprime>=0)
qb = (-gamma_a*sqrt(vp+phi+1e-6))/vt-(nq-1)/nq*qi;
else
qb = -vgprime/vt;
// eq. 76 (qox removed since it is assumed to be zero)
qg = -qi-qb;
if (E0!=0) begin
// eq. 61
beta0prime = beta0*(1+COX/(E0*`EPSSI)*qb0);
// eq. 62
beta = beta0prime/(1+COX/(E0*`EPSSI)*vt*abs(qb+nau*qi));
end else begin
// eq. 63
vpprime = 0.5*(vp+sqrt(`SQR(vp)+2*`SQR(vt)));
// eq. 64
beta = beta0/(1+THETA*vpprime);
end // else: !if(e0!=0)
// eq. 65
is = 2*n*beta*`SQR(vt);
// $strobe("beta0=%e beta0prime=%e beta=%e E0=%e qb0=%e qb=%e qi=%e",beta0,beta0prime,beta,E0,qb0,qb,qi);
// eq. 66
ids = is*(ifwd-irprime);
// eq. 67
vib = vd-vs-IBN*2*vdss;
// eq. 68
if (vib>0) begin
`expl((-ibb*lc)/vib,tmp)
idb = ids*IBA/ibb*vib*tmp;
end else
idb = 0;
// $strobe("ids=%e idb=%e",ids,idb);
if (mode>1) begin
if (isat_s>0) begin
if (-vs>vexp_s)
ibsj = IMAX+gexp_s*(-vs-vexp_s);
else begin
`expl(-vs/vt,tmp);
ibsj = isat_s*(tmp-1);
end
end else
ibsj = 0;
if (isat_d>0) begin
if (-vd>vexp_d)
ibdj = IMAX+gexp_d*(-vd-vexp_d);
else begin
`expl(-vd/vt,tmp);
ibdj = isat_d*(tmp-1);
end
end else
ibdj = 0;
end else begin // if (mode>1)
if (isat_s>0) begin
if (-vd>vexp_s)
ibsj = IMAX+gexp_s*(-vd-vexp_s);
else begin
`expl(-vd/vt,tmp);
ibsj = isat_s*(tmp-1);
end
end else
ibsj = 0;
if (isat_d>0) begin
if (-vs>vexp_d)
ibdj = IMAX+gexp_d*(-vs-vexp_d);
else begin
`expl(-vs/vt,tmp);
ibdj = isat_d*(tmp-1);
end
end else
ibdj = 0;
end // else: !if(mode>1)
qdt = coxt*vt*qd;
qst = coxt*vt*qs;
qgt = coxt*vt*qg;
qbt = coxt*vt*qb;
cbs = 0;
cbd = 0;
if (cbs0>0) begin
if (MTYPE*V(b,si)>FC*PB)
cbs = cbs+cbs0/pow(1-FC,MJ)*(1+MJ*(MTYPE*V(b,si)-PB*FC))/(PB*(1-FC));
else
cbs = cbs+cbs0/pow(1-MTYPE*V(b,si),MJ);
end
if (cbd0>0) begin
if (MTYPE*V(b,di)>FC*PB)
cbd = cbd+cbd0/pow(1-FC,MJ)*(1+MJ*(MTYPE*V(b,di)-PB*FC))/(PB*(1-FC));
else
cbd = cbd+cbd0/pow(1-MTYPE*V(b,di),MJ);
end
if (cbs0sw>0) begin
if (MTYPE*V(b,si)>FCSW*PBSW)
cbs = cbs+cbs0sw/pow(1-FCSW,MJSW)*(1+MJSW*(MTYPE*V(b,si)-PBSW*FCSW))/(PBSW*(1-FCSW));
else
cbs = cbs+cbs0sw/pow(1-MTYPE*V(b,si),MJSW);
end
if (cbd0sw>0) begin
if (MTYPE*V(b,di)>FCSW*PBSW)
cbd = cbd+cbd0sw/pow(1-FCSW,MJSW)*(1+MJSW*(MTYPE*V(b,di)-PBSW*FCSW))/(PBSW*(1-FCSW));
else
cbd = cbd+cbd0sw/pow(1-MTYPE*V(b,di),MJSW);
end
end //Bias-dependent model evaluation
begin //Define branch sources
I(di,si) <+ MTYPE*mode*ids;
if (mode>0) begin
I(di,b) <+ MTYPE*idb;
I(di,g) <+ MTYPE*ddt(qdt);
I(si,g) <+ MTYPE*ddt(qst);
end else begin
I(si,b) <+ MTYPE*idb;
I(si,g) <+ MTYPE*ddt(qdt);
I(di,g) <+ MTYPE*ddt(qst);
end // else: !if(mode>0)
I(b,si) <+ MTYPE*ibsj;
I(b,di) <+ MTYPE*ibdj;
I(b,g) <+ MTYPE*ddt(qbt);
I(g,si) <+ cgso*ddt(V(g,si));
I(g,di) <+ cgdo*ddt(V(g,di));
I(g,b) <+ cgbo*ddt(V(g,b));
if (RD>0)
I(d,di) <+ V(d,di)/rd;
else
V(d,di) <+ 0.0;
if (RS>0)
I(s,si) <+ V(s,si)/rs;
else
V(s,si) <+ 0.0;
I(b,si) <+ cbs*ddt(V(b,si));
I(b,di) <+ cbd*ddt(V(b,di));
end // begin
// `NOISE begin //Define noise sources
//
// end // noise
end //analog
endmodule