![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to acans.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [CSRG BSD Unix] / 3BSD / cmd / spice|
Return to acans.f CVS log | Up to [CSRG BSD Unix] / 3BSD / cmd / spice |
1.1 root 1: subroutine acan
2: implicit double precision (a-h,o-z)
3: c
4: c
5: c this routine drives the small-signal analyses.
6: c
7: common /tabinf/ ielmnt,isbckt,nsbckt,iunsat,nunsat,itemps,numtem,
8: 1 isens,nsens,ifour,nfour,ifield,icode,idelim,icolum,insize,
9: 2 junode,lsbkpt,numbkp,iorder,jmnode,iur,iuc,ilc,ilr,numoff,isr,
10: 3 nmoffc,iseq,iseq1,neqn,nodevs,ndiag,iswap,iequa,macins,lvnim1,
11: 4 lx0,lvn,lynl,lyu,lyl,lx1,lx2,lx3,lx4,lx5,lx6,lx7,ld0,ld1,ltd,
12: 5 imynl,imvn,lcvn,loutpt,nsnod,nsmat,nsval,icnod,icmat,icval
13: common /miscel/ atime,aprog(3),adate,atitle(10),defl,defw,defad,
14: 1 defas,rstats(50),iwidth,lwidth,nopage
15: common /cirdat/ locate(50),jelcnt(50),nunods,ncnods,numnod,nstop,
16: 1 nut,nlt,nxtrm,ndist,ntlin,ibr,numvs
17: common /status/ omega,time,delta,delold(7),ag(7),vt,xni,egfet,
18: 1 xmu,mode,modedc,icalc,initf,method,iord,maxord,noncon,iterno,
19: 2 itemno,nosolv,ipostp,iscrch
20: common /flags/ iprnta,iprntl,iprntm,iprntn,iprnto,limtim,limpts,
21: 1 lvlcod,lvltim,itl1,itl2,itl3,itl4,itl5,igoof,nogo,keof
22: common /knstnt/ twopi,xlog2,xlog10,root2,rad,boltz,charge,ctok,
23: 1 gmin,reltol,abstol,vntol,trtol,chgtol,eps0,epssil,epsox
24: common /ac/ fstart,fstop,fincr,skw2,refprl,spw2,jacflg,idfreq,
25: 1 inoise,nosprt,nosout,nosin,idist,idprt
26: common /cje/ maxtim,itime,icost
27: common /blank/ value(1000)
28: integer nodplc(64)
29: complex*16 cvalue(32)
30: equivalence (value(1),nodplc(1),cvalue(1))
31: data aendor /9.87654321d0/
32: call second(t1)
33: c.. post-processor initialization
34: if(ipostp.eq.0) go to 1
35: numcur=jelcnt(9)
36: numpos=nunods+numcur
37: call getm16(ibuff,numpos)
38: numpos=numpos*4
39: if(numcur.eq.0) go to 1
40: loc=locate(9)
41: loccur=nodplc(loc+6)-1
42: c
43: c allocate storage
44: c
45: 1 call getm8(ndiag,2*nstop)
46: call getm8(lvn,nstop)
47: call getm8(imvn,nstop)
48: call getm16(lcvn,nstop)
49: call getm8(lynl,nstop+nut+nlt)
50: call getm8(imynl,nstop+nut+nlt)
51: if (idist.ne.0) call dinit
52: nandd=0
53: if (inoise.eq.0) go to 10
54: if (idist.eq.0) go to 10
55: nandd=1
56: call getm16(lvnctc,nstop)
57: 10 call getm16(loutpt,0)
58: call crunch
59: lyu=lynl+nstop
60: lyl=lyu+nut
61: numout=jelcnt(43)+jelcnt(44)+jelcnt(45)+1
62: c if(ipostp.ne.0) call pheadr(atitle)
63: icalc=0
64: freq=fstart
65: c
66: c load y matrix and c vector, solve for v vector
67: c
68: 100 call getcje
69: if ((maxtim-itime).le.limtim) go to 900
70: omega=twopi*freq
71: call acload
72: call acdcmp
73: call acsol
74: if (igoof.eq.0) go to 200
75: write (6,121) igoof,freq
76: 121 format('0warning: underflow ',i4,' time(s) in ac analysis at freq
77: 1 = ',1pd9.3,' hz')
78: igoof=0
79: c
80: c store outputs
81: c
82: 200 call extmem(loutpt,numout)
83: loco=loutpt+icalc*numout
84: icalc=icalc+1
85: cvalue(loco+1)=dcmplx(freq,omega)
86: loc=locate(43)
87: 310 if (loc.eq.0) go to 350
88: if (nodplc(loc+5).ne.0) go to 320
89: node1=nodplc(loc+2)
90: node2=nodplc(loc+3)
91: iseq=nodplc(loc+4)
92: cvalue(loco+iseq)=cvalue(lcvn+node1)-cvalue(lcvn+node2)
93: loc=nodplc(loc)
94: go to 310
95: 320 iptr=nodplc(loc+2)
96: iptr=nodplc(iptr+6)
97: iseq=nodplc(loc+4)
98: cvalue(loco+iseq)=cvalue(lcvn+iptr)
99: loc=nodplc(loc)
100: go to 310
101: 350 if(ipostp.eq.0) go to 400
102: cvalue(ibuff+1)=dcmplx(freq,0.0d0)
103: call copy16(cvalue(lcvn+2),cvalue(ibuff+2),nunods-1)
104: if(numcur.ne.0) call copy16(cvalue(lcvn+loccur+1),
105: 1 cvalue(ibuff+nunods+1),numcur)
106: call dblsgl(cvalue(ibuff+1),numpos)
107: c call fwrite(cvalue(ibuff+1),numpos)
108: c
109: c noise and distortion analyses
110: c
111: 400 if (nandd.eq.0) go to 410
112: call copy16(cvalue(lcvn+1),cvalue(lvnctc+1),nstop)
113: 410 if (inoise.ne.0) call noise(loco)
114: if (nandd.eq.0) go to 420
115: call copy16(cvalue(lvnctc+1),cvalue(lcvn+1),nstop)
116: 420 if (idist.ne.0) call disto(loco)
117: c
118: c increment frequency
119: c
120: if (icalc.ge.jacflg) go to 1000
121: if (idfreq.ge.3) go to 510
122: freq=freq*fincr
123: go to 100
124: 510 freq=freq+fincr
125: go to 100
126: c
127: c finished
128: c
129: 900 write (6,901)
130: 901 format('0*error*: cpu time limit exceeded ... analysis stopped'/)
131: nogo=1
132: 1000 if(ipostp.eq.0) go to 1010
133: cvalue(ibuff+1)=aendor
134: c call fwrite(cvalue(ibuff+1),numpos)
135: if(ipostp.ne.0) call clrmem(ibuff)
136: 1010 call clrmem(lvnim1)
137: call clrmem(lx0)
138: call clrmem(lvn)
139: call clrmem(imvn)
140: call clrmem(lcvn)
141: call clrmem(imynl)
142: call clrmem(lynl)
143: call clrmem(ndiag)
144: if (idist.eq.0) go to 1020
145: call clrmem(ld0)
146: call clrmem(ld1)
147: 1020 if (nandd.eq.0) go to 1040
148: call clrmem(lvnctc)
149: 1040 call second(t2)
150: rstats(7)=rstats(7)+t2-t1
151: rstats(8)=rstats(8)+icalc
152: return
153: end
154: subroutine cdiv(xr,xi,yr,yi,cr,ci)
155: c.. ok if cr and ci are really xr and xi or yr and yi
156: implicit double precision (a-h,o-z)
157: xrtemp=xr
158: xitemp=xi
159: yrtemp=yr
160: yitemp=yi
161: amag2=yrtemp*yrtemp+yitemp*yitemp
162: cr=(xrtemp*yrtemp+xitemp*yitemp)/amag2
163: ci=(xitemp*yrtemp-xrtemp*yitemp)/amag2
164: return
165: end
166: subroutine cmult(xr,xi,yr,yi,cr,ci)
167: c.. ok if cr and ci are really xr and xi or yr and yi
168: implicit double precision (a-h,o-z)
169: xrtemp=xr
170: xitemp=xi
171: yrtemp=yr
172: yitemp=yi
173: cr=xrtemp*yrtemp-xitemp*yitemp
174: ci=xitemp*yrtemp+xrtemp*yitemp
175: return
176: end
177: subroutine dblsgl(carray,nwds)
178: c.. note that carray really contains complex*16
179: complex*8 carray(1)
180: complex*8 ctemp8(2),ctemp3
181: complex*16 ctemp
182: equivalence (ctemp,ctemp8(1))
183: c
184: c.. gather up numpos 16-bit words from dbl-complex to
185: c.. make up sgl-complex
186: c
187: num=nwds/4
188: do 10 i=1,num
189: ndex=2*i-1
190: ctemp8(1)=carray(ndex)
191: ctemp8(2)=carray(ndex+1)
192: ctemp3=ctemp
193: carray(i)=ctemp3
194: 10 continue
195: return
196: end
197: subroutine acdcmp
198: implicit double precision (a-h,o-z)
199: c
200: c this routine performs an lu factorization of the circuit equation
201: c coefficient matrix.
202: c
203: common /tabinf/ ielmnt,isbckt,nsbckt,iunsat,nunsat,itemps,numtem,
204: 1 isens,nsens,ifour,nfour,ifield,icode,idelim,icolum,insize,
205: 2 junode,lsbkpt,numbkp,iorder,jmnode,iur,iuc,ilc,ilr,numoff,isr,
206: 3 nmoffc,iseq,iseq1,neqn,nodevs,ndiag,iswap,iequa,macins,lvnim1,
207: 4 lx0,lvn,lynl,lyu,lyl,lx1,lx2,lx3,lx4,lx5,lx6,lx7,ld0,ld1,ltd,
208: 5 imynl,imvn,lcvn,loutpt,nsnod,nsmat,nsval,icnod,icmat,icval
209: common /cirdat/ locate(50),jelcnt(50),nunods,ncnods,numnod,nstop,
210: 1 nut,nlt,nxtrm,ndist,ntlin,ibr,numvs
211: common /flags/ iprnta,iprntl,iprntm,iprntn,iprnto,limtim,limpts,
212: 1 lvlcod,lvltim,itl1,itl2,itl3,itl4,itl5,igoof,nogo,keof
213: common /knstnt/ twopi,xlog2,xlog10,root2,rad,boltz,charge,ctok,
214: 1 gmin,reltol,abstol,vntol,trtol,chgtol,eps0,epssil,epsox
215: common /blank/ value(1000)
216: integer nodplc(64)
217: complex*16 cvalue(32)
218: equivalence (value(1),nodplc(1),cvalue(1))
219: c
220: c
221: do 100 i=2,nstop
222: io=nodplc(iorder+i)
223: gdiag=dabs(value(lynl+io))+dabs(value(imynl+io))
224: if (gdiag.ge.gmin) go to 10
225: value(lynl+io)=gmin
226: value(imynl+io)=0.0d0
227: igoof=igoof+1
228: 10 jstart=nodplc(ilc+i)
229: jstop=nodplc(ilc+i+1)-1
230: if (jstart.gt.jstop) go to 100
231: do 90 j=jstart,jstop
232: call cdiv(value(lyl+j),value(imynl+nstop+nut+j),value(lynl+io),
233: 1 value(imynl+io),value(lyl+j),value(imynl+nstop+nut+j))
234: icol=nodplc(ilr+j)
235: kstart=nodplc(iur+i)
236: kstop=nodplc(iur+i+1)-1
237: if (kstart.gt.kstop) go to 90
238: do 80 k=kstart,kstop
239: irow=nodplc(iuc+k)
240: if (icol-irow) 20,60,40
241: c
242: c find (icol,irow) matrix term (upper triangle)
243: c
244: 20 l=nodplc(iur+icol+1)
245: 30 l=l-1
246: if (nodplc(iuc+l).ne.irow) go to 30
247: ispot=lyu+l
248: ispot2=imynl+nstop+l
249: go to 70
250: c
251: c find (icol,irow) matrix term (lower triangle)
252: c
253: 40 l=nodplc(ilc+irow+1)
254: 50 l=l-1
255: if (nodplc(ilr+l).ne.icol) go to 50
256: ispot=lyl+l
257: ispot2=imynl+nstop+nut+l
258: go to 70
259: c
260: c find (icol,irow) matrix term (diagonal)
261: c
262: 60 ispot=lynl+nodplc(iorder+irow)
263: ispot2=imynl+nodplc(iorder+irow)
264: c
265: 70 call cmult(value(lyl+j),value(imynl+nstop+nut+j),
266: 1 value(lyu+k),value(imynl+nstop+k),xreal,ximag)
267: value(ispot)=value(ispot)-xreal
268: value(ispot2)=value(ispot2)-ximag
269: 80 continue
270: 90 continue
271: 100 continue
272: return
273: end
274: subroutine acsol
275: implicit double precision (a-h,o-z)
276: c
277: c this routine solves the circuit equations by performing a forward
278: c and backward substitution using the previously-computed lu factors.
279: c
280: common /tabinf/ ielmnt,isbckt,nsbckt,iunsat,nunsat,itemps,numtem,
281: 1 isens,nsens,ifour,nfour,ifield,icode,idelim,icolum,insize,
282: 2 junode,lsbkpt,numbkp,iorder,jmnode,iur,iuc,ilc,ilr,numoff,isr,
283: 3 nmoffc,iseq,iseq1,neqn,nodevs,ndiag,iswap,iequa,macins,lvnim1,
284: 4 lx0,lvn,lynl,lyu,lyl,lx1,lx2,lx3,lx4,lx5,lx6,lx7,ld0,ld1,ltd,
285: 5 imynl,imvn,lcvn,loutpt,nsnod,nsmat,nsval,icnod,icmat,icval
286: common /cirdat/ locate(50),jelcnt(50),nunods,ncnods,numnod,nstop,
287: 1 nut,nlt,nxtrm,ndist,ntlin,ibr,numvs
288: common /blank/ value(1000)
289: integer nodplc(64)
290: complex*16 cvalue(32)
291: equivalence (value(1),nodplc(1),cvalue(1))
292: c
293: c forward substitution
294: c
295: do 20 i=2,nstop
296: jstart=nodplc(ilc+i)
297: jstop=nodplc(ilc+i+1)-1
298: if (jstart.gt.jstop) go to 20
299: io=nodplc(iorder+i)
300: if (value(lvn+io).ne.0.0d0) go to 5
301: if (value(imvn+io).eq.0.0d0) go to 20
302: 5 do 10 j=jstart,jstop
303: jo=nodplc(ilr+j)
304: jo=nodplc(iorder+jo)
305: call cmult(value(lyl+j),value(imynl+nstop+nut+j),value(lvn+io),
306: 1 value(imvn+io),xreal,ximag)
307: value(lvn+jo)=value(lvn+jo)-xreal
308: value(imvn+jo)=value(imvn+jo)-ximag
309: 10 continue
310: 20 continue
311: c
312: c back substitution
313: c
314: k=nstop+1
315: do 50 i=2,nstop
316: k=k-1
317: io=nodplc(iorder+k)
318: jstart=nodplc(iur+k)
319: jstop=nodplc(iur+k+1)-1
320: if (jstart.gt.jstop) go to 40
321: do 30 j=jstart,jstop
322: jo=nodplc(iuc+j)
323: jo=nodplc(iorder+jo)
324: call cmult(value(lyu+j),value(imynl+nstop+j),value(lvn+jo),
325: 1 value(imvn+jo),xreal,ximag)
326: value(lvn+io)=value(lvn+io)-xreal
327: value(imvn+io)=value(imvn+io)-ximag
328: 30 continue
329: 40 call cdiv(value(lvn+io),value(imvn+io),value(lynl+io),
330: 1 value(imynl+io),value(lvn+io),value(imvn+io))
331: 50 continue
332: do 100 i=2,nstop
333: cvalue(lcvn+i)=dcmplx(value(lvn+i),value(imvn+i))
334: 100 continue
335: cvalue(lcvn+1)=dcmplx(0.0d0,0.0d0)
336: return
337: end
338: subroutine acload
339: implicit double precision (a-h,o-z)
340: c
341: c this routine zeroes-out and then loads the complex coefficient
342: c matrix.
343: c
344: common /tabinf/ ielmnt,isbckt,nsbckt,iunsat,nunsat,itemps,numtem,
345: 1 isens,nsens,ifour,nfour,ifield,icode,idelim,icolum,insize,
346: 2 junode,lsbkpt,numbkp,iorder,jmnode,iur,iuc,ilc,ilr,numoff,isr,
347: 3 nmoffc,iseq,iseq1,neqn,nodevs,ndiag,iswap,iequa,macins,lvnim1,
348: 4 lx0,lvn,lynl,lyu,lyl,lx1,lx2,lx3,lx4,lx5,lx6,lx7,ld0,ld1,ltd,
349: 5 imynl,imvn,lcvn,loutpt,nsnod,nsmat,nsval,icnod,icmat,icval
350: common /cirdat/ locate(50),jelcnt(50),nunods,ncnods,numnod,nstop,
351: 1 nut,nlt,nxtrm,ndist,ntlin,ibr,numvs
352: common /status/ omega,time,delta,delold(7),ag(7),vt,xni,egfet,
353: 1 xmu,mode,modedc,icalc,initf,method,iord,maxord,noncon,iterno,
354: 2 itemno,nosolv,ipostp,iscrch
355: common /blank/ value(1000)
356: integer nodplc(64)
357: complex*16 cvalue(32)
358: equivalence (value(1),nodplc(1),cvalue(1))
359: c
360: c
361: complex*16 cval
362: c
363: c zero y matrix and current vector
364: c
365: call zero8(value(lvn+1),nstop)
366: call zero8(value(imvn+1),nstop)
367: call zero8(value(lynl+1),nstop+nut+nlt)
368: call zero8(value(imynl+1),nstop+nut+nlt)
369: c
370: c resistors
371: c
372: loc=locate(1)
373: 20 if (loc.eq.0) go to 30
374: locv=nodplc(loc+1)
375: val=value(locv+1)
376: locy=lynl+nodplc(loc+6)
377: value(locy)=value(locy)+val
378: locy=lynl+nodplc(loc+7)
379: value(locy)=value(locy)+val
380: locy=lynl+nodplc(loc+4)
381: value(locy)=value(locy)-val
382: locy=lynl+nodplc(loc+5)
383: value(locy)=value(locy)-val
384: loc=nodplc(loc)
385: go to 20
386: c
387: c capacitors
388: c
389: 30 loc=locate(2)
390: 40 if (loc.eq.0) go to 50
391: locv=nodplc(loc+1)
392: val=omega*value(locv+1)
393: locyi=imynl+nodplc(loc+10)
394: value(locyi)=value(locyi)+val
395: locyi=imynl+nodplc(loc+11)
396: value(locyi)=value(locyi)+val
397: locyi=imynl+nodplc(loc+5)
398: value(locyi)=value(locyi)-val
399: locyi=imynl+nodplc(loc+6)
400: value(locyi)=value(locyi)-val
401: loc=nodplc(loc)
402: go to 40
403: c
404: c inductors
405: c
406: 50 loc=locate(3)
407: 60 if (loc.eq.0) go to 70
408: locv=nodplc(loc+1)
409: val=omega*value(locv+1)
410: locyi=imynl+nodplc(loc+13)
411: locy=lynl+nodplc(loc+13)
412: value(locy)=0.0d0
413: value(locyi)=-val
414: locy=lynl+nodplc(loc+6)
415: locyi=imynl+nodplc(loc+6)
416: value(locy)=1.0d0
417: value(locyi)=0.0d0
418: locy=lynl+nodplc(loc+7)
419: locyi=imynl+nodplc(loc+7)
420: value(locy)=-1.0d0
421: value(locyi)=0.0d0
422: locy=lynl+nodplc(loc+8)
423: locyi=imynl+nodplc(loc+8)
424: value(locy)=1.0d0
425: value(locyi)=0.0d0
426: locy=lynl+nodplc(loc+9)
427: locyi=imynl+nodplc(loc+9)
428: value(locy)=-1.0d0
429: value(locyi)=0.0d0
430: loc=nodplc(loc)
431: go to 60
432: c
433: c mutual inductors
434: c
435: 70 loc=locate(4)
436: 80 if (loc.eq.0) go to 90
437: locv=nodplc(loc+1)
438: val=omega*value(locv+1)
439: locy=lynl+nodplc(loc+4)
440: locyi=imynl+nodplc(loc+4)
441: value(locy)=0.0d0
442: value(locyi)=-val
443: locy=lynl+nodplc(loc+5)
444: locyi=imynl+nodplc(loc+5)
445: value(locy)=0.0d0
446: value(locyi)=-val
447: loc=nodplc(loc)
448: go to 80
449: c
450: c nonlinear voltage controlled current sources
451: c
452: 90 loc=locate(5)
453: 95 if (loc.eq.0) go to 100
454: ndim=nodplc(loc+4)
455: lmat=nodplc(loc+7)
456: loct=lx0+nodplc(loc+12)+2
457: do 97 i=1,ndim
458: val=value(loct)
459: loct=loct+2
460: locy=lynl+nodplc(lmat+1)
461: value(locy)=value(locy)+val
462: locy=lynl+nodplc(lmat+2)
463: value(locy)=value(locy)-val
464: locy=lynl+nodplc(lmat+3)
465: value(locy)=value(locy)-val
466: locy=lynl+nodplc(lmat+4)
467: value(locy)=value(locy)+val
468: lmat=lmat+4
469: 97 continue
470: loc=nodplc(loc)
471: go to 95
472: c
473: c nonlinear voltage controlled voltage sources
474: c
475: 100 loc=locate(6)
476: 105 if (loc.eq.0) go to 110
477: ndim=nodplc(loc+4)
478: lmat=nodplc(loc+8)
479: loct=lx0+nodplc(loc+13)+3
480: locy=lynl+nodplc(lmat+1)
481: locyi=imynl+nodplc(lmat+1)
482: value(locy)=+1.0d0
483: value(locyi)=0.0d0
484: locy=lynl+nodplc(lmat+2)
485: locyi=imynl+nodplc(lmat+2)
486: value(locy)=-1.0d0
487: value(locyi)=0.0d0
488: locy=lynl+nodplc(lmat+3)
489: locyi=imynl+nodplc(lmat+3)
490: value(locy)=+1.0d0
491: value(locyi)=0.0d0
492: locy=lynl+nodplc(lmat+4)
493: locyi=imynl+nodplc(lmat+4)
494: value(locy)=-1.0d0
495: value(locyi)=0.0d0
496: lmat=lmat+4
497: do 107 i=1,ndim
498: val=value(loct)
499: loct=loct+2
500: locy=lynl+nodplc(lmat+1)
501: value(locy)=value(locy)-val
502: locy=lynl+nodplc(lmat+2)
503: value(locy)=value(locy)+val
504: lmat=lmat+2
505: 107 continue
506: loc=nodplc(loc)
507: go to 105
508: c
509: c nonlinear current controlled current sources
510: c
511: 110 loc=locate(7)
512: 115 if (loc.eq.0) go to 120
513: ndim=nodplc(loc+4)
514: lmat=nodplc(loc+7)
515: loct=lx0+nodplc(loc+12)+2
516: do 117 i=1,ndim
517: val=value(loct)
518: loct=loct+2
519: locy=lynl+nodplc(lmat+1)
520: locyi=imynl+nodplc(lmat+1)
521: value(locy)=+val
522: value(locyi)=0.0d0
523: locy=lynl+nodplc(lmat+2)
524: locyi=imynl+nodplc(lmat+2)
525: value(locy)=-val
526: value(locyi)=0.0d0
527: lmat=lmat+2
528: 117 continue
529: loc=nodplc(loc)
530: go to 115
531: c
532: c nonlinear current controlled voltage sources
533: c
534: 120 loc=locate(8)
535: 125 if (loc.eq.0) go to 140
536: ndim=nodplc(loc+4)
537: lmat=nodplc(loc+8)
538: loct=lx0+nodplc(loc+13)+3
539: locy=lynl+nodplc(lmat+1)
540: locyi=imynl+nodplc(lmat+1)
541: value(locy)=+1.0d0
542: value(locyi)=0.0d0
543: locy=lynl+nodplc(lmat+2)
544: locyi=imynl+nodplc(lmat+2)
545: value(locy)=-1.0d0
546: value(locyi)=0.0d0
547: locy=lynl+nodplc(lmat+3)
548: locyi=imynl+nodplc(lmat+3)
549: value(locy)=+1.0d0
550: value(locyi)=0.0d0
551: locy=lynl+nodplc(lmat+4)
552: locyi=imynl+nodplc(lmat+4)
553: value(locy)=-1.0d0
554: value(locyi)=0.0d0
555: lmat=lmat+4
556: do 127 i=1,ndim
557: val=value(loct)
558: loct=loct+2
559: locy=lynl+nodplc(lmat+i)
560: value(locy)=value(locy)-val
561: 127 continue
562: loc=nodplc(loc)
563: go to 125
564: c
565: c voltage sources
566: c
567: 140 loc=locate(9)
568: 150 if (loc.eq.0) go to 160
569: locv=nodplc(loc+1)
570: iptr=nodplc(loc+6)
571: value(lvn+iptr)=value(locv+2)
572: value(imvn+iptr)=value(locv+3)
573: locy=lynl+nodplc(loc+7)
574: value(locy)=value(locy)+1.0d0
575: locy=lynl+nodplc(loc+8)
576: value(locy)=value(locy)-1.0d0
577: locy=lynl+nodplc(loc+9)
578: value(locy)=value(locy)+1.0d0
579: locy=lynl+nodplc(loc+10)
580: value(locy)=value(locy)-1.0d0
581: loc=nodplc(loc)
582: go to 150
583: c
584: c current sources
585: c
586: 160 loc=locate(10)
587: 170 if (loc.eq.0) go to 200
588: locv=nodplc(loc+1)
589: node1=nodplc(loc+2)
590: node2=nodplc(loc+3)
591: value(lvn+node1)=value(lvn+node1)-value(locv+2)
592: value(imvn+node1)=value(imvn+node1)-value(locv+3)
593: value(lvn+node2)=value(lvn+node2)+value(locv+2)
594: value(imvn+node2)=value(imvn+node2)+value(locv+3)
595: loc=nodplc(loc)
596: go to 170
597: c
598: c diodes
599: c
600: 200 loc=locate(11)
601: 210 if (loc.eq.0) go to 250
602: locv=nodplc(loc+1)
603: area=value(locv+1)
604: locm=nodplc(loc+5)
605: locm=nodplc(locm+1)
606: loct=lx0+nodplc(loc+11)
607: gspr=value(locm+2)*area
608: geq=value(loct+2)
609: xceq=value(loct+4)*omega
610: locy=lynl+nodplc(loc+13)
611: value(locy)=value(locy)+gspr
612: locy=lynl+nodplc(loc+14)
613: locyi=imynl+nodplc(loc+14)
614: value(locy)=value(locy)+geq
615: value(locyi)=value(locyi)+xceq
616: locy=lynl+nodplc(loc+15)
617: locyi=imynl+nodplc(loc+15)
618: value(locy)=value(locy)+geq+gspr
619: value(locyi)=value(locyi)+xceq
620: locy=lynl+nodplc(loc+7)
621: value(locy)=value(locy)-gspr
622: locy=lynl+nodplc(loc+8)
623: locyi=imynl+nodplc(loc+8)
624: value(locy)=value(locy)-geq
625: value(locyi)=value(locyi)-xceq
626: locy=lynl+nodplc(loc+9)
627: value(locy)=value(locy)-gspr
628: locy=lynl+nodplc(loc+10)
629: locyi=imynl+nodplc(loc+10)
630: value(locy)=value(locy)-geq
631: value(locyi)=value(locyi)-xceq
632: loc=nodplc(loc)
633: go to 210
634: c
635: c bjts
636: c
637: 250 loc=locate(12)
638: 260 if (loc.eq.0) go to 300
639: locv=nodplc(loc+1)
640: area=value(locv+1)
641: locm=nodplc(loc+8)
642: locm=nodplc(locm+1)
643: loct=lx0+nodplc(loc+22)
644: gcpr=value(locm+20)*area
645: gepr=value(locm+19)*area
646: gpi=value(loct+4)
647: gmu=value(loct+5)
648: gm=value(loct+6)
649: go=value(loct+7)
650: xgm=0.0d0
651: td=value(locm+28)
652: if(td.eq.0.0d0) go to 270
653: arg=td*omega
654: gm=gm+go
655: xgm=-gm*dsin(arg)
656: gm=gm*dcos(arg)-go
657: 270 gx=value(loct+16)
658: xcpi=value(loct+9)*omega
659: xcmu=value(loct+11)*omega
660: xcbx=value(loct+15)*omega
661: xccs=value(loct+13)*omega
662: xcmcb=value(loct+17)*omega
663: locy=lynl+nodplc(loc+24)
664: value(locy)=value(locy)+gcpr
665: locy=lynl+nodplc(loc+25)
666: locyi=imynl+nodplc(loc+25)
667: value(locy)=value(locy)+gx
668: value(locyi)=value(locyi)+xcbx
669: locy=lynl+nodplc(loc+26)
670: value(locy)=value(locy)+gepr
671: locy=lynl+nodplc(loc+27)
672: locyi=imynl+nodplc(loc+27)
673: value(locy)=value(locy)+gmu+go+gcpr
674: value(locyi)=value(locyi)+xcmu+xccs+xcbx
675: locy=lynl+nodplc(loc+28)
676: locyi=imynl+nodplc(loc+28)
677: value(locy)=value(locy)+gx+gpi+gmu
678: value(locyi)=value(locyi)+xcpi+xcmu+xcmcb
679: locy=lynl+nodplc(loc+29)
680: locyi=imynl+nodplc(loc+29)
681: value(locy)=value(locy)+gpi+gepr+gm+go
682: value(locyi)=value(locyi)+xcpi+xgm
683: locy=lynl+nodplc(loc+10)
684: value(locy)=value(locy)-gcpr
685: locy=lynl+nodplc(loc+11)
686: value(locy)=value(locy)-gx
687: locy=lynl+nodplc(loc+12)
688: value(locy)=value(locy)-gepr
689: locy=lynl+nodplc(loc+13)
690: value(locy)=value(locy)-gcpr
691: locy=lynl+nodplc(loc+14)
692: locyi=imynl+nodplc(loc+14)
693: value(locy)=value(locy)-gmu+gm
694: value(locyi)=value(locyi)-xcmu+xgm
695: locy=lynl+nodplc(loc+15)
696: locyi=imynl+nodplc(loc+15)
697: value(locy)=value(locy)-gm-go
698: value(locyi)=value(locyi)-xgm
699: locy=lynl+nodplc(loc+16)
700: value(locy)=value(locy)-gx
701: locy=lynl+nodplc(loc+17)
702: locyi=imynl+nodplc(loc+17)
703: value(locy)=value(locy)-gmu
704: value(locyi)=value(locyi)-xcmu-xcmcb
705: locy=lynl+nodplc(loc+18)
706: locyi=imynl+nodplc(loc+18)
707: value(locy)=value(locy)-gpi
708: value(locyi)=value(locyi)-xcpi
709: locy=lynl+nodplc(loc+19)
710: value(locy)=value(locy)-gepr
711: locy=lynl+nodplc(loc+20)
712: locyi=imynl+nodplc(loc+20)
713: value(locy)=value(locy)-go
714: value(locyi)=value(locyi)+xcmcb
715: locy=lynl+nodplc(loc+21)
716: locyi=imynl+nodplc(loc+21)
717: value(locy)=value(locy)-gpi-gm
718: value(locyi)=value(locyi)-xcpi-xgm-xcmcb
719: locyi=imynl+nodplc(loc+31)
720: value(locyi)=value(locyi)+xccs
721: locyi=imynl+nodplc(loc+32)
722: value(locyi)=value(locyi)-xccs
723: locyi=imynl+nodplc(loc+33)
724: value(locyi)=value(locyi)-xccs
725: locyi=imynl+nodplc(loc+34)
726: value(locyi)=value(locyi)-xcbx
727: locyi=imynl+nodplc(loc+35)
728: value(locyi)=value(locyi)-xcbx
729: loc=nodplc(loc)
730: go to 260
731: c
732: c jfets
733: c
734: 300 loc=locate(13)
735: 310 if (loc.eq.0) go to 350
736: locv=nodplc(loc+1)
737: area=value(locv+1)
738: locm=nodplc(loc+7)
739: locm=nodplc(locm+1)
740: loct=lx0+nodplc(loc+19)
741: gdpr=value(locm+4)*area
742: gspr=value(locm+5)*area
743: gm=value(loct+5)
744: gds=value(loct+6)
745: ggs=value(loct+7)
746: xgs=value(loct+9)*omega
747: ggd=value(loct+8)
748: xgd=value(loct+11)*omega
749: locy=lynl+nodplc(loc+20)
750: value(locy)=value(locy)+gdpr
751: locy=lynl+nodplc(loc+21)
752: locyi=imynl+nodplc(loc+21)
753: value(locy)=value(locy)+ggd+ggs
754: value(locyi)=value(locyi)+xgd+xgs
755: locy=lynl+nodplc(loc+22)
756: value(locy)=value(locy)+gspr
757: locy=lynl+nodplc(loc+23)
758: locyi=imynl+nodplc(loc+23)
759: value(locy)=value(locy)+gdpr+gds+ggd
760: value(locyi)=value(locyi)+xgd
761: locy=lynl+nodplc(loc+24)
762: locyi=imynl+nodplc(loc+24)
763: value(locy)=value(locy)+gspr+gds+gm+ggs
764: value(locyi)=value(locyi)+xgs
765: locy=lynl+nodplc(loc+9)
766: value(locy)=value(locy)-gdpr
767: locy=lynl+nodplc(loc+10)
768: locyi=imynl+nodplc(loc+10)
769: value(locy)=value(locy)-ggd
770: value(locyi)=value(locyi)-xgd
771: locy=lynl+nodplc(loc+11)
772: locyi=imynl+nodplc(loc+11)
773: value(locy)=value(locy)-ggs
774: value(locyi)=value(locyi)-xgs
775: locy=lynl+nodplc(loc+12)
776: value(locy)=value(locy)-gspr
777: locy=lynl+nodplc(loc+13)
778: value(locy)=value(locy)-gdpr
779: locy=lynl+nodplc(loc+14)
780: locyi=imynl+nodplc(loc+14)
781: value(locy)=value(locy)-ggd+gm
782: value(locyi)=value(locyi)-xgd
783: locy=lynl+nodplc(loc+15)
784: value(locy)=value(locy)-gds-gm
785: locy=lynl+nodplc(loc+16)
786: locyi=imynl+nodplc(loc+16)
787: value(locy)=value(locy)-ggs-gm
788: value(locyi)=value(locyi)-xgs
789: locy=lynl+nodplc(loc+17)
790: value(locy)=value(locy)-gspr
791: locy=lynl+nodplc(loc+18)
792: value(locy)=value(locy)-gds
793: loc=nodplc(loc)
794: go to 310
795: c
796: c mosfets
797: c
798: 350 loc=locate(14)
799: 360 if (loc.eq.0) go to 400
800: locv=nodplc(loc+1)
801: locm=nodplc(loc+8)
802: itype=nodplc(locm+2)
803: locm=nodplc(locm+1)
804: loct=lx0+nodplc(loc+26)
805: gdpr=value(locv+11)
806: gspr=value(locv+12)
807: if(itype.eq.0) go to 380
808: xl=value(locv+1)-2.0d0*value(locm+20)
809: xw=value(locv+2)-2.0d0*value(locm+36)
810: covlgs=value(locm+8)*xw
811: covlgd=value(locm+9)*xw
812: covlgb=value(locm+10)*xl
813: didvg=value(loct)
814: didvd=value(loct+1)
815: didvs=value(loct+2)
816: geqbd=value(loct+3)
817: geqbs=value(loct+5)
818: ccgg=value(loct+8)
819: ccgd=value(loct+9)
820: ccgs=value(loct+10)
821: ccbg=value(loct+11)
822: ccbd=value(loct+12)
823: ccbs=value(loct+13)
824: capbd=value(loct+14)
825: capbs=value(loct+15)
826: xccg2=-0.5d0*(ccgg+ccbg)*omega
827: xccd2=-0.5d0*(ccgd+ccbd)*omega
828: xccs2=-0.5d0*(ccgs+ccbs)*omega
829: xccdg=xccg2-covlgd*omega
830: xccdd=xccd2+(capbd+covlgd)*omega
831: xccds=xccs2
832: xccsg=xccg2-covlgs*omega
833: xccsd=xccd2
834: xccss=xccs2+(capbs+covlgs)*omega
835: xccgg=(ccgg+covlgd+covlgs+covlgb)*omega
836: xccgd=(ccgd-covlgd)*omega
837: xccgs=(ccgs-covlgs)*omega
838: xccbg=(ccbg-covlgb)*omega
839: xccbd=xccbd+(ccbd-capbd)*omega
840: xccbs=xccbs+(ccbs-capbs)*omega
841: gccdg=didvg
842: gccdd=geqbd+didvd
843: gccds=didvs
844: gccsg=-didvg
845: gccsd=-didvd
846: gccss=geqbs-didvs
847: gccbd=-geqbd
848: gccbs=-geqbs
849: locyi=imynl+nodplc(loc+28)
850: value(locyi)=value(locyi)+xccgg
851: locyi=imynl+nodplc(loc+30)
852: value(locyi)=value(locyi)-xccbg-xccbd-xccbs
853: locyi=imynl+nodplc(loc+31)
854: value(locyi)=value(locyi)+xccdd
855: locyi=imynl+nodplc(loc+32)
856: value(locyi)=value(locyi)+xccss
857: locyi=imynl+nodplc(loc+11)
858: value(locyi)=value(locyi)-xccgg-xccgd-xccgs
859: locyi=imynl+nodplc(loc+12)
860: value(locyi)=value(locyi)+xccgd
861: locyi=imynl+nodplc(loc+13)
862: value(locyi)=value(locyi)+xccgs
863: locyi=imynl+nodplc(loc+15)
864: value(locyi)=value(locyi)+xccbg
865: locyi=imynl+nodplc(loc+16)
866: value(locyi)=value(locyi)+xccbd
867: locyi=imynl+nodplc(loc+17)
868: value(locyi)=value(locyi)+xccbs
869: locyi=imynl+nodplc(loc+19)
870: value(locyi)=value(locyi)+xccdg
871: locyi=imynl+nodplc(loc+20)
872: value(locyi)=value(locyi)-xccdg-xccdd-xccds
873: locyi=imynl+nodplc(loc+21)
874: value(locyi)=value(locyi)+xccds
875: locyi=imynl+nodplc(loc+22)
876: value(locyi)=value(locyi)+xccsg
877: locyi=imynl+nodplc(loc+24)
878: value(locyi)=value(locyi)-xccsg-xccsd-xccss
879: locyi=imynl+nodplc(loc+25)
880: value(locyi)=value(locyi)+xccsd
881: locy=lynl+nodplc(loc+27)
882: value(locy)=value(locy)+gdpr
883: locy=lynl+nodplc(loc+29)
884: value(locy)=value(locy)+gspr
885: locy=lynl+nodplc(loc+30)
886: value(locy)=value(locy)-gccbd-gccbs
887: locy=lynl+nodplc(loc+31)
888: value(locy)=value(locy)+gdpr+gccdd
889: locy=lynl+nodplc(loc+32)
890: value(locy)=value(locy)+gspr+gccss
891: locy=lynl+nodplc(loc+10)
892: value(locy)=value(locy)-gdpr
893: locy=lynl+nodplc(loc+14)
894: value(locy)=value(locy)-gspr
895: locy=lynl+nodplc(loc+16)
896: value(locy)=value(locy)+gccbd
897: locy=lynl+nodplc(loc+17)
898: value(locy)=value(locy)+gccbs
899: locy=lynl+nodplc(loc+18)
900: value(locy)=value(locy)-gdpr
901: locy=lynl+nodplc(loc+19)
902: value(locy)=value(locy)+gccdg
903: locy=lynl+nodplc(loc+20)
904: value(locy)=value(locy)-gccdg-gccdd-gccds
905: locy=lynl+nodplc(loc+21)
906: value(locy)=value(locy)+gccds
907: locy=lynl+nodplc(loc+22)
908: value(locy)=value(locy)+gccsg
909: locy=lynl+nodplc(loc+23)
910: value(locy)=value(locy)-gspr
911: locy=lynl+nodplc(loc+24)
912: value(locy)=value(locy)-gccsg-gccsd-gccss
913: locy=lynl+nodplc(loc+25)
914: value(locy)=value(locy)+gccsd
915: loc=nodplc(loc)
916: go to 360
917: c... ga-as fets
918: 380 continue
919: devmod=value(locv+8)
920: ggd=value(loct+6)
921: gm=0.0d0
922: gds=0.0d0
923: gmj1=value(loct+7)
924: gdb=value(loct+8)
925: ggs=value(loct+9)
926: xcds=value(loct+10)*omega
927: gsb=value(loct+11)
928: xcgs=value(loct+12)*omega
929: gmj2=value(loct+13)
930: xcgd=value(loct+14)*omega
931: xcgb=value(loct+16)*omega
932: c write(6,1001) ggd,gm,gds,gmj1,gdb,ggs,gsb,gmj2
933: 1001 format(' ggd gm gds gmj1 gdb ggs gsb gmj2'/,1x,1p8e9.1)
934: if(devmod.gt.0.0d0) go to 385
935: gmrev=gm
936: gmnrm=0.0d0
937: gm=0.0d0
938: gmj1r=gmj1
939: gmj1=0.0d0
940: gmj1n=0.0d0
941: gmj2n=0.0d0
942: gmj2r=0.0d0
943: go to 390
944: 385 gmnrm=gm
945: gmrev=0.0d0
946: gm=0.0d0
947: gmj2n=gmj2
948: gmj2r=0.0d0
949: gmj2=0.0d0
950: gmj1r=0.0d0
951: gmj1n=0.0d0
952: 390 locy=lynl+nodplc(loc+27)
953: value(locy)=value(locy)+gdpr
954: locy=lynl+nodplc(loc+28)
955: locyi=imynl+nodplc(loc+28)
956: value(locy)=value(locy)+ggd+ggs
957: value(locyi)=value(locyi)+xcgd+xcgs+xcgb
958: locy=lynl+nodplc(loc+29)
959: value(locy)=value(locy)+gspr
960: locy=lynl+nodplc(loc+30)
961: locyi=imynl+nodplc(loc+30)
962: value(locy)=value(locy)+gdb+gsb+gmj1+gmj2
963: value(locyi)=value(locyi)+xcgb
964: locy=lynl+nodplc(loc+31)
965: locyi=imynl+nodplc(loc+31)
966: value(locy)=value(locy)+gdpr+gds+gdb+ggd-gmrev-gmj1r
967: value(locyi)=value(locyi)+xcds+xcgd
968: locy=lynl+nodplc(loc+32)
969: locyi=imynl+nodplc(loc+32)
970: value(locy)=value(locy)+gspr+gds+gsb+ggs+gmnrm-gmj2n
971: value(locyi)=value(locyi)+xcds+xcgs
972: locy=lynl+nodplc(loc+10)
973: value(locy)=value(locy)-gdpr
974: locyi=imynl+nodplc(loc+11)
975: value(locyi)=value(locyi)-xcgb
976: locy=lynl+nodplc(loc+12)
977: locyi=imynl+nodplc(loc+12)
978: value(locy)=value(locy)-ggd
979: value(locyi)=value(locyi)-xcgd
980: locy=lynl+nodplc(loc+13)
981: locyi=imynl+nodplc(loc+13)
982: value(locy)=value(locy)-ggs
983: value(locyi)=value(locyi)-xcgs
984: locy=lynl+nodplc(loc+14)
985: value(locy)=value(locy)-gspr
986: locy=lynl+nodplc(loc+15)
987: locyi=imynl+nodplc(loc+15)
988: value(locy)=value(locy)-gmj2n-gmj1n-gmj1r-gmj2r-gmj1-gmj2
989: value(locyi)=value(locyi)-xcgb
990: locy=lynl+nodplc(loc+16)
991: value(locy)=value(locy)-gdb+gmj2r+gmj1r
992: locy=lynl+nodplc(loc+17)
993: value(locy)=value(locy)-gsb+gmj1n+gmj2n
994: locy=lynl+nodplc(loc+18)
995: value(locy)=value(locy)-gdpr
996: locy=lynl+nodplc(loc+19)
997: locyi=imynl+nodplc(loc+19)
998: value(locy)=value(locy)-ggd+gmnrm+gmrev+gmj1r+gmj1n+gmj1+gm
999: value(locyi)=value(locyi)-xcgd
1000: locy=lynl+nodplc(loc+20)
1001: value(locy)=value(locy)-gdb-gmj1-gm
1002: locy=lynl+nodplc(loc+21)
1003: locyi=imynl+nodplc(loc+21)
1004: value(locy)=value(locy)-gds-gmnrm-gmj1n
1005: value(locyi)=value(locyi)-xcds
1006: locy=lynl+nodplc(loc+22)
1007: locyi=imynl+nodplc(loc+22)
1008: value(locy)=value(locy)-ggs-gmnrm-gmrev+gmj2r+gmj2n+gmj2-gm
1009: value(locyi)=value(locyi)-xcgs
1010: locy=lynl+nodplc(loc+23)
1011: value(locy)=value(locy)-gspr
1012: locy=lynl+nodplc(loc+24)
1013: value(locy)=value(locy)-gsb-gmj2+gm
1014: locy=lynl+nodplc(loc+25)
1015: locyi=imynl+nodplc(loc+25)
1016: value(locy)=value(locy)-gds+gmrev-gmj2r
1017: value(locyi)=value(locyi)-xcds
1018: loc=nodplc(loc)
1019: go to 360
1020: c
1021: c transmission lines
1022: c
1023: 400 loc=locate(17)
1024: 410 if (loc.eq.0) go to 1000
1025: locv=nodplc(loc+1)
1026: z0=value(locv+1)
1027: y0=1.0d0/z0
1028: td=value(locv+2)
1029: arg=-omega*td
1030: rval=dcos(arg)
1031: xval=dsin(arg)
1032: locy=lynl+nodplc(loc+10)
1033: value(locy)=value(locy)+y0
1034: locy=lynl+nodplc(loc+11)
1035: locyi=imynl+nodplc(loc+11)
1036: value(locy)=-y0
1037: value(locyi)=0.0d0
1038: locy=lynl+nodplc(loc+12)
1039: locyi=imynl+nodplc(loc+12)
1040: value(locy)=-1.0d0
1041: value(locyi)=0.0d0
1042: locy=lynl+nodplc(loc+13)
1043: value(locy)=value(locy)+y0
1044: locy=lynl+nodplc(loc+14)
1045: locyi=imynl+nodplc(loc+14)
1046: value(locy)=-1.0d0
1047: value(locyi)=0.0d0
1048: locy=lynl+nodplc(loc+15)
1049: locyi=imynl+nodplc(loc+15)
1050: value(locy)=-y0
1051: value(locyi)=0.0d0
1052: locy=lynl+nodplc(loc+16)
1053: locyi=imynl+nodplc(loc+16)
1054: value(locy)=+y0
1055: value(locyi)=0.0d0
1056: locy=lynl+nodplc(loc+17)
1057: locyi=imynl+nodplc(loc+17)
1058: value(locy)=+1.0d0
1059: value(locyi)=0.0d0
1060: locy=lynl+nodplc(loc+18)
1061: locyi=imynl+nodplc(loc+18)
1062: value(locy)=+y0
1063: value(locyi)=0.0d0
1064: locy=lynl+nodplc(loc+19)
1065: locyi=imynl+nodplc(loc+19)
1066: value(locy)=+1.0d0
1067: value(locyi)=0.0d0
1068: locy=lynl+nodplc(loc+20)
1069: locyi=imynl+nodplc(loc+20)
1070: value(locy)=-1.0d0
1071: value(locyi)=0.0d0
1072: locy=lynl+nodplc(loc+21)
1073: locyi=imynl+nodplc(loc+21)
1074: value(locy)=-rval
1075: value(locyi)=-xval
1076: locy=lynl+nodplc(loc+22)
1077: locyi=imynl+nodplc(loc+22)
1078: value(locy)=+rval
1079: value(locyi)=+xval
1080: locy=lynl+nodplc(loc+23)
1081: locyi=imynl+nodplc(loc+23)
1082: value(locy)=+1.0d0
1083: value(locyi)=0.0d0
1084: locy=lynl+nodplc(loc+24)
1085: locyi=imynl+nodplc(loc+24)
1086: value(locy)=-rval*z0
1087: value(locyi)=-xval*z0
1088: locy=lynl+nodplc(loc+25)
1089: locyi=imynl+nodplc(loc+25)
1090: value(locy)=-rval
1091: value(locyi)=-xval
1092: locy=lynl+nodplc(loc+26)
1093: locyi=imynl+nodplc(loc+26)
1094: value(locy)=+rval
1095: value(locyi)=+xval
1096: locy=lynl+nodplc(loc+27)
1097: locyi=imynl+nodplc(loc+27)
1098: value(locy)=-1.0d0
1099: value(locyi)=0.0d0
1100: locy=lynl+nodplc(loc+28)
1101: locyi=imynl+nodplc(loc+28)
1102: value(locy)=+1.0d0
1103: value(locyi)=0.0d0
1104: locy=lynl+nodplc(loc+29)
1105: locyi=imynl+nodplc(loc+29)
1106: value(locy)=-rval*z0
1107: value(locyi)=-xval*z0
1108: locy=lynl+nodplc(loc+31)
1109: locyi=imynl+nodplc(loc+31)
1110: value(locy)=-y0
1111: value(locyi)=0.0d0
1112: locy=lynl+nodplc(loc+32)
1113: locyi=imynl+nodplc(loc+32)
1114: value(locy)=-y0
1115: value(locyi)=0.0d0
1116: loc=nodplc(loc)
1117: go to 410
1118: c
1119: c reorder right-hand side
1120: c
1121: 1000 do 1010 i=2,nstop
1122: j=nodplc(iswap+i)
1123: value(ndiag+i)=value(lvn+j)
1124: value(ndiag+i+nstop)=value(imvn+j)
1125: 1010 continue
1126: call copy8(value(ndiag+1),value(lvn+1),nstop)
1127: call copy8(value(ndiag+nstop+1),value(imvn+1),nstop)
1128: c
1129: c finished
1130: c
1131: return
1132: end
1133: subroutine noise(loco)
1134: implicit double precision (a-h,o-z)
1135: c
1136: c this routine computes the noise due to various circuit elements.
1137: c
1138: common /tabinf/ ielmnt,isbckt,nsbckt,iunsat,nunsat,itemps,numtem,
1139: 1 isens,nsens,ifour,nfour,ifield,icode,idelim,icolum,insize,
1140: 2 junode,lsbkpt,numbkp,iorder,jmnode,iur,iuc,ilc,ilr,numoff,isr,
1141: 3 nmoffc,iseq,iseq1,neqn,nodevs,ndiag,iswap,iequa,macins,lvnim1,
1142: 4 lx0,lvn,lynl,lyu,lyl,lx1,lx2,lx3,lx4,lx5,lx6,lx7,ld0,ld1,ltd,
1143: 5 imynl,imvn,lcvn,loutpt,nsnod,nsmat,nsval,icnod,icmat,icval
1144: common /cirdat/ locate(50),jelcnt(50),nunods,ncnods,numnod,nstop,
1145: 1 nut,nlt,nxtrm,ndist,ntlin,ibr,numvs
1146: common /status/ omega,time,delta,delold(7),ag(7),vt,xni,egfet,
1147: 1 xmu,mode,modedc,icalc,initf,method,iord,maxord,noncon,iterno,
1148: 2 itemno,nosolv,ipostp,iscrch
1149: common /miscel/ atime,aprog(3),adate,atitle(10),defl,defw,defad,
1150: 1 defas,rstats(50),iwidth,lwidth,nopage
1151: common /knstnt/ twopi,xlog2,xlog10,root2,rad,boltz,charge,ctok,
1152: 1 gmin,reltol,abstol,vntol,trtol,chgtol,eps0,epssil,epsox
1153: common /ac/ fstart,fstop,fincr,skw2,refprl,spw2,jacflg,idfreq,
1154: 1 inoise,nosprt,nosout,nosin,idist,idprt
1155: common /blank/ value(1000)
1156: integer nodplc(64)
1157: complex*16 cvalue(32)
1158: equivalence (value(1),nodplc(1),cvalue(1))
1159: c
1160: c
1161: dimension vno1(12),vno2(12),vno3(12),vno4(12),vno5(12),vno6(12)
1162: dimension vntot(12),anam(12),string(5)
1163: dimension titln(4),v(2)
1164: dimension afmt1(3),afmt2(3)
1165: complex*16 cval,c(1)
1166: equivalence (c(1),v(1),cval)
1167: equivalence (v(1),vreal),(v(2),vimag)
1168: data titln / 8hnoise an, 8halysis , 8h , 8h /
1169: data alsrb,alsrc,alsre,alsrs,alsrd / 2hrb,2hrc,2hre,2hrs,2hrd /
1170: data alsib,alsic,alsid,alsfn / 2hib,2hic,2hid,2hfn /
1171: data alstot / 5htotal /
1172: data aslash,ablnk / 1h/, 1h /
1173: data afmt1 /8h(////,11,8hx, (2x,,8ha8)) /
1174: data afmt2 /8h(1h0,a8,,8h1p d10.,8h3) /
1175: kntr=12
1176: if(lwidth.le.80) kntr=7
1177: ipos=11
1178: call move(afmt1,ipos,ablnk,1,2)
1179: call alfnum(kntr,afmt1,ipos)
1180: ipos=11
1181: call move(afmt2,ipos,ablnk,1,2)
1182: call alfnum(kntr,afmt2,ipos)
1183: nprnt=0
1184: freq=omega/twopi
1185: if (icalc.ge.2) go to 10
1186: fourkt=4.0d0*charge*vt
1187: twoq=2.0d0*charge
1188: noposo=nodplc(nosout+2)
1189: nonego=nodplc(nosout+3)
1190: kntlim=lwidth/11
1191: nkntr=1
1192: 10 if (nosprt.eq.0) go to 30
1193: if (nkntr.gt.icalc) go to 30
1194: nprnt=1
1195: nkntr=nkntr+nosprt
1196: call title(0,lwidth,1,titln)
1197: write (6,16) freq
1198: 16 format('0 frequency = ',1pd10.3,' hz'/)
1199: c
1200: c obtain adjoint circuit solution
1201: c
1202: 30 vnrms=0.0d0
1203: cval=cvalue(lcvn+noposo)-cvalue(lcvn+nonego)
1204: vout=dsqrt(vreal*vreal+vimag*vimag)
1205: vout=dmax1(vout,1.0d-20)
1206: call zero8(value(lvn+1),nstop)
1207: call zero8(value(imvn+1),nstop)
1208: value(lvn+noposo)=-1.0d0
1209: value(lvn+nonego)=+1.0d0
1210: call acasol
1211: c
1212: c resistors
1213: c
1214: if (jelcnt(1).eq.0) go to 200
1215: ititle=0
1216: 91 format(//'0**** resistor squared noise voltages (sq v/hz)')
1217: 100 loc=locate(1)
1218: kntr=0
1219: 110 if (loc.eq.0) go to 130
1220: kntr=kntr+1
1221: locv=nodplc(loc+1)
1222: anam(kntr)=value(locv)
1223: node1=nodplc(loc+2)
1224: node2=nodplc(loc+3)
1225: cval=cvalue(lcvn+node1)-cvalue(lcvn+node2)
1226: vntot(kntr)=(vreal*vreal+vimag*vimag)*fourkt*value(locv+1)
1227: vnrms=vnrms+vntot(kntr)
1228: if (kntr.ge.kntlim) go to 140
1229: 120 loc=nodplc(loc)
1230: go to 110
1231: 130 if (kntr.eq.0) go to 200
1232: 140 if (nprnt.eq.0) go to 160
1233: if (ititle.eq.0) write (6,91)
1234: ititle=1
1235: write (6,afmt1) (anam(i),i=1,kntr)
1236: write (6,afmt2) alstot,(vntot(i),i=1,kntr)
1237: 160 kntr=0
1238: if (loc.ne.0) go to 120
1239: c
1240: c diodes
1241: c
1242: 200 if (jelcnt(11).eq.0) go to 300
1243: ititle=0
1244: 201 format(//'0**** diode squared noise voltages (sq v/hz)')
1245: 210 loc=locate(11)
1246: kntr=0
1247: 220 if (loc.eq.0) go to 240
1248: kntr=kntr+1
1249: locv=nodplc(loc+1)
1250: anam(kntr)=value(locv)
1251: node1=nodplc(loc+2)
1252: node2=nodplc(loc+3)
1253: node3=nodplc(loc+4)
1254: locm=nodplc(loc+5)
1255: locm=nodplc(locm+1)
1256: loct=nodplc(loc+11)
1257: area=value(locv+1)
1258: fnk=value(locm+10)
1259: fna=value(locm+11)
1260: c
1261: c ohmic resistance
1262: c
1263: cval=cvalue(lcvn+node1)-cvalue(lcvn+node3)
1264: vno1(kntr)=(vreal*vreal+vimag*vimag)*fourkt*value(locm+2)*area
1265: c
1266: c junction shot noise and flicker noise
1267: c
1268: cval=cvalue(lcvn+node3)-cvalue(lcvn+node2)
1269: vtemp=vreal*vreal+vimag*vimag
1270: arg=dmax1(dabs(value(lx0+loct+1)),1.0d-20)
1271: vno2(kntr)=vtemp*twoq*arg
1272: vno3(kntr)=vtemp*fnk*dexp(fna*dlog(arg))/freq
1273: vntot(kntr)=vno1(kntr)+vno2(kntr)+vno3(kntr)
1274: vnrms=vnrms+vntot(kntr)
1275: if (kntr.ge.kntlim) go to 250
1276: 230 loc=nodplc(loc)
1277: go to 220
1278: 240 if (kntr.eq.0) go to 300
1279: 250 if (nprnt.eq.0) go to 260
1280: if (ititle.eq.0) write (6,201)
1281: ititle=1
1282: write (6,afmt1) (anam(i),i=1,kntr)
1283: write (6,afmt2) alsrs,(vno1(i),i=1,kntr)
1284: write (6,afmt2) alsid,(vno2(i),i=1,kntr)
1285: write (6,afmt2) alsfn,(vno3(i),i=1,kntr)
1286: write (6,afmt2) alstot,(vntot(i),i=1,kntr)
1287: 260 kntr=0
1288: if (loc.ne.0) go to 230
1289: c
1290: c bipolar junction transistors
1291: c
1292: 300 if (jelcnt(12).eq.0) go to 400
1293: ititle=0
1294: 301 format(//'0**** transistor squared noise voltages (sq v/hz)')
1295: 310 loc=locate(12)
1296: kntr=0
1297: 320 if (loc.eq.0) go to 340
1298: kntr=kntr+1
1299: locv=nodplc(loc+1)
1300: anam(kntr)=value(locv)
1301: node1=nodplc(loc+2)
1302: node2=nodplc(loc+3)
1303: node3=nodplc(loc+4)
1304: node4=nodplc(loc+5)
1305: node5=nodplc(loc+6)
1306: node6=nodplc(loc+7)
1307: locm=nodplc(loc+8)
1308: locm=nodplc(locm+1)
1309: loct=nodplc(loc+22)
1310: area=value(locv+1)
1311: fnk=value(locm+44)
1312: fna=value(locm+45)
1313: c
1314: c extrinsic resistances
1315: c
1316: c... base resistance
1317: cval=cvalue(lcvn+node2)-cvalue(lcvn+node5)
1318: vno1(kntr)=(vreal*vreal+vimag*vimag)*fourkt*value(lx0+loct+16)
1319: 1 *area
1320: c... collector resistance
1321: cval=cvalue(lcvn+node1)-cvalue(lcvn+node4)
1322: vno2(kntr)=(vreal*vreal+vimag*vimag)*fourkt*value(locm+20)*area
1323: c... emitter resistance
1324: cval=cvalue(lcvn+node3)-cvalue(lcvn+node6)
1325: vno3(kntr)=(vreal*vreal+vimag*vimag)*fourkt*value(locm+19)*area
1326: c
1327: c base current shot noise and flicker noise
1328: c
1329: cval=cvalue(lcvn+node5)-cvalue(lcvn+node6)
1330: vtemp=vreal*vreal+vimag*vimag
1331: arg=dmax1(dabs(value(lx0+loct+3)),1.0d-20)
1332: vno4(kntr)=vtemp*twoq*arg
1333: vno5(kntr)=vtemp*fnk*dexp(fna*dlog(arg))/freq
1334: c
1335: c collector current shot noise
1336: c
1337: cval=cvalue(lcvn+node4)-cvalue(lcvn+node6)
1338: vno6(kntr)=(vreal*vreal+vimag*vimag)*twoq*dabs(value(lx0+loct+2))
1339: vntot(kntr)=vno1(kntr)+vno2(kntr)+vno3(kntr)+vno4(kntr)+vno5(kntr)
1340: 1 +vno6(kntr)
1341: vnrms=vnrms+vntot(kntr)
1342: if (kntr.ge.kntlim) go to 350
1343: 330 loc=nodplc(loc)
1344: go to 320
1345: 340 if (kntr.eq.0) go to 400
1346: 350 if (nprnt.eq.0) go to 360
1347: if (ititle.eq.0) write (6,301)
1348: ititle=1
1349: write (6,afmt1) (anam(i),i=1,kntr)
1350: write (6,afmt2) alsrb,(vno1(i),i=1,kntr)
1351: write (6,afmt2) alsrc,(vno2(i),i=1,kntr)
1352: write (6,afmt2) alsre,(vno3(i),i=1,kntr)
1353: write (6,afmt2) alsib,(vno4(i),i=1,kntr)
1354: write (6,afmt2) alsic,(vno6(i),i=1,kntr)
1355: write (6,afmt2) alsfn,(vno5(i),i=1,kntr)
1356: write (6,afmt2) alstot,(vntot(i),i=1,kntr)
1357: 360 kntr=0
1358: if (loc.ne.0) go to 330
1359: c
1360: c jfets
1361: c
1362: 400 if (jelcnt(13).eq.0) go to 500
1363: ititle=0
1364: 401 format(//'0**** jfet squared noise voltages (sq v/hz)')
1365: 410 loc=locate(13)
1366: kntr=0
1367: 420 if (loc.eq.0) go to 440
1368: kntr=kntr+1
1369: locv=nodplc(loc+1)
1370: anam(kntr)=value(locv)
1371: node1=nodplc(loc+2)
1372: node2=nodplc(loc+3)
1373: node3=nodplc(loc+4)
1374: node4=nodplc(loc+5)
1375: node5=nodplc(loc+6)
1376: locm=nodplc(loc+7)
1377: locm=nodplc(locm+1)
1378: loct=nodplc(loc+19)
1379: area=value(locv+1)
1380: fnk=value(locm+10)
1381: fna=value(locm+11)
1382: c
1383: c extrinsic resistances
1384: c
1385: c... drain resistance
1386: cval=cvalue(lcvn+node1)-cvalue(lcvn+node4)
1387: vno1(kntr)=(vreal*vreal+vimag*vimag)*fourkt*value(locm+4)*area
1388: c... source resistance
1389: cval=cvalue(lcvn+node3)-cvalue(lcvn+node5)
1390: vno2(kntr)=(vreal*vreal+vimag*vimag)*fourkt*value(locm+5)*area
1391: c
1392: c drain current shot noise and flicker noise
1393: c
1394: cval=cvalue(lcvn+node4)-cvalue(lcvn+node5)
1395: vtemp=vreal*vreal+vimag*vimag
1396: vno3(kntr)=vtemp*fourkt*2.0d0*dabs(value(lx0+loct+5))/3.0d0
1397: arg=dmax1(dabs(value(lx0+loct+3)),1.0d-20)
1398: vno4(kntr)=vtemp*fnk*dexp(fna*dlog(arg))/freq
1399: vntot(kntr)=vno1(kntr)+vno2(kntr)+vno3(kntr)+vno4(kntr)
1400: vnrms=vnrms+vntot(kntr)
1401: if (kntr.ge.kntlim) go to 450
1402: 430 loc=nodplc(loc)
1403: go to 420
1404: 440 if (kntr.eq.0) go to 500
1405: 450 if (nprnt.eq.0) go to 460
1406: if (ititle.eq.0) write (6,401)
1407: ititle=1
1408: write (6,afmt1) (anam(i),i=1,kntr)
1409: write (6,afmt2) alsrd,(vno1(i),i=1,kntr)
1410: write (6,afmt2) alsrs,(vno2(i),i=1,kntr)
1411: write (6,afmt2) alsid,(vno3(i),i=1,kntr)
1412: write (6,afmt2) alsfn,(vno4(i),i=1,kntr)
1413: write (6,afmt2) alstot,(vntot(i),i=1,kntr)
1414: 460 kntr=0
1415: if (loc.ne.0) go to 430
1416: c
1417: c mosfets
1418: c
1419: 500 if (jelcnt(14).eq.0) go to 600
1420: ititle=0
1421: 501 format(//'0**** mosfet squared noise voltages (sq v/hz)')
1422: 510 loc=locate(14)
1423: kntr=0
1424: 520 if (loc.eq.0) go to 540
1425: kntr=kntr+1
1426: locv=nodplc(loc+1)
1427: anam(kntr)=value(locv)
1428: node1=nodplc(loc+2)
1429: node2=nodplc(loc+3)
1430: node3=nodplc(loc+4)
1431: node4=nodplc(loc+5)
1432: node5=nodplc(loc+6)
1433: node6=nodplc(loc+7)
1434: locm=nodplc(loc+8)
1435: itype=nodplc(locm+2)
1436: loct=nodplc(loc+26)
1437: locm=nodplc(locm+1)
1438: if(itype.eq.0) go to 522
1439: xl=value(locv+1)-2.0d0*value(locm+20)
1440: xw=value(locm+2)-2.0d0*value(locm+36)
1441: cox=value(locm+13)*xl*xw
1442: fnk=value(locm+27)
1443: fna=value(locm+28)
1444: c
1445: c extrinsic resistances
1446: c
1447: c... drain resistance
1448: 522 cval=cvalue(lcvn+node1)-cvalue(lcvn+node5)
1449: vno1(kntr)=(vreal*vreal+vimag*vimag)*fourkt*value(locv+11)
1450: c... source resistance
1451: cval=cvalue(lcvn+node3)-cvalue(lcvn+node6)
1452: vno2(kntr)=(vreal*vreal+vimag*vimag)*fourkt*value(locv+12)
1453: c
1454: c drain current shot noise and flicker noise
1455: c
1456: cval=cvalue(lcvn+node5)-cvalue(lcvn+node6)
1457: vtemp=vreal*vreal+vimag*vimag
1458: gm=value(lx0+loct+7)
1459: arg=dmax1(dabs(value(lx0+loct+4)),1.0d-20)
1460: if(itype.ne.0) go to 524
1461: modeop=value(locv+8)
1462: if(modeop.le.0) gm=value(lx0+loct+13)
1463: if(value(locm+10).ne.0.0d0) gm=value(locm+10)
1464: xnexp=value(locm+11)
1465: fnk=value(locm+12)
1466: fna=value(locm+13)
1467: vno4(kntr)=vtemp*fnk*dexp(fna*dlog(arg))/(freq**xnexp)
1468: 524 vno3(kntr)=vtemp*fourkt*dabs(gm)/1.5d0
1469: if(itype.eq.0) go to 525
1470: vno4(kntr)=vtemp*fnk*dexp(fna*dlog(arg))/(freq*cox)
1471: 525 vntot(kntr)=vno1(kntr)+vno2(kntr)+vno3(kntr)+vno4(kntr)
1472: vnrms=vnrms+vntot(kntr)
1473: if (kntr.ge.kntlim) go to 550
1474: 530 loc=nodplc(loc)
1475: go to 520
1476: 540 if (kntr.eq.0) go to 600
1477: 550 if (nprnt.eq.0) go to 560
1478: if (ititle.eq.0) write (6,501)
1479: ititle=1
1480: write (6,afmt1) (anam(i),i=1,kntr)
1481: write (6,afmt2) alsrd,(vno1(i),i=1,kntr)
1482: write (6,afmt2) alsrs,(vno2(i),i=1,kntr)
1483: write (6,afmt2) alsid,(vno3(i),i=1,kntr)
1484: write (6,afmt2) alsfn,(vno4(i),i=1,kntr)
1485: write (6,afmt2) alstot,(vntot(i),i=1,kntr)
1486: 560 kntr=0
1487: if (loc.ne.0) go to 530
1488: c
1489: c compute equivalent input noise voltage
1490: c
1491: 600 vnout=dsqrt(vnrms)
1492: vnin=vnout/vout
1493: if (nprnt.eq.0) go to 620
1494: do 610 i=1,5
1495: string(i)=ablnk
1496: 610 continue
1497: ioutyp=1
1498: ipos=1
1499: call outnam(nosout,ioutyp,string,ipos)
1500: call move(string,ipos,aslash,1,1)
1501: ipos=ipos+1
1502: locv=nodplc(nosin+1)
1503: anam1=value(locv)
1504: call move(string,ipos,anam1,1,8)
1505: write (6,611) vnrms,vnout,string,vout,anam1,vnin
1506: 611 format(////,
1507: 1 '0**** total output noise voltage',9x,'= ',1pd10.3,' sq v/hz'/,
1508: 2 1h0,40x,'= ',d10.3,' v/rt hz'/,
1509: 3 '0 transfer function value:',
1510: 4 1h0,7x,4a8,a1,'= ',d10.3,/,
1511: 5 '0 equivalent input noise at ',a8,' = ',d10.3,' /rt hz')
1512: c
1513: c save noise outputs
1514: c
1515: 620 loc=locate(44)
1516: 630 if (loc.eq.0) go to 1000
1517: iseq=nodplc(loc+4)
1518: if (nodplc(loc+5).ne.2) go to 640
1519: cvalue(loco+iseq)=vnout
1520: go to 650
1521: 640 cvalue(loco+iseq)=vnin
1522: 650 loc=nodplc(loc)
1523: go to 630
1524: c
1525: c finished
1526: c
1527: 1000 return
1528: end
1529: subroutine acasol
1530: implicit double precision (a-h,o-z)
1531: c
1532: c this routine evaluates the response of the adjoint circuit by
1533: c doing a forward/backward substitution step using the transpose of the
1534: c circuit equation coefficient matrix.
1535: c
1536: common /tabinf/ ielmnt,isbckt,nsbckt,iunsat,nunsat,itemps,numtem,
1537: 1 isens,nsens,ifour,nfour,ifield,icode,idelim,icolum,insize,
1538: 2 junode,lsbkpt,numbkp,iorder,jmnode,iur,iuc,ilc,ilr,numoff,isr,
1539: 3 nmoffc,iseq,iseq1,neqn,nodevs,ndiag,iswap,iequa,macins,lvnim1,
1540: 4 lx0,lvn,lynl,lyu,lyl,lx1,lx2,lx3,lx4,lx5,lx6,lx7,ld0,ld1,ltd,
1541: 5 imynl,imvn,lcvn,loutpt,nsnod,nsmat,nsval,icnod,icmat,icval
1542: common /cirdat/ locate(50),jelcnt(50),nunods,ncnods,numnod,nstop,
1543: 1 nut,nlt,nxtrm,ndist,ntlin,ibr,numvs
1544: common /blank/ value(1000)
1545: integer nodplc(64)
1546: complex*16 cvalue(32)
1547: equivalence (value(1),nodplc(1),cvalue(1))
1548: c
1549: c evaluates adjoint response by doing forward/backward substitution on
1550: c the transpose of the y matrix
1551: c
1552: c
1553: c forward substitution
1554: c
1555: do 20 i=2,nstop
1556: io=nodplc(iorder+i)
1557: call cdiv(value(lvn+io),value(imvn+io),value(lynl+io),
1558: 1 value(imynl+io),value(lvn+io),value(imvn+io))
1559: jstart=nodplc(iur+i)
1560: jstop=nodplc(iur+i+1)-1
1561: if (jstart.gt.jstop) go to 20
1562: if (value(lvn+io).ne.0.0d0) go to 5
1563: if (value(imvn+io).eq.0.0d0) go to 20
1564: 5 do 10 j=jstart,jstop
1565: jo=nodplc(iuc+j)
1566: jo=nodplc(iorder+jo)
1567: call cmult(value(lyu+j),value(imynl+nstop+j),value(lvn+io),
1568: 1 value(imvn+io),xreal,ximag)
1569: value(lvn+jo)=value(lvn+jo)-xreal
1570: value(imvn+jo)=value(imvn+jo)-ximag
1571: 10 continue
1572: 20 continue
1573: c
1574: c backward substitution
1575: c
1576: k=nstop+1
1577: do 40 i=2,nstop
1578: k=k-1
1579: io=nodplc(iorder+k)
1580: jstart=nodplc(ilc+k)
1581: jstop=nodplc(ilc+k+1)-1
1582: if (jstart.gt.jstop) go to 40
1583: do 30 j=jstart,jstop
1584: jo=nodplc(ilr+j)
1585: jo=nodplc(iorder+jo)
1586: call cmult(value(lyl+j),value(imynl+nstop+nut+j),value(lvn+jo),
1587: 1 value(imvn+jo),xreal,ximag)
1588: value(lvn+io)=value(lvn+io)-xreal
1589: value(imvn+io)=value(imvn+io)-ximag
1590: 30 continue
1591: 40 continue
1592: c
1593: c reorder right-hand side
1594: c
1595: do 50 i=2,nstop
1596: j=nodplc(iswap+i)
1597: value(ndiag+i)=value(lvn+j)
1598: value(ndiag+i+nstop)=value(imvn+j)
1599: 50 continue
1600: call copy8(value(ndiag+1),value(lvn+1),nstop)
1601: call copy8(value(ndiag+nstop+1),value(imvn+1),nstop)
1602: do 120 i=2,nstop
1603: cvalue(lcvn+i)=dcmplx(value(lvn+i),value(imvn+i))
1604: 120 continue
1605: cvalue(lcvn+1)=dcmplx(0.0d0,0.0d0)
1606: c
1607: c finished
1608: c
1609: return
1610: end
1611: subroutine dinit
1612: implicit double precision (a-h,o-z)
1613: c
1614: c this routine performs storage-allocation and one-time computation
1615: c needed to do the small-signal distortion analysis.
1616: c
1617: common /tabinf/ ielmnt,isbckt,nsbckt,iunsat,nunsat,itemps,numtem,
1618: 1 isens,nsens,ifour,nfour,ifield,icode,idelim,icolum,insize,
1619: 2 junode,lsbkpt,numbkp,iorder,jmnode,iur,iuc,ilc,ilr,numoff,isr,
1620: 3 nmoffc,iseq,iseq1,neqn,nodevs,ndiag,iswap,iequa,macins,lvnim1,
1621: 4 lx0,lvn,lynl,lyu,lyl,lx1,lx2,lx3,lx4,lx5,lx6,lx7,ld0,ld1,ltd,
1622: 5 imynl,imvn,lcvn,loutpt,nsnod,nsmat,nsval,icnod,icmat,icval
1623: common /cirdat/ locate(50),jelcnt(50),nunods,ncnods,numnod,nstop,
1624: 1 nut,nlt,nxtrm,ndist,ntlin,ibr,numvs
1625: common /status/ omega,time,delta,delold(7),ag(7),vt,xni,egfet,
1626: 1 xmu,mode,modedc,icalc,initf,method,iord,maxord,noncon,iterno,
1627: 2 itemno,nosolv,ipostp,iscrch
1628: common /flags/ iprnta,iprntl,iprntm,iprntn,iprnto,limtim,limpts,
1629: 1 lvlcod,lvltim,itl1,itl2,itl3,itl4,itl5,igoof,nogo,keof
1630: common /blank/ value(1000)
1631: integer nodplc(64)
1632: complex*16 cvalue(32)
1633: equivalence (value(1),nodplc(1),cvalue(1))
1634: c
1635: c
1636: call getm8(ld0,ndist)
1637: call getm16(ld1,5*nstop)
1638: c
1639: c bipolar junction transistors
1640: c
1641: loc=locate(12)
1642: 100 if (loc.eq.0) go to 200
1643: locv=nodplc(loc+1)
1644: area=value(locv+1)
1645: locm=nodplc(loc+8)
1646: locm=nodplc(locm+1)
1647: loct=lx0+nodplc(loc+22)
1648: locd=ld0+nodplc(loc+23)
1649: csat=value(locm+1)*area
1650: ova=value(locm+4)
1651: tf=value(locm+24)
1652: tr=value(locm+33)
1653: czbe=value(locm+21)*area
1654: czbc=value(locm+29)*area
1655: pe=value(locm+22)
1656: xme=value(locm+23)
1657: pc=value(locm+30)
1658: xmc=value(locm+31)
1659: fcpe=value(locm+46)
1660: fcpc=value(locm+50)
1661: vbe=value(loct)
1662: vbc=value(loct+1)
1663: gpi=value(loct+4)
1664: go=value(loct+7)
1665: gm=value(loct+6)
1666: gmu=value(loct+5)
1667: if (vbe.gt.0.0d0) go to 110
1668: evbe=1.0d0
1669: cbe=csat*vbe/vt
1670: go to 120
1671: 110 evbe=dexp(vbe/vt)
1672: cbe=csat*(evbe-1.0d0)
1673: 120 if (vbc.gt.0.0d0) go to 130
1674: evbc=1.0d0
1675: cbc=csat*vbc/vt
1676: arg=1.0d0-vbc/pc
1677: go to 140
1678: 130 evbc=dexp(vbc/vt)
1679: cbc=csat*(evbc-1.0d0)
1680: 140 if (vbe.ge.fcpe) go to 150
1681: arg=1.0d0-vbe/pe
1682: sarg=dexp(xme*dlog(arg))
1683: cjeo=czbe/sarg
1684: argbe=pe-vbe
1685: cje1=xme*cjeo/argbe
1686: cje2=xme*(1.0d0+xme)*cje1/argbe
1687: go to 160
1688: 150 denom=dexp((1.0d0+xme)*dlog(1.0d0-fcpe))
1689: cjeo=czbe*(1.0d0-fcpe*(1.0d0+xme)+xme*vbe/pe)/denom
1690: cje1=czbe*xme/(denom*pe)
1691: cje2=0.0d0
1692: 160 if (vbc.ge.fcpc) go to 170
1693: arg=1.0d0-vbc/pc
1694: sarg=dexp(xmc*dlog(arg))
1695: cjco=czbc/sarg
1696: argbc=pc-vbc
1697: cjc1=xmc*cjco/argbc
1698: cjc2=xmc*(1.0d0+xmc)*cjc1/argbc
1699: go to 180
1700: 170 denom=dexp((1.0d0+xmc)*dlog(1.0d0-fcpc))
1701: cjco=czbc*(1.0d0-fcpc*(1.0d0+xmc)+xmc*vbc/pc)/denom
1702: cjc1=czbc*xmc/(denom*pc)
1703: cjc2=0.0d0
1704: 180 twovt=vt+vt
1705: go2=(-go+csat*(evbe+evbc)*ova)/twovt
1706: gmo2=(cbe+csat)*ova/vt-2.0d0*go2
1707: gm2=(gm+go)/twovt-gmo2-go2
1708: gmu2=gmu/twovt
1709: if (vbc.le.0.0d0) gmu2=0.0d0
1710: gpi2=gpi/twovt
1711: if (vbe.le.0.0d0) gpi2=0.0d0
1712: cbo=tf*csat*evbe/vt
1713: cbor=tr*csat*evbc/vt
1714: cb1=cbo/vt
1715: cb1r=cbor/vt
1716: trivt=3.0d0*vt
1717: go3=-(go2+(cbc+csat)*ova/twovt)/trivt
1718: gmo23=-3.0d0*go3
1719: gm2o3=-gmo23+(cbe+csat)*ova/(vt*twovt)
1720: gm3=(gm2-(cbe-cbc)*ova/twovt)/trivt
1721: gmu3=gmu2/trivt
1722: gpi3=gpi2/trivt
1723: cb2=cb1/twovt
1724: cb2r=cb1r/twovt
1725: value(locd)=cje1
1726: value(locd+1)=cje2
1727: value(locd+2)=cjc1
1728: value(locd+3)=cjc2
1729: value(locd+4)=go2
1730: value(locd+5)=gmo2
1731: value(locd+6)=gm2
1732: value(locd+7)=gmu2
1733: value(locd+8)=gpi2
1734: value(locd+9)=cbo
1735: value(locd+10)=cbor
1736: value(locd+11)=cb1
1737: value(locd+12)=cb1r
1738: value(locd+13)=go3
1739: value(locd+14)=gmo23
1740: value(locd+15)=gm2o3
1741: value(locd+16)=gm3
1742: value(locd+17)=gmu3
1743: value(locd+18)=gpi3
1744: value(locd+19)=cb2
1745: value(locd+20)=cb2r
1746: loc=nodplc(loc)
1747: go to 100
1748: c
1749: c diodes
1750: c
1751: 200 loc=locate(11)
1752: 210 if (loc.eq.0) go to 300
1753: locv=nodplc(loc+1)
1754: area=value(locv+1)
1755: locm=nodplc(loc+5)
1756: locm=nodplc(locm+1)
1757: loct=lx0+nodplc(loc+11)
1758: locd=ld0+nodplc(loc+12)
1759: csat=value(locm+1)*area
1760: vte=value(locm+3)*vt
1761: tau=value(locm+4)
1762: czero=value(locm+5)*area
1763: phib=value(locm+6)
1764: xm=value(locm+7)
1765: fcpb=value(locm+12)
1766: vd=value(loct)
1767: geq=value(loct+2)
1768: evd=1.0d0
1769: if (vd.ge.0.0d0) evd=dexp(vd/vte)
1770: if (vd.ge.fcpb) go to 220
1771: arg=1.0d0-vd/phib
1772: sarg=dexp(xm*dlog(arg))
1773: cdjo=czero/sarg
1774: argd=phib-vd
1775: cdj1=xm*czero/argd
1776: cdj2=xm*(1.0d0+xm)*cdj1/argd
1777: go to 230
1778: 220 denom=dexp((1.0d0+xm)*dlog(1.0d0-fcpb))
1779: cdjo=czero*(1.0d0-fcpb*(1.0d0+xm)+xm*vd/phib)/denom
1780: cdj1=czero*xm/(denom*phib)
1781: cdj2=0.0d0
1782: cdj2=0.0d0
1783: 230 cdbo=tau*csat*evd/vte
1784: cdb1=cdbo/vte
1785: twovte=2.0d0*vte
1786: geq2=geq/twovte
1787: if (vd.le.0.0d0) geq2=0.0d0
1788: trivte=3.0d0*vte
1789: geq3=geq2/trivte
1790: cdb2=cdb1/twovte
1791: value(locd)=cdj1
1792: value(locd+1)=cdj2
1793: value(locd+2)=cdbo
1794: value(locd+3)=cdb1
1795: value(locd+4)=geq2
1796: value(locd+5)=geq3
1797: value(locd+6)=cdb2
1798: loc=nodplc(loc)
1799: go to 210
1800: c
1801: c finished
1802: c
1803: 300 return
1804: end
1805: subroutine disto(loco)
1806: implicit double precision (a-h,o-z)
1807: c
1808: c this routine performs the small-signal distortion analysis.
1809: c
1810: common /tabinf/ ielmnt,isbckt,nsbckt,iunsat,nunsat,itemps,numtem,
1811: 1 isens,nsens,ifour,nfour,ifield,icode,idelim,icolum,insize,
1812: 2 junode,lsbkpt,numbkp,iorder,jmnode,iur,iuc,ilc,ilr,numoff,isr,
1813: 3 nmoffc,iseq,iseq1,neqn,nodevs,ndiag,iswap,iequa,macins,lvnim1,
1814: 4 lx0,lvn,lynl,lyu,lyl,lx1,lx2,lx3,lx4,lx5,lx6,lx7,ld0,ld1,ltd,
1815: 5 imynl,imvn,lcvn,loutpt,nsnod,nsmat,nsval,icnod,icmat,icval
1816: common /miscel/ atime,aprog(3),adate,atitle(10),defl,defw,defad,
1817: 1 defas,rstats(50),iwidth,lwidth,nopage
1818: common /cirdat/ locate(50),jelcnt(50),nunods,ncnods,numnod,nstop,
1819: 1 nut,nlt,nxtrm,ndist,ntlin,ibr,numvs
1820: common /status/ omega,time,delta,delold(7),ag(7),vt,xni,egfet,
1821: 1 xmu,mode,modedc,icalc,initf,method,iord,maxord,noncon,iterno,
1822: 2 itemno,nosolv,ipostp,iscrch
1823: common /knstnt/ twopi,xlog2,xlog10,root2,rad,boltz,charge,ctok,
1824: 1 gmin,reltol,abstol,vntol,trtol,chgtol,eps0,epssil,epsox
1825: common /flags/ iprnta,iprntl,iprntm,iprntn,iprnto,limtim,limpts,
1826: 1 lvlcod,lvltim,itl1,itl2,itl3,itl4,itl5,igoof,nogo,keof
1827: common /ac/ fstart,fstop,fincr,skw2,refprl,spw2,jacflg,idfreq,
1828: 1 inoise,nosprt,nosout,nosin,idist,idprt
1829: common /blank/ value(1000)
1830: integer nodplc(64)
1831: complex*16 cvalue(32)
1832: equivalence (value(1),nodplc(1),cvalue(1))
1833: c
1834: c
1835: complex*16 difvn1,difvn2,difvn3,difvi1,difvi2,difvi3,dsgo2,dsgm2,
1836: 1 dsgmu2,dsgpi2,dscb1,dscb1r,dscje1,dscjc1,disto1,disto2,disto3,
1837: 2 dsgmo2,dgm2o3,dgmo23,bew,cew,bcw,be2w,ce2w,bc2w,bew2,cew2,
1838: 3 bcw2,bew12,cew12,bcw12,dscdb1,dscdj1,dsg2,cvabe,cvabc,cvace,
1839: 4 cvout,cvdist
1840: dimension distit(4)
1841: dimension vdo(2,12)
1842: complex*16 cvdo(12)
1843: equivalence (cvdo(1),vdo(1,1))
1844: data distit / 8hdistorti, 8hon analy, 8hsis , 8h /
1845: icvw1=ld1
1846: icv2w1=icvw1+nstop
1847: icvw2=icv2w1+nstop
1848: icvw12=icvw2+nstop
1849: icvadj=icvw12+nstop
1850: iprnt=0
1851: if (icalc.ge.2) go to 10
1852: idnp=nodplc(idist+2)
1853: idnn=nodplc(idist+3)
1854: locv=nodplc(idist+1)
1855: rload=1.0d0/value(locv+1)
1856: kntr=1
1857: 10 if (idprt.eq.0) go to 30
1858: if (kntr.gt.icalc) go to 30
1859: iprnt=1
1860: kntr=kntr+idprt
1861: call title(0,lwidth,1,distit)
1862: 30 freq1=dreal(cvalue(loco+1))
1863: freq2=skw2*freq1
1864: call copy16(cvalue(lcvn+1),cvalue(icvw1+1),nstop)
1865: cvout=cvalue(icvw1+idnp)-cvalue(icvw1+idnn)
1866: call magphs(cvout,omag,ophase)
1867: c
1868: c begin the distortion analysis
1869: c
1870: do 1000 kdisto=1,7
1871: cvdist=dcmplx(0.0d0,0.0d0)
1872: go to (1000,110,120,130,140,160,170),kdisto
1873: 110 freqd=2.0d0*freq1
1874: arg=dsqrt(2.0d0*rload*refprl)/(omag*omag)
1875: if (iprnt.eq.0) go to 200
1876: write (6,111) freq1,freqd,omag,ophase
1877: 111 format (///5x,'2nd harmonic distortion',30x,'freq1 = ',1pd9.2,
1878: 1 ' hz'//5x,'distortion frequency ',d9.2,' hz',16x,
1879: 2 'mag ',d9.3,3x,'phs ',0pf7.2)
1880: go to 200
1881: 120 freqd=3.0d0*freq1
1882: arg=2.0d0*rload*refprl/(omag*omag*omag)
1883: if (iprnt.eq.0) go to 200
1884: write (6,121) freq1,freqd,omag,ophase
1885: 121 format (1h1,4x,'3rd harmonic distortion',30x,'freq1 = ',1pd9.2,
1886: 1 ' hz'//5x,'distortion frequency ',d9.2,' hz',16x,
1887: 2 'mag ',d9.3,3x,'phs ',0pf7.2)
1888: go to 200
1889: 130 freqd=freq2
1890: go to 200
1891: 140 freqd=freq1-freq2
1892: arg=dsqrt(2.0d0*rload*refprl)*spw2/(omag*omag)
1893: if (iprnt.eq.0) go to 200
1894: write (6,151) freq1,freq2,freqd,omag,ophase,ow2mag,ow2phs
1895: 151 format (1h1,4x,'2nd order intermodulation difference component',
1896: 1 7x,'freq1 = ',1pd9.2,' hz',15x,'freq2 = ',d9.2,' hz'//
1897: 2 5x,'distortion frequency ',d9.2,' hz',16x,'mag ',
1898: 3 d9.3,3x,'phs ',0pf7.2,9x,'mag ',1pd9.3,3x,'phs ',0pf7.2)
1899: go to 200
1900: 160 freqd=freq1+freq2
1901: arg=dsqrt(2.0d0*rload*refprl)*spw2/(omag*omag)
1902: if (iprnt.eq.0) go to 200
1903: write (6,161) freq1,freq2,freqd,omag,ophase,ow2mag,ow2phs
1904: 161 format (1h1,4x,'2nd order intermodulation sum component',
1905: 1 14x,'freq1 = ',1pd9.2,' hz',15x,'freq2 = ',d9.2,' hz'//
1906: 2 5x,'distortion frequency ',d9.2,' hz',16x,'mag ',
1907: 3 d9.3,3x,'phs ',0pf7.2,9x,'mag ',1pd9.3,3x,'phs ',0pf7.2)
1908: go to 200
1909: 170 freqd=2.0d0*freq1-freq2
1910: arg=2.0d0*rload*refprl*spw2/(omag*omag*omag)
1911: if (iprnt.eq.0) go to 200
1912: write (6,171) freq1,freq2,freqd,omag,ophase,ow2mag,ow2phs
1913: 171 format (1h1,4x,'3rd order intermodulation difference component',
1914: 1 7x,'freq1 = ',1pd9.2,' hz',15x,'freq2 = ',d9.2,' hz'//
1915: 2 5x,'distortion frequency ',d9.2,' hz',16x,'mag ',
1916: 3 d9.3,3x,'phs ',0pf7.2,9x,'mag ',1pd9.3,3x,'phs ',0pf7.2)
1917: c
1918: c load and decompose y matrix
1919: c
1920: 200 omega=twopi*freqd
1921: igoof=0
1922: call acload
1923: call acdcmp
1924: if (igoof.eq.0) go to 220
1925: write (6,211) igoof,freqd
1926: 211 format('0warning: underflow ',i4,' time(s) in distortion analysis
1927: 1 at freq = ',1pd9.3,' hz')
1928: igoof=0
1929: 220 if (kdisto.eq.4) go to 710
1930: c
1931: c obtain adjoint solution
1932: c
1933: call zero8(value(lvn+1),nstop)
1934: call zero8(value(imvn+1),nstop)
1935: value(lvn+idnp)=-1.0d0
1936: value(lvn+idnn)=+1.0d0
1937: call acasol
1938: call copy16(cvalue(lcvn+1),cvalue(icvadj+1),nstop)
1939: call zero8(value(lvn+1),nstop)
1940: call zero8(value(imvn+1),nstop)
1941: c
1942: c bjts
1943: c
1944: if (jelcnt(12).eq.0) go to 500
1945: ititle=0
1946: 301 format (////1x,'bjt distortion components'//1x,'name',11x,'gm',
1947: 1 8x,'gpi',7x,'go',8x,'gmu',6x,'gmo2',7x,'cb',8x,'cbr',7x,'cje',
1948: 2 7x,'cjc',6x,'total')
1949: 311 format (////1x,'bjt distortion components'//1x,'name',11x,'gm',
1950: 1 8x,'gpi',7x,'go',8x,'gmu',6x,'gmo2',7x,'cb',8x,'cbr',7x,'cje',
1951: 2 7x,'cjc',6x,'gm203',5x,'gmo23',5x,'total')
1952: 320 loc=locate(12)
1953: 330 if (loc.eq.0) go to 500
1954: locv=nodplc(loc+1)
1955: loct=lx0+nodplc(loc+22)
1956: locd=ld0+nodplc(loc+23)
1957: node1=nodplc(loc+5)
1958: node2=nodplc(loc+6)
1959: node3=nodplc(loc+7)
1960: cje1=value(locd)
1961: cje2=value(locd+1)
1962: cjc1=value(locd+2)
1963: cjc2=value(locd+3)
1964: go2=value(locd+4)
1965: gmo2=value(locd+5)
1966: gm2=value(locd+6)
1967: gmu2=value(locd+7)
1968: gpi2=value(locd+8)
1969: cb1=value(locd+11)
1970: cb1r=value(locd+12)
1971: go3=value(locd+13)
1972: gmo23=value(locd+14)
1973: gm2o3=value(locd+15)
1974: gm3=value(locd+16)
1975: gmu3=value(locd+17)
1976: gpi3=value(locd+18)
1977: cb2=value(locd+19)
1978: cb2r=value(locd+20)
1979: bew=cvalue(icvw1+node2)-cvalue(icvw1+node3)
1980: cew=cvalue(icvw1+node1)-cvalue(icvw1+node3)
1981: bcw=cvalue(icvw1+node2)-cvalue(icvw1+node1)
1982: if (kdisto.eq.2) go to 370
1983: be2w=cvalue(icv2w1+node2)-cvalue(icv2w1+node3)
1984: ce2w=cvalue(icv2w1+node1)-cvalue(icv2w1+node3)
1985: bc2w=cvalue(icv2w1+node2)-cvalue(icv2w1+node1)
1986: if (kdisto.eq.3) go to 380
1987: bew2=cvalue(icvw2+node2)-cvalue(icvw2+node3)
1988: cew2=cvalue(icvw2+node1)-cvalue(icvw2+node3)
1989: bcw2=cvalue(icvw2+node2)-cvalue(icvw2+node1)
1990: if (kdisto.eq.5) go to 390
1991: if (kdisto.eq.6) go to 400
1992: bew12=cvalue(icvw12+node2)-cvalue(icvw12+node3)
1993: cew12=cvalue(icvw12+node1)-cvalue(icvw12+node3)
1994: bcw12=cvalue(icvw12+node2)-cvalue(icvw12+node1)
1995: go to 410
1996: c
1997: c calculate hd2 current generators
1998: c
1999: 370 difvn1=0.5d0*cew*cew
2000: difvn2=0.5d0*bew*bew
2001: difvn3=0.5d0*bcw*bcw
2002: dsgmo2=gmo2*0.5d0*bew*cew
2003: go to 420
2004: c
2005: c calculate hd3 current generators
2006: c
2007: 380 difvi1=0.50d0*cew*ce2w
2008: difvn1=0.25d0*cew*cew*cew
2009: difvi2=0.50d0*bew*be2w
2010: difvn2=0.25d0*bew*bew*bew
2011: difvi3=0.50d0*bcw*bc2w
2012: difvn3=0.25d0*bcw*bcw*bcw
2013: dsgmo2=gmo2*(bew*ce2w+be2w*cew)*0.5d0
2014: go to 430
2015: c
2016: c calculate im2d current generators
2017: c
2018: 390 difvn1=cew*dconjg(cew2)
2019: difvn2=bew*dconjg(bew2)
2020: difvn3=bcw*dconjg(bcw2)
2021: dsgmo2=gmo2*0.5d0*(bew*dconjg(cew2)+cew*dconjg(bew2))
2022: go to 420
2023: c
2024: c calculate im2s current generators
2025: c
2026: 400 difvn1=cew*cew2
2027: difvn2=bew*bew2
2028: difvn3=bcw*bcw2
2029: dsgmo2=gmo2*0.5d0*(bew*cew2+bew2*cew)
2030: go to 420
2031: c
2032: c calculate im3 current generators
2033: c
2034: 410 difvi1=0.5d0*(ce2w*dconjg(cew2)+cew*cew12)
2035: difvi2=0.5d0*(be2w*dconjg(bew2)+bew*bew12)
2036: difvi3=0.5d0*(bc2w*dconjg(bcw2)+bcw*bcw12)
2037: difvn1=cew*cew*dconjg(cew2)*0.75d0
2038: difvn2=bew*bew*dconjg(bew2)*0.75d0
2039: difvn3=bcw*bcw*dconjg(bcw2)*0.75d0
2040: dsgmo2=gmo2*0.5d0*(dconjg(bew2)*ce2w+bew*cew12+dconjg(cew2)*be2w+
2041: 1 cew*bew12)
2042: go to 430
2043: c
2044: 420 dsgo2=go2*difvn1
2045: dsgm2=gm2*difvn2
2046: dsgmu2=gmu2*difvn3
2047: dsgpi2=gpi2*difvn2
2048: dscb1=0.5d0*cb1*omega*dcmplx(-dimag(difvn2),dreal(difvn2))
2049: dscb1r=0.5d0*cb1r*omega*dcmplx(-dimag(difvn3),dreal(difvn3))
2050: dscje1=0.5d0*cje1*omega*dcmplx(-dimag(difvn2),dreal(difvn2))
2051: dscjc1=0.5d0*cjc1*omega*dcmplx(-dimag(difvn3),dreal(difvn3))
2052: go to 440
2053: c
2054: 430 dsgo2=2.0d0*go2*difvi1+go3*difvn1
2055: dsgm2=2.0d0*gm2*difvi2+gm3*difvn2
2056: dsgmu2=2.0d0*gmu2*difvi3+gmu3*difvn3
2057: dsgpi2=2.0d0*gpi2*difvi2+gpi3*difvn2
2058: dscb1=omega*(cb1*difvi2+cb2*difvn2/3.0d0)
2059: dscb1=dcmplx(-dimag(dscb1),dreal(dscb1))
2060: dscb1r=omega*(cb1r*difvi3+cb2r*difvn3/3.0d0)
2061: dscb1r=dcmplx(-dimag(dscb1r),dreal(dscb1r))
2062: dscje1=omega*(cje1*difvi2+cje2*difvn2/3.0d0)
2063: dscje1=dcmplx(-dimag(dscje1),dreal(dscje1))
2064: dscjc1=omega*(cjc1*difvi3+cjc2*difvn3/3.0d0)
2065: dscjc1=dcmplx(-dimag(dscjc1),dreal(dscjc1))
2066: c
2067: c determine contribution of each distortion source
2068: c
2069: 440 cvabe=cvalue(icvadj+node2)-cvalue(icvadj+node3)
2070: cvabc=cvalue(icvadj+node2)-cvalue(icvadj+node1)
2071: cvace=cvalue(icvadj+node1)-cvalue(icvadj+node3)
2072: disto1=dsgm2+dsgo2+dsgmo2
2073: disto2=dsgpi2+dscb1+dscje1
2074: disto3=dsgmu2+dscb1r+dscjc1
2075: cvdo(1)=dsgm2*cvace*arg
2076: cvdo(2)=dsgpi2*cvabe*arg
2077: cvdo(3)=dsgo2*cvace*arg
2078: cvdo(4)=dsgmu2*cvabc*arg
2079: cvdo(5)=dsgmo2*cvace*arg
2080: cvdo(6)=dscb1*cvabe*arg
2081: cvdo(7)=dscb1r*cvabc*arg
2082: cvdo(8)=dscje1*cvabe*arg
2083: cvdo(9)=dscjc1*cvabc*arg
2084: if (kdisto.eq.3) go to 450
2085: if (kdisto.eq.7) go to 460
2086: cvdo(10)=cvdo(1)+cvdo(2)+cvdo(3)+cvdo(4)+cvdo(5)+cvdo(6)+cvdo(7)+
2087: 1 cvdo(8)+cvdo(9)
2088: cvdist=cvdist+cvdo(10)
2089: if (iprnt.eq.0) go to 480
2090: do 445 j=1,10
2091: call magphs(cvdo(j),xmag,xphs)
2092: cvdo(j)=dcmplx(xmag,xphs)
2093: 445 continue
2094: if (ititle.eq.0) write (6,301)
2095: ititle=1
2096: write (6,446) value(locv),(vdo(1,j),j=1,10)
2097: 446 format(1h0,a8,'mag',1p12d10.3)
2098: write (6,447) (vdo(2,j),j=1,10)
2099: 447 format(9x,'phs',12(1x,f7.2,2x))
2100: go to 480
2101: 450 dgm2o3=gm2o3*cew*bew*bew*0.25d0
2102: dgmo23=gmo23*bew*cew*cew*0.25d0
2103: go to 470
2104: 460 dgm2o3=gm2o3*(0.5d0*bew*dconjg(bew2)*cew+0.25d0*bew*bew*
2105: 1 dconjg(cew2))
2106: dgmo23=gmo23*(0.5d0*cew*dconjg(cew2)*bew+0.25d0*cew*cew*
2107: 1 dconjg(bew2))
2108: 470 disto1=disto1+dgm2o3+dgmo23
2109: cvdo(10)=dgm2o3*cvace*arg
2110: cvdo(11)=dgmo23*cvace*arg
2111: cvdo(12)=cvdo(1)+cvdo(2)+cvdo(3)+cvdo(4)+cvdo(5)+cvdo(6)+cvdo(7)+
2112: 1 cvdo(8)+cvdo(9)+cvdo(10)+cvdo(11)
2113: cvdist=cvdist+cvdo(12)
2114: if (iprnt.eq.0) go to 480
2115: do 475 j=1,12
2116: call magphs(cvdo(j),xmag,xphs)
2117: cvdo(j)=dcmplx(xmag,xphs)
2118: 475 continue
2119: if (ititle.eq.0) write (6,311)
2120: ititle=1
2121: write (6,446) value(locv),(vdo(1,j),j=1,12)
2122: write (6,447) (vdo(2,j),j=1,12)
2123: 480 value(lvn+node1)=value(lvn+node1)
2124: 1 -dreal(disto1-disto3)
2125: value(lvn+node2)=value(lvn+node2)
2126: 1 -dreal(disto2+disto3)
2127: value(lvn+node3)=value(lvn+node3)
2128: 1 +dreal(disto1+disto2)
2129: value(imvn+node1)=value(imvn+node1)
2130: 1 -dimag(disto1-disto3)
2131: value(imvn+node2)=value(imvn+node2)
2132: 1 -dimag(disto2+disto3)
2133: value(imvn+node3)=value(imvn+node3)
2134: 1 +dimag(disto1+disto2)
2135: loc=nodplc(loc)
2136: go to 330
2137: c
2138: c junction diodes
2139: c
2140: 500 if (jelcnt(11).eq.0) go to 700
2141: ititle=0
2142: 501 format (////1x,'diode distortion components'//1x,'name',
2143: 1 11x,'geq',7x,'cb',8x,'cj',7x,'total')
2144: 510 loc=locate(11)
2145: 520 if (loc.eq.0) go to 700
2146: locv=nodplc(loc+1)
2147: node1=nodplc(loc+2)
2148: node2=nodplc(loc+3)
2149: node3=nodplc(loc+4)
2150: locm=nodplc(loc+5)
2151: locm=nodplc(locm+1)
2152: loct=lx0+nodplc(loc+11)
2153: locd=ld0+nodplc(loc+12)
2154: cdj1=value(locd)
2155: cdj2=value(locd+1)
2156: cdb1=value(locd+3)
2157: geq2=value(locd+4)
2158: geq3=value(locd+5)
2159: cdb2=value(locd+6)
2160: bew=cvalue(icvw1+node3)-cvalue(icvw1+node2)
2161: if (kdisto.eq.2) go to 540
2162: be2w=cvalue(icv2w1+node3)-cvalue(icv2w1+node2)
2163: if (kdisto.eq.3) go to 550
2164: bew2=cvalue(icvw2+node3)-cvalue(icvw2+node2)
2165: if (kdisto.eq.5) go to 560
2166: if (kdisto.eq.6) go to 570
2167: bew12=cvalue(icvw12+node3)-cvalue(icvw12+node2)
2168: go to 580
2169: c
2170: c calculate hd2 current generators
2171: c
2172: 540 difvn1=0.5d0*bew*bew
2173: go to 590
2174: c
2175: c calculate hd3 current generators
2176: c
2177: 550 difvi1=0.5d0*bew*be2w
2178: difvn1=0.25d0*bew*bew*bew
2179: go to 600
2180: c
2181: c calculate im2d current generators
2182: c
2183: 560 difvn1=bew*dconjg(bew2)
2184: go to 590
2185: c
2186: c calculate im2s current generators
2187: c
2188: 570 difvn1=bew*bew2
2189: go to 590
2190: c
2191: c calculate im3 current generators
2192: c
2193: 580 difvi1=0.5d0*(be2w*dconjg(bew2)+bew*bew12)
2194: difvn1=bew*bew*dconjg(bew2)*0.75d0
2195: go to 600
2196: 590 dsg2=geq2*difvn1
2197: dscdb1=0.5d0*cdb1*omega*dcmplx(-dimag(difvn1),dreal(difvn1))
2198: dscdj1=0.5d0*cdj1*omega*dcmplx(-dimag(difvn1),dreal(difvn1))
2199: go to 610
2200: c
2201: 600 dsg2=2.0d0*geq2*difvi1+geq3*difvn1
2202: dscdb1=omega*(cdb1*difvi1+cdb2*difvn1/3.0d0)
2203: dscdb1=dcmplx(-dimag(dscdb1),dreal(dscdb1))
2204: dscdj1=omega*(cdj1*difvi1+cdj2*difvn1/3.0d0)
2205: dscdj1=dcmplx(-dimag(dscdj1),dreal(dscdj1))
2206: c
2207: c determine contribution of each distortion source
2208: c
2209: 610 cvabe=cvalue(icvadj+node3)-cvalue(icvadj+node2)
2210: disto1=dsg2+dscdb1+dscdj1
2211: cvdo(1)=dsg2*cvabe*arg
2212: cvdo(2)=dscdb1*cvabe*arg
2213: cvdo(3)=dscdj1*cvabe*arg
2214: cvdo(4)=cvdo(1)+cvdo(2)+cvdo(3)
2215: cvdist=cvdist+cvdo(4)
2216: if (iprnt.eq.0) go to 680
2217: do 670 j=1,4
2218: call magphs(cvdo(j),xmag,xphs)
2219: cvdo(j)=dcmplx(xmag,xphs)
2220: 670 continue
2221: if (ititle.eq.0) write (6,501)
2222: ititle=1
2223: write (6,446) value(locv),(vdo(1,j),j=1,4)
2224: write (6,447) (vdo(2,j),j=1,4)
2225: 680 value(lvn+node2)=value(lvn+node2)+dreal(disto1)
2226: value(lvn+node3)=value(lvn+node3)-dreal(disto1)
2227: value(imvn+node2)=value(imvn+node2)+dimag(disto1)
2228: value(imvn+node3)=value(imvn+node3)-dimag(disto1)
2229: loc=nodplc(loc)
2230: go to 520
2231: c
2232: c obtain total distortion solution if necessary
2233: c
2234: 700 go to (1000,710,790,710,710,840,860),kdisto
2235: 710 call acsol
2236: c
2237: c store solution, print and store answers
2238: c
2239: 760 go to (1000,770,790,800,820,840,860),kdisto
2240: 770 call copy16(cvalue(lcvn+1),cvalue(icv2w1+1),nstop)
2241: call magphs(cvdist,o2mag,o2phs)
2242: if (iprnt.eq.0) go to 900
2243: o2log=20.0d0*dlog10(o2mag)
2244: write (6,781) o2mag,o2phs,o2log
2245: 781 format (///5x,'hd2 magnitude ',1pd10.3,5x,'phase ',0pf7.2,
2246: 1 5x,'= ',f7.2,' db')
2247: go to 900
2248: 790 call magphs(cvdist,o3mag,o3phs)
2249: if (iprnt.eq.0) go to 900
2250: o3log=20.0d0*dlog10(o3mag)
2251: write (6,791) o3mag,o3phs,o3log
2252: 791 format (///5x,'hd3 magnitude ',1pd10.3,5x,'phase ',0pf7.2,
2253: 1 5x,'= ',f7.2,' db')
2254: go to 900
2255: 800 call copy16(cvalue(lcvn+1),cvalue(icvw2+1),nstop)
2256: cvout=cvalue(icvw2+idnp)-cvalue(icvw2+idnn)
2257: call magphs(cvout,ow2mag,ow2phs)
2258: go to 1000
2259: 820 call copy16(cvalue(lcvn+1),cvalue(icvw12+1),nstop)
2260: 840 call magphs(cvdist,o12mag,o12phs)
2261: if (iprnt.eq.0) go to 900
2262: o12log=20.0d0*dlog10(o12mag)
2263: if (kdisto.eq.6) go to 850
2264: write (6,841) o12mag,o12phs,o12log
2265: 841 format (///5x,'im2d magnitude ',1pd10.3,5x,'phase ',0pf7.2,
2266: 1 5x,'= ',f7.2,' db')
2267: go to 900
2268: 850 write (6,851) o12mag,o12phs,o12log
2269: 851 format (///5x,'im2s magnitude ',1pd10.3,5x,'phase ',0pf7.2,
2270: 1 5x,'= ',f7.2,' db')
2271: go to 900
2272: 860 call magphs(cvdist,o21mag,o21phs)
2273: if (iprnt.eq.0) go to 900
2274: o21log=20.0d0*dlog10(o21mag)
2275: write (6,861) o21mag,o21phs,o21log
2276: 861 format (///5x,'im3 magnitude ',1pd10.3,5x,'phase ',0pf7.2,
2277: 1 5x,'= ',f7.2,' db')
2278: cma=dabs(4.0d0*o21mag*dcos((o21phs-ophase)/rad))
2279: cma=dmax1(cma,1.0d-20)
2280: cmp=dabs(4.0d0*o21mag*dsin((o21phs-ophase)/rad))
2281: cmp=dmax1(cmp,1.0d-20)
2282: cmalog=20.0d0*dlog10(cma)
2283: cmplog=20.0d0*dlog10(cmp)
2284: write (6,866)
2285: 866 format (////5x,'approximate cross modulation components')
2286: write (6,871) cma,cmalog
2287: 871 format (/5x,'cma magnitude ',1pd10.3,24x,'= ',0pf7.2,' db')
2288: write (6,881) cmp,cmplog
2289: 881 format (/5x,'cmp magnitude ',1pd10.3,24x,'= ',0pf7.2,' db')
2290: c
2291: c save distortion outputs
2292: c
2293: 900 iflag=kdisto+2
2294: if (iflag.ge.7) iflag=iflag-1
2295: loc=locate(45)
2296: 910 if (loc.eq.0) go to 1000
2297: if (nodplc(loc+5).ne.iflag) go to 920
2298: iseq=nodplc(loc+4)
2299: cvalue(loco+iseq)=cvdist
2300: 920 loc=nodplc(loc)
2301: go to 910
2302: 1000 continue
2303: c
2304: c finished
2305: c
2306: 2000 return
2307: end
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.