|
|
1.1 root 1: subroutine dcop
2: implicit double precision (a-h,o-z)
3: c
4: c
5: c this routine prints out the operating points of the nonlinear
6: c circuit elements.
7: c
8: common /tabinf/ ielmnt,isbckt,nsbckt,iunsat,nunsat,itemps,numtem,
9: 1 isens,nsens,ifour,nfour,ifield,icode,idelim,icolum,insize,
10: 2 junode,lsbkpt,numbkp,iorder,jmnode,iur,iuc,ilc,ilr,numoff,isr,
11: 3 nmoffc,iseq,iseq1,neqn,nodevs,ndiag,iswap,iequa,macins,lvnim1,
12: 4 lx0,lvn,lynl,lyu,lyl,lx1,lx2,lx3,lx4,lx5,lx6,lx7,ld0,ld1,ltd,
13: 5 imynl,imvn,lcvn,loutpt,nsnod,nsmat,nsval,icnod,icmat,icval
14: common /cirdat/ locate(50),jelcnt(50),nunods,ncnods,numnod,nstop,
15: 1 nut,nlt,nxtrm,ndist,ntlin,ibr,numvs
16: common /status/ omega,time,delta,delold(7),ag(7),vt,xni,egfet,
17: 1 xmu,mode,modedc,icalc,initf,method,iord,maxord,noncon,iterno,
18: 2 itemno,nosolv,ipostp,iscrch
19: common /knstnt/ twopi,xlog2,xlog10,root2,rad,boltz,charge,ctok,
20: 1 gmin,reltol,abstol,vntol,trtol,chgtol,eps0,epssil,epsox
21: common /miscel/ atime,aprog(3),adate,atitle(10),defl,defw,defad,
22: 1 defas,rstats(50),iwidth,lwidth,nopage
23: common /dc/ tcstar(2),tcstop(2),tcincr(2),icvflg,itcelm(2),kssop,
24: 1 kinel,kidin,kovar,kidout
25: common /ac/ fstart,fstop,fincr,skw2,refprl,spw2,jacflg,idfreq,
26: 1 inoise,nosprt,nosout,nosin,idist,idprt
27: common /blank/ value(1000)
28: integer nodplc(64)
29: complex*16 cvalue(32)
30: equivalence (value(1),nodplc(1),cvalue(1))
31: c
32: c
33: dimension optitl(4)
34: dimension anam(12),av1(12),ai1(12),req(12)
35: dimension amod(12),vd(12),cap(12)
36: dimension cb(12),cc(12),vbe(12),vbc(12),vce(12),rpi(12),
37: 1 ro(12),cpi(12),cmu(12),betadc(12),betaac(12),ft(12),
38: 2 ccs(12),cbx(12),rx(12)
39: dimension cg(12),vgs(12),vds(12),gds(12),vbs(12),cbd(12),cbs(12),
40: 2 cgsov(12),cgdov(12),cgbov(12),vth(12),vdsat(12),cd(12),gm(12),
41: 3 ccgg(12),ccgd(12),ccgs(12),ccbg(12),ccbd(12),ccbs(12),
42: 4 gmb(12)
43: dimension cgs(12),cgd(12),cgb(12),cds(12)
44: equivalence(cb(1),cg(1)),(cc(1),vgs(1)),(vbe(1),vds(1)),
45: 1(vbc(1),gds(1)),(vce(1),vbs(1)),(rpi(1),cbd(1)),
46: 2(ro(1),cbs(1)),(cpi(1),cgsov(1)),(cmu(1),cgdov(1)),
47: 3(betadc(1),cgbov(1)),(betaac(1),vth(1)),(ft(1),vdsat(1)),
48: 4(ccs(1),cd(1)),(cbx(1),ccgg(1)),(rx(1),ccgd(1))
49: equivalence(vd(1),cg(1)),(cap(1),vgs(1)),(av1(1),vds(1)),
50: 1 (ai1(1),gds(1)),(req(1),vbs(1))
51: equivalence (cgs(1),ccgg(1)),(cgd(1),ccgd(1)),(cgb(1),ccgs(1)),
52: 1 (cds(1),ccbg(1))
53: dimension afmt1(3),afmt2(2),afmt3(3),afmt4(3)
54: data optitl / 8hoperatin, 8hg point , 8hinformat, 8hion /
55: data av,avd,avbe,avbc,avce,avgs,avds,avbs / 1hv,2hvd,3hvbe,3hvbc,
56: 1 3hvce,3hvgs,3hvds,3hvbs /
57: data acntrv,acntri,asrcv,asrci,atrang,atranr,avgain,aigain /
58: 1 8hv-contrl, 8hi-contrl, 8hv-source, 8hi-source,
59: 2 8htrans-g , 8htrans-r , 8hv gain , 8hi gain /
60: data ai,aid,aib,aic,aig / 1hi,2hid,2hib,2hic,2hig /
61: data areq,arpi,aro / 3hreq,3hrpi,2hro /
62: data acap,acpi,acmu,acgs,acgd,acbd,acbs / 3hcap,3hcpi,3hcmu,3hcgs,
63: 1 3hcgd,3hcbd,3hcbs /
64: data acgsov,acgdov,acgbov /6hcgsovl,6hcgdovl,6hcgbovl/
65: data accgg,accgd,accgs,accbg,accbd,accbs /7hdqgdvgb,7hdqgdvdb,
66: 1 7hdqgdvsb,7hdqbdvgb,7hdqbdvdb,7hdqbdvsb/
67: data acgb,acds / 3hcgb,3hcds /
68: data avth, avdsat / 3hvth, 5hvdsat /
69: data agm,agds / 2hgm,3hgds /
70: data agmb / 4hgmb /
71: data accs,acbx,arx /3hccs,3hcbx,2hrx/
72: data abetad,abetaa / 6hbetadc,6hbetaac /
73: data aft / 2hft /
74: c
75: data ablnk /1h /
76: data afmt1 /8h(//1h0,1,8h0x, (2x,8h,a8)) /
77: data afmt2 /8h(1h ,a8,,8h f10.3)/
78: data afmt3 /8h(1h ,a8,,8h1p d10.,8h2) /
79: data afmt4 /8h('0model,8h ', (,8h2x,a8)) /
80: c
81: c.. fix-up the format statements
82: c
83: kntr=12
84: if(lwidth.le.80) kntr=7
85: ipos=12
86: call move(afmt1,ipos,ablnk,1,2)
87: call alfnum(kntr,afmt1,ipos)
88: ipos=9
89: call move(afmt2,ipos,ablnk,1,2)
90: call alfnum(kntr,afmt2,ipos)
91: ipos=11
92: call move(afmt3,ipos,ablnk,1,2)
93: call alfnum(kntr,afmt3,ipos)
94: ipos=14
95: call move(afmt4,ipos,ablnk,1,2)
96: call alfnum(kntr,afmt4,ipos)
97: c
98: c compute voltage source currents and power dissipation
99: c
100: call second(t1)
101: if ((mode.eq.1).and.(modedc.eq.2).and.(nosolv.ne.0)) go to 700
102: power=0.0d0
103: if (jelcnt(9).eq.0) go to 50
104: ititle=0
105: 11 format (////5x,'voltage source currents'//5x,'name',
106: 1 7x,'current'/)
107: loc=locate(9)
108: 20 if (loc.eq.0) go to 50
109: locv=nodplc(loc+1)
110: iptr=nodplc(loc+6)
111: creal=value(lvnim1+iptr)
112: power=power-creal*value(locv+1)
113: if (ititle.eq.0) write (6,11)
114: ititle=1
115: write (6,21) value(locv),creal
116: 21 format (/5x,a8,1x,1pd10.3)
117: 30 loc=nodplc(loc)
118: go to 20
119: 50 loc=locate(10)
120: 60 if (loc.eq.0) go to 90
121: locv=nodplc(loc+1)
122: node1=nodplc(loc+2)
123: node2=nodplc(loc+3)
124: power=power-value(locv+1)
125: 1 *(value(lvnim1+node1)-value(lvnim1+node2))
126: loc=nodplc(loc)
127: go to 60
128: 90 write (6,91) power
129: 91 format (//5x,'total power dissipation ',1pd9.2,' watts')
130: c
131: c small signal device parameters
132: c
133: numdev=jelcnt(5)+jelcnt(6)+jelcnt(7)+jelcnt(8)+jelcnt(11)
134: 1 +jelcnt(12)+jelcnt(13)+jelcnt(14)
135: if (numdev.eq.0) go to 600
136: call title(0,lwidth,1,optitl)
137: kntlim=lwidth/11
138: c
139: c nonlinear voltage controlled current sources
140: c
141: if (jelcnt(5).eq.0) go to 175
142: ititle=0
143: 111 format(1h0,/,'0**** voltage-controlled current sources')
144: loc=locate(5)
145: kntr=0
146: 120 if (loc.eq.0) go to 140
147: kntr=kntr+1
148: locv=nodplc(loc+1)
149: loct=lx0+nodplc(loc+12)
150: anam(kntr)=value(locv)
151: ai1(kntr)=value(loct)
152: if (kntr.ge.kntlim) go to 150
153: 130 loc=nodplc(loc)
154: go to 120
155: 140 if (kntr.eq.0) go to 175
156: 150 if (ititle.eq.0) write (6,111)
157: ititle=1
158: write (6,afmt1) (anam(i),i=1,kntr)
159: write (6,afmt3) asrci,(ai1(i),i=1,kntr)
160: kntr=0
161: if (loc.ne.0) go to 130
162: c
163: c nonlinear voltage controlled voltage sources
164: c
165: 175 if (jelcnt(6).eq.0) go to 186
166: ititle=0
167: 176 format(1h0,/,'0**** voltage-controlled voltage sources')
168: loc=locate(6)
169: kntr=0
170: 178 if (loc.eq.0) go to 182
171: kntr=kntr+1
172: locv=nodplc(loc+1)
173: loct=lx0+nodplc(loc+13)
174: anam(kntr)=value(locv)
175: av1(kntr)=value(loct)
176: ai1(kntr)=value(loct+1)
177: if (kntr.ge.kntlim) go to 184
178: 180 loc=nodplc(loc)
179: go to 178
180: 182 if (kntr.eq.0) go to 186
181: 184 if (ititle.eq.0) write (6,176)
182: ititle=1
183: write (6,afmt1) (anam(i),i=1,kntr)
184: write (6,afmt2) asrcv,(av1(i),i=1,kntr)
185: write (6,afmt3) asrci,(ai1(i),i=1,kntr)
186: kntr=0
187: if (loc.ne.0) go to 180
188: c
189: c nonlinear current controlled current sources
190: c
191: 186 if (jelcnt(7).eq.0) go to 196
192: ititle=0
193: 187 format(1h0,/,'0**** current-controlled current sources')
194: loc=locate(7)
195: kntr=0
196: 188 if (loc.eq.0) go to 192
197: kntr=kntr+1
198: locv=nodplc(loc+1)
199: loct=lx0+nodplc(loc+12)
200: anam(kntr)=value(locv)
201: ai1(kntr)=value(loct)
202: if (kntr.ge.kntlim) go to 194
203: 190 loc=nodplc(loc)
204: go to 188
205: 192 if (kntr.eq.0) go to 196
206: 194 if (ititle.eq.0) write (6,187)
207: ititle=1
208: write (6,afmt1) (anam(i),i=1,kntr)
209: write (6,afmt3) asrci,(ai1(i),i=1,kntr)
210: kntr=0
211: if (loc.ne.0) go to 190
212: c
213: c nonlinear current controlled voltage sources
214: c
215: 196 if (jelcnt(8).eq.0) go to 210
216: ititle=0
217: 197 format(1h0,/,'0**** current-controlled voltage sources')
218: loc=locate(8)
219: kntr=0
220: 198 if (loc.eq.0) go to 202
221: kntr=kntr+1
222: locv=nodplc(loc+1)
223: loct=lx0+nodplc(loc+13)
224: anam(kntr)=value(locv)
225: av1(kntr)=value(loct)
226: ai1(kntr)=value(loct+1)
227: if (kntr.ge.kntlim) go to 204
228: 200 loc=nodplc(loc)
229: go to 198
230: 202 if (kntr.eq.0) go to 210
231: 204 if (ititle.eq.0) write (6,197)
232: ititle=1
233: write (6,afmt1) (anam(i),i=1,kntr)
234: write (6,afmt2) asrcv,(av1(i),i=1,kntr)
235: write (6,afmt3) asrci,(ai1(i),i=1,kntr)
236: kntr=0
237: if (loc.ne.0) go to 200
238: c
239: c diodes
240: c
241: 210 if (jelcnt(11).eq.0) go to 300
242: ititle=0
243: 211 format(1h0,/,'0**** diodes')
244: loc=locate(11)
245: kntr=0
246: 220 if (loc.eq.0) go to 240
247: kntr=kntr+1
248: locv=nodplc(loc+1)
249: node1=nodplc(loc+2)
250: node2=nodplc(loc+3)
251: locm=nodplc(loc+5)
252: locm=nodplc(locm+1)
253: loct=lx0+nodplc(loc+11)
254: anam(kntr)=value(locv)
255: amod(kntr)=value(locm)
256: cd(kntr)=value(loct+1)
257: vd(kntr)=value(lvnim1+node1)-value(lvnim1+node2)
258: if (modedc.ne.1) go to 225
259: req(kntr)=1.0d0/value(loct+2)
260: cap(kntr)=value(loct+4)
261: 225 if (kntr.ge.kntlim) go to 250
262: 230 loc=nodplc(loc)
263: go to 220
264: 240 if (kntr.eq.0) go to 300
265: 250 if (ititle.eq.0) write (6,211)
266: ititle=1
267: write (6,afmt1) (anam(i),i=1,kntr)
268: write (6,afmt4) (amod(i),i=1,kntr)
269: write (6,afmt3) aid,(cd(i),i=1,kntr)
270: write (6,afmt2) avd,(vd(i),i=1,kntr)
271: if (modedc.ne.1) go to 260
272: write (6,afmt3) areq,(req(i),i=1,kntr)
273: write (6,afmt3) acap,(cap(i),i=1,kntr)
274: 260 kntr=0
275: if (loc.ne.0) go to 230
276: c
277: c bipolar junction transistors
278: c
279: 300 if (jelcnt(12).eq.0) go to 400
280: ititle=0
281: 301 format(1h0,/,'0**** bipolar junction transistors')
282: loc=locate(12)
283: kntr=0
284: 320 if (loc.eq.0) go to 340
285: kntr=kntr+1
286: locv=nodplc(loc+1)
287: node1=nodplc(loc+2)
288: node2=nodplc(loc+3)
289: node3=nodplc(loc+4)
290: locm=nodplc(loc+8)
291: type=nodplc(locm+2)
292: locm=nodplc(locm+1)
293: loct=lx0+nodplc(loc+22)
294: anam(kntr)=value(locv)
295: amod(kntr)=value(locm)
296: cb(kntr)=type*value(loct+3)
297: cc(kntr)=type*value(loct+2)
298: vbe(kntr)=value(lvnim1+node2)-value(lvnim1+node3)
299: vbc(kntr)=value(lvnim1+node2)-value(lvnim1+node1)
300: vce(kntr)=vbe(kntr)-vbc(kntr)
301: betadc(kntr)=cc(kntr)/dsign(dmax1(dabs(cb(kntr)),1.0d-20),
302: 1 cb(kntr))
303: if (modedc.ne.1) go to 325
304: rx(kntr)=0.0d0
305: if(value(loct+16).ne.0.0d0) rx(kntr)=1.0d0/value(loct+16)
306: ccs(kntr)=value(loct+13)
307: cbx(kntr)=value(loct+15)
308: rpi(kntr)=1.0d0/value(loct+4)
309: gm(kntr)=value(loct+6)
310: ro(kntr)=1.0d0/value(loct+7)
311: cpi(kntr)=value(loct+9)
312: cmu(kntr)=value(loct+11)
313: betaac(kntr)=gm(kntr)*rpi(kntr)
314: ft(kntr)=gm(kntr)/(twopi*dmax1(cpi(kntr)+cmu(kntr)+cbx(kntr),
315: 1 1.0d-20))
316: 325 if (kntr.ge.kntlim) go to 350
317: 330 loc=nodplc(loc)
318: go to 320
319: 340 if (kntr.eq.0) go to 400
320: 350 if (ititle.eq.0) write (6,301)
321: ititle=1
322: write (6,afmt1) (anam(i),i=1,kntr)
323: write (6,afmt4) (amod(i),i=1,kntr)
324: write (6,afmt3) aib,(cb(i),i=1,kntr)
325: write (6,afmt3) aic,(cc(i),i=1,kntr)
326: write (6,afmt2) avbe,(vbe(i),i=1,kntr)
327: write (6,afmt2) avbc,(vbc(i),i=1,kntr)
328: write (6,afmt2) avce,(vce(i),i=1,kntr)
329: write (6,afmt2) abetad,(betadc(i),i=1,kntr)
330: if (modedc.ne.1) go to 360
331: write (6,afmt3) agm,(gm(i),i=1,kntr)
332: write (6,afmt3) arpi,(rpi(i),i=1,kntr)
333: write(6,afmt3) arx,(rx(i),i=1,kntr)
334: write (6,afmt3) aro,(ro(i),i=1,kntr)
335: write (6,afmt3) acpi,(cpi(i),i=1,kntr)
336: write (6,afmt3) acmu,(cmu(i),i=1,kntr)
337: write(6,afmt3) acbx,(cbx(i),i=1,kntr)
338: write(6,afmt3) accs,(ccs(i),i=1,kntr)
339: write (6,afmt2) abetaa,(betaac(i),i=1,kntr)
340: write (6,afmt3) aft,(ft(i),i=1,kntr)
341: 360 kntr=0
342: if (loc.ne.0) go to 330
343: c
344: c jfets
345: c
346: 400 if (jelcnt(13).eq.0) go to 500
347: ititle=0
348: 401 format(1h0,/,'0**** jfets')
349: loc=locate(13)
350: kntr=0
351: 420 if (loc.eq.0) go to 440
352: kntr=kntr+1
353: locv=nodplc(loc+1)
354: node1=nodplc(loc+2)
355: node2=nodplc(loc+3)
356: node3=nodplc(loc+4)
357: locm=nodplc(loc+7)
358: type=nodplc(locm+2)
359: locm=nodplc(locm+1)
360: loct=lx0+nodplc(loc+19)
361: anam(kntr)=value(locv)
362: amod(kntr)=value(locm)
363: cd(kntr)=type*(value(loct+3)-value(loct+4))
364: vgs(kntr)=value(lvnim1+node2)-value(lvnim1+node3)
365: vds(kntr)=value(lvnim1+node1)-value(lvnim1+node3)
366: if (modedc.ne.1) go to 425
367: gm(kntr)=value(loct+5)
368: gds(kntr)=value(loct+6)
369: cgs(kntr)=value(loct+9)
370: cgd(kntr)=value(loct+11)
371: 425 if (kntr.ge.kntlim) go to 450
372: 430 loc=nodplc(loc)
373: go to 420
374: 440 if (kntr.eq.0) go to 500
375: 450 if (ititle.eq.0) write (6,401)
376: ititle=1
377: write (6,afmt1) (anam(i),i=1,kntr)
378: write (6,afmt4) (amod(i),i=1,kntr)
379: write (6,afmt3) aid,(cd(i),i=1,kntr)
380: write (6,afmt2) avgs,(vgs(i),i=1,kntr)
381: write (6,afmt2) avds,(vds(i),i=1,kntr)
382: if (modedc.ne.1) go to 460
383: write (6,afmt3) agm,(gm(i),i=1,kntr)
384: write (6,afmt3) agds,(gds(i),i=1,kntr)
385: write (6,afmt3) acgs,(cgs(i),i=1,kntr)
386: write (6,afmt3) acgd,(cgd(i),i=1,kntr)
387: 460 kntr=0
388: if (loc.ne.0) go to 430
389: c
390: c mosfets
391: c
392: 500 if (jelcnt(14).eq.0) go to 600
393: ititle=0
394: 501 format(1h0,/,'0**** mosfets')
395: loc=locate(14)
396: kntr=0
397: 520 if (loc.eq.0) go to 540
398: kntr=kntr+1
399: locv=nodplc(loc+1)
400: node1=nodplc(loc+2)
401: node2=nodplc(loc+3)
402: node3=nodplc(loc+4)
403: node4=nodplc(loc+5)
404: node5=nodplc(loc+6)
405: node6=nodplc(loc+7)
406: locm=nodplc(loc+8)
407: type=nodplc(locm+2)
408: locm=nodplc(locm+1)
409: loct=lx0+nodplc(loc+26)
410: anam(kntr)=value(locv)
411: amod(kntr)=value(locm)
412: if(type.eq.0.0d0) go to 522
413: cd(kntr)=type*value(loct+4)
414: vgs(kntr)=value(lvnim1+node2)-value(lvnim1+node3)
415: vds(kntr)=value(lvnim1+node1)-value(lvnim1+node3)
416: vbs(kntr)=value(lvnim1+node4)-value(lvnim1+node3)
417: if (modedc.ne.1) go to 525
418: xl=value(locv+1)-2.0d0*value(locm+20)
419: xw=value(locv+2)-2.0d0*value(locm+36)
420: covlgs=value(locm+8)*xw
421: covlgd=value(locm+9)*xw
422: covlgb=value(locm+10)*xl
423: devmod=value(locv+8)
424: vdsat(kntr)=value(locv+10)
425: vth(kntr)=value(locv+9)+type*(value(lvnim1+node4)-
426: 1 value(lvnim1+node6))
427: gds(kntr)=value(loct+1)
428: gm(kntr)=value(loct)
429: gmb(kntr)=-(value(loct)+value(loct+1)+value(loct+2))
430: if(devmod.gt.0.0d0) go to 521
431: gds(kntr)=value(loct+2)
432: vth(kntr)=value(locv+9)+type*(value(lvnim1+node4)-
433: 1 value(lvnim1+node5))
434: 521 cbd(kntr)=value(loct+14)
435: cbs(kntr)=value(loct+15)
436: cgsov(kntr)=covlgs
437: cgdov(kntr)=covlgd
438: cgbov(kntr)=covlgb
439: ccgg(kntr)=value(loct+8)
440: ccgd(kntr)=value(loct+9)
441: ccgs(kntr)=value(loct+10)
442: ccbg(kntr)=value(loct+11)
443: ccbd(kntr)=value(loct+12)
444: ccbs(kntr)=value(loct+13)
445: go to 525
446: c... special case for ga-as
447: 522 cd(kntr)=value(loct+4)
448: cg(kntr)=value(loct+5)
449: vgs(kntr)=value(lvnim1+node2)-value(lvnim1+node3)
450: vds(kntr)=value(lvnim1+node1)-value(lvnim1+node3)
451: vbs(kntr)=value(lvnim1+node4)-value(lvnim1+node3)
452: if(modedc.ne.1) go to 525
453: modeop=value(locv+8)
454: gm(kntr)=value(loct+7)
455: gds(kntr)=(value(loct+8)*value(loct+11))
456: 1 /(value(loct+8)+value(loct+11))
457: if(modeop.le.0) gm(kntr)=value(loct+13)
458: cds(kntr)=value(loct+10)
459: cgs(kntr)=value(loct+12)
460: cgd(kntr)=value(loct+14)
461: cgb(kntr)=value(loct+16)
462: 525 if (kntr.ge.kntlim) go to 550
463: 530 loc=nodplc(loc)
464: go to 520
465: 540 if (kntr.eq.0) go to 600
466: 550 if (ititle.eq.0) write (6,501)
467: ititle=1
468: write (6,afmt1) (anam(i),i=1,kntr)
469: write (6,afmt4) (amod(i),i=1,kntr)
470: if(type.eq.0.0d0) go to 555
471: write (6,afmt3) aid,(cd(i),i=1,kntr)
472: write (6,afmt2) avgs,(vgs(i),i=1,kntr)
473: write (6,afmt2) avds,(vds(i),i=1,kntr)
474: write (6,afmt2) avbs,(vbs(i),i=1,kntr)
475: if (modedc.ne.1) go to 560
476: write (6,afmt2) avth,(vth(i),i=1,kntr)
477: write (6,afmt2) avdsat,(vdsat(i),i=1,kntr)
478: write (6,afmt3) agm,(gm(i),i=1,kntr)
479: write (6,afmt3) agds,(gds(i),i=1,kntr)
480: write (6,afmt3) agmb,(gmb(i),i=1,kntr)
481: write (6,afmt3) acbd,(cbd(i),i=1,kntr)
482: write (6,afmt3) acbs,(cbs(i),i=1,kntr)
483: write (6,afmt3) acgsov,(cgsov(i),i=1,kntr)
484: write (6,afmt3) acgdov,(cgdov(i),i=1,kntr)
485: write (6,afmt3) acgbov,(cgbov(i),i=1,kntr)
486: write(6,551)
487: 551 format(' derivatives of gate (dqgdvx) and bulk (dqbdvx) charges')
488: write (6,afmt3) accgg,(ccgg(i),i=1,kntr)
489: write (6,afmt3) accgd,(ccgd(i),i=1,kntr)
490: write (6,afmt3) accgs,(ccgs(i),i=1,kntr)
491: write (6,afmt3) accbg,(ccbg(i),i=1,kntr)
492: write (6,afmt3) accbd,(ccbd(i),i=1,kntr)
493: write (6,afmt3) accbs,(ccbs(i),i=1,kntr)
494: go to 560
495: 555 write(6,afmt3) aid,(cd(i),i=1,kntr)
496: write(6,afmt3) aig,(cg(i),i=1,kntr)
497: write (6,afmt2) avgs,(vgs(i),i=1,kntr)
498: write (6,afmt2) avds,(vds(i),i=1,kntr)
499: write (6,afmt2) avbs,(vbs(i),i=1,kntr)
500: if (modedc.ne.1) go to 560
501: write (6,afmt3) agm,(gm(i),i=1,kntr)
502: write (6,afmt3) agds,(gds(i),i=1,kntr)
503: write (6,afmt3) acgs,(cgs(i),i=1,kntr)
504: write (6,afmt3) acgd,(cgd(i),i=1,kntr)
505: write (6,afmt3) acgb,(cgb(i),i=1,kntr)
506: write (6,afmt3) acds,(cds(i),i=1,kntr)
507: 560 kntr=0
508: if (loc.ne.0) go to 530
509: 600 if (modedc.ne.1) go to 700
510: if (kinel.eq.0) go to 610
511: call sstf
512: 610 if (nsens.eq.0) go to 700
513: call sencal
514: c
515: c finished
516: c
517: 700 if (modedc.eq.2) go to 710
518: if (jacflg.ne.0) go to 705
519: call clrmem(lvnim1)
520: call clrmem(lx0)
521: 705 call clrmem(lvn)
522: call clrmem(ndiag)
523: 710 call second(t2)
524: rstats(5)=rstats(5)+t2-t1
525: return
526: end
527: subroutine sstf
528: implicit double precision (a-h,o-z)
529: c
530: c this routine computes the value of the small-signal transfer
531: c function specified by the user.
532: c
533: common /tabinf/ ielmnt,isbckt,nsbckt,iunsat,nunsat,itemps,numtem,
534: 1 isens,nsens,ifour,nfour,ifield,icode,idelim,icolum,insize,
535: 2 junode,lsbkpt,numbkp,iorder,jmnode,iur,iuc,ilc,ilr,numoff,isr,
536: 3 nmoffc,iseq,iseq1,neqn,nodevs,ndiag,iswap,iequa,macins,lvnim1,
537: 4 lx0,lvn,lynl,lyu,lyl,lx1,lx2,lx3,lx4,lx5,lx6,lx7,ld0,ld1,ltd,
538: 5 imynl,imvn,lcvn,loutpt,nsnod,nsmat,nsval,icnod,icmat,icval
539: common /cirdat/ locate(50),jelcnt(50),nunods,ncnods,numnod,nstop,
540: 1 nut,nlt,nxtrm,ndist,ntlin,ibr,numvs
541: common /flags/ iprnta,iprntl,iprntm,iprntn,iprnto,limtim,limpts,
542: 1 lvlcod,lvltim,itl1,itl2,itl3,itl4,itl5,igoof,nogo,keof
543: common /dc/ tcstar(2),tcstop(2),tcincr(2),icvflg,itcelm(2),kssop,
544: 1 kinel,kidin,kovar,kidout
545: common /blank/ value(1000)
546: integer nodplc(64)
547: complex*16 cvalue(32)
548: equivalence (value(1),nodplc(1),cvalue(1))
549: c
550: c
551: dimension string(5),save(3)
552: data aslash, ablnk / 1h/, 1h /
553: c
554: c setup current vector for input resistance and transfer function
555: c
556: call zero8(value(lvn+1),nstop)
557: if (kidin.eq.10) go to 5
558: c... voltage source input
559: iptri=nodplc(kinel+6)
560: value(lvn+iptri)=+1.0d0
561: go to 10
562: c... current source input
563: 5 noposi=nodplc(kinel+2)
564: nonegi=nodplc(kinel+3)
565: value(lvn+noposi)=-1.0d0
566: value(lvn+nonegi)=+1.0d0
567: c
568: c lu decompose and solve the system of circuit equations
569: c
570: c... reorder the right-hand side
571: 10 do 15 i=2,nstop
572: j=nodplc(iswap+i)
573: value(ndiag+i)=value(lvn+j)
574: 15 continue
575: call copy8(value(ndiag+1),value(lvn+1),nstop)
576: 20 call dcdcmp
577: call dcsol
578: value(lvn+1)=0.0d0
579: c
580: c evaluate transfer function
581: c
582: if (nodplc(kovar+5).ne.0) go to 30
583: c... voltage output
584: noposo=nodplc(kovar+2)
585: nonego=nodplc(kovar+3)
586: trfn=value(lvn+noposo)-value(lvn+nonego)
587: go to 40
588: c... current output (through voltage source)
589: 30 iptro=nodplc(kovar+2)
590: iptro=nodplc(iptro+6)
591: trfn=value(lvn+iptro)
592: c
593: c evaluate input resistance
594: c
595: 40 if (kidin.eq.9) go to 50
596: c... current source input
597: zin=value(lvn+nonegi)-value(lvn+noposi)
598: go to 70
599: c... voltage source input
600: 50 creal=value(lvn+iptri)
601: if (dabs(creal).ge.1.0d-20) go to 60
602: zin=1.0d20
603: go to 70
604: 60 zin=-1.0d0/creal
605: c
606: c setup current vector for output resistance
607: c
608: 70 call zero8(value(lvn+1),nstop)
609: if (nodplc(kovar+5).ne.0) go to 80
610: c... voltage output
611: value(lvn+noposo)=-1.0d0
612: value(lvn+nonego)=+1.0d0
613: go to 90
614: 80 if (nodplc(kovar+2).ne.kinel) go to 85
615: zout=zin
616: go to 200
617: c... current output (through voltage source)
618: 85 value(lvn+iptro)=+1.0d0
619: c
620: c perform new forward and backward substitution
621: c
622: c... reorder the right-hand side
623: 90 do 95 i=2,nstop
624: j=nodplc(iswap+i)
625: value(ndiag+i)=value(lvn+j)
626: 95 continue
627: call copy8(value(ndiag+1),value(lvn+1),nstop)
628: call dcsol
629: value(lvn+1)=0.0d0
630: c
631: c evaluate output resistance
632: c
633: 100 if (nodplc(kovar+5).ne.0) go to 110
634: c... voltage output
635: zout=value(lvn+nonego)-value(lvn+noposo)
636: go to 200
637: c... current output (through voltage source)
638: 110 creal=value(lvn+iptro)
639: if (dabs(creal).ge.1.0d-20) go to 120
640: zout=1.0d20
641: go to 200
642: 120 zout=-1.0d0/creal
643: c
644: c print results
645: c
646: 200 do 210 i=1,5
647: string(i)=ablnk
648: 210 continue
649: ipos=1
650: call outnam(kovar,1,string,ipos)
651: call copy8(string,save,3)
652: call move(string,ipos,aslash,1,1)
653: ipos=ipos+1
654: locv=nodplc(kinel+1)
655: anam=value(locv)
656: call move(string,ipos,anam,1,8)
657: write (6,231) string,trfn,anam,zin,save,zout
658: 231 format(////,'0**** small-signal characteristics'//,
659: 1 1h0,5x,5a8,3h = ,1pd10.3,/,
660: 2 1h0,5x,'input resistance at ',a8,12x,3h = ,d10.3,/,
661: 3 1h0,5x,'output resistance at ',2a8,a3,3h = ,d10.3)
662: return
663: end
664: subroutine sencal
665: implicit double precision (a-h,o-z)
666: c
667: c this routine computes the dc sensitivities of circuit elements
668: c with respect to user specified outputs.
669: c
670: common /tabinf/ ielmnt,isbckt,nsbckt,iunsat,nunsat,itemps,numtem,
671: 1 isens,nsens,ifour,nfour,ifield,icode,idelim,icolum,insize,
672: 2 junode,lsbkpt,numbkp,iorder,jmnode,iur,iuc,ilc,ilr,numoff,isr,
673: 3 nmoffc,iseq,iseq1,neqn,nodevs,ndiag,iswap,iequa,macins,lvnim1,
674: 4 lx0,lvn,lynl,lyu,lyl,lx1,lx2,lx3,lx4,lx5,lx6,lx7,ld0,ld1,ltd,
675: 5 imynl,imvn,lcvn,loutpt,nsnod,nsmat,nsval,icnod,icmat,icval
676: common /cirdat/ locate(50),jelcnt(50),nunods,ncnods,numnod,nstop,
677: 1 nut,nlt,nxtrm,ndist,ntlin,ibr,numvs
678: common /status/ omega,time,delta,delold(7),ag(7),vt,xni,egfet,
679: 1 xmu,mode,modedc,icalc,initf,method,iord,maxord,noncon,iterno,
680: 2 itemno,nosolv,ipostp,iscrch
681: common /miscel/ atime,aprog(3),adate,atitle(10),defl,defw,defad,
682: 1 defas,rstats(50),iwidth,lwidth,nopage
683: common /flags/ iprnta,iprntl,iprntm,iprntn,iprnto,limtim,limpts,
684: 1 lvlcod,lvltim,itl1,itl2,itl3,itl4,itl5,igoof,nogo,keof
685: common /dc/ tcstar(2),tcstop(2),tcincr(2),icvflg,itcelm(2),kssop,
686: 1 kinel,kidin,kovar,kidout
687: common /blank/ value(1000)
688: integer nodplc(64)
689: complex*16 cvalue(32)
690: equivalence (value(1),nodplc(1),cvalue(1))
691: c
692: c
693: dimension string(5),sentit(4)
694: data alsrs,alsis,alsn,alsrb,alsrc,alsre / 2hrs,2his,1hn,2hrb,2hrc,
695: 1 2hre /
696: data alsbf,alsc2,alsbr,alsc4,alsne,alsnc,alsik,alsikr,alsva,alsvb
697: 1 / 2hbf,3hjle,2hbr,3hjlc,3hnle,3hnlc,3hjbf,3hjbr,3hvbf,3hvbr/
698: data alsjs /2hjs/
699: data sentit / 8hdc sensi, 8htivity a, 8hnalysis , 8h /
700: data ablnk / 1h /
701: c
702: c
703: if (kinel.ne.0) go to 8
704: 4 call dcdcmp
705: c
706: c
707: 8 do 1000 n=1,nsens
708: c
709: c prepare adjoint excitation vector
710: c
711: call zero8(value(lvn+1),nstop)
712: locs=nodplc(isens+n)
713: ioutyp=nodplc(locs+5)
714: if (ioutyp.ne.0) go to 10
715: c... voltage output
716: ivolts=1
717: noposo=nodplc(locs+2)
718: nonego=nodplc(locs+3)
719: value(lvn+noposo)=-1.0d0
720: value(lvn+nonego)=+1.0d0
721: go to 20
722: c... current output (through voltage source)
723: 10 iptro=nodplc(locs+2)
724: ivolts=0
725: iptro=nodplc(iptro+6)
726: value(lvn+iptro)=-1.0d0
727: c
728: c obtain adjoint solution by doing forward/backward substitution on
729: c the transpose of the y matrix
730: c
731: 20 call asol
732: value(lvn+1)=0.0d0
733: c
734: c real solution in lvnim1; adjoint solution in lvn ...
735: c
736: call title(0,lwidth,1,sentit)
737: ipos=1
738: call outnam(locs,1,string,ipos)
739: call move(string,ipos,ablnk,1,7)
740: jstop=(ipos+6)/8
741: write (6,36) (string(j),j=1,jstop)
742: 36 format('0dc sensitivities of output ',5a8)
743: if(ivolts.ne.0) write (6,41)
744: if(ivolts.eq.0) write(6,42)
745: 41 format(1h0,8x,'element',9x,'element',7x,'element',7x,'normalized'/
746: 1 10x,'name',12x,'value',6x,'sensitivity sensitivity'/35x,
747: 2 ' (volts/unit) (volts/percent)'/)
748: 42 format(1h0,8x,'element',9x,'element',7x,'element',7x,'normalized'/
749: 1 10x,'name',12x,'value',6x,'sensitivity sensitivity'/35x,
750: 2 ' (amps/unit) (amps/percent)'/)
751: c
752: c resistors
753: c
754: loc=locate(1)
755: 100 if (loc.eq.0) go to 110
756: locv=nodplc(loc+1)
757: node1=nodplc(loc+2)
758: node2=nodplc(loc+3)
759: val=1.0d0/value(locv+1)
760: sens=-(value(lvnim1+node1)-value(lvnim1+node2))*
761: 1 (value(lvn +node1)-value(lvn +node2))/(val*val)
762: sensn=val*sens/100.0d0
763: write (6,101) value(locv),val,sens,sensn
764: 101 format(10x,a8,4x,1pd10.3,5x,d10.3,5x,d10.3)
765: 105 loc=nodplc(loc)
766: go to 100
767: c
768: c voltage sources
769: c
770: 110 loc=locate(9)
771: 140 if (loc.eq.0) go to 150
772: locv=nodplc(loc+1)
773: val=value(locv+1)
774: iptrv=nodplc(loc+6)
775: sens=-value(lvn+iptrv)
776: sensn=val*sens/100.0d0
777: write (6,101) value(locv),val,sens,sensn
778: 145 loc=nodplc(loc)
779: go to 140
780: c
781: c current sources
782: c
783: 150 loc=locate(10)
784: 160 if (loc.eq.0) go to 170
785: locv=nodplc(loc+1)
786: node1=nodplc(loc+2)
787: node2=nodplc(loc+3)
788: val=value(locv+1)
789: sens=value(lvn+node1)-value(lvn+node2)
790: sensn=val*sens/100.0d0
791: write (6,101) value(locv),val,sens,sensn
792: 165 loc=nodplc(loc)
793: go to 160
794: c
795: c diodes
796: c
797: 170 loc=locate(11)
798: 180 if (loc.eq.0) go to 210
799: locv=nodplc(loc+1)
800: write (6,181) value(locv)
801: 181 format(1x,a8)
802: node1=nodplc(loc+2)
803: node2=nodplc(loc+3)
804: node3=nodplc(loc+4)
805: locm=nodplc(loc+5)
806: locm=nodplc(locm+1)
807: area=value(locv+1)
808: c
809: c series resistance (rs)
810: c
811: val=value(locm+2)*area
812: if (val.ne.0.0d0) go to 190
813: write (6,186) alsrs
814: 186 format(10x,a8,5x,2h0.,13x,2h0.,13x,2h0.)
815: go to 200
816: 190 val=1.0d0/val
817: sens=-(value(lvnim1+node1)-value(lvnim1+node3))*
818: 1 (value(lvn +node1)-value(lvn +node3))/(val*val)
819: sensn=val*sens/100.0d0
820: write (6,101) alsrs,val,sens,sensn
821: c
822: c intrinsic parameters
823: c
824: 200 csat=value(locm+1)*area
825: xn=value(locm+3)
826: vbe=value(lvnim1+node3)-value(lvnim1+node2)
827: vte=xn*vt
828: evbe=dexp(vbe/vte)
829: vabe=value(lvn+node3)-value(lvn+node2)
830: c
831: c saturation current (is)
832: c
833: sens=vabe*(evbe-1.0d0)
834: sensn=csat*sens/100.0d0
835: write (6,101) alsis,csat,sens,sensn
836: c
837: c ideality factor (n)
838: c
839: sens=-vabe*(csat/xn)*(vbe/vte)*evbe
840: if (dabs(sens).lt.1.0d-50) sens=0.0d0
841: sensn=xn*sens/100.0d0
842: write (6,101) alsn,xn,sens,sensn
843: 205 loc=nodplc(loc)
844: go to 180
845: c
846: c bipolar junction transistors
847: c
848: 210 loc=locate(12)
849: 220 if (loc.eq.0) go to 1000
850: locv=nodplc(loc+1)
851: write (6,181) value(locv)
852: node1=nodplc(loc+2)
853: node2=nodplc(loc+3)
854: node3=nodplc(loc+4)
855: node4=nodplc(loc+5)
856: node5=nodplc(loc+6)
857: node6=nodplc(loc+7)
858: locm=nodplc(loc+8)
859: type=nodplc(locm+2)
860: locm=nodplc(locm+1)
861: loct=lx0+nodplc(loc+22)
862: area=value(locv+1)
863: c
864: c base resistance (rb)
865: c
866: val=value(loct+16)
867: if (val.ne.0.0d0) go to 230
868: write (6,186) alsrb
869: go to 240
870: 230 val=1.0d0/val
871: sens=-(value(lvnim1+node2)-value(lvnim1+node5))*
872: 1 (value(lvn +node2)-value(lvn +node5))/(val*val)
873: sensn=val*sens/100.0d0
874: write (6,101) alsrb,val,sens,sensn
875: c
876: c collector resistance (rc)
877: c
878: 240 val=value(locm+20)*area
879: if (val.ne.0.0d0) go to 250
880: write (6,186) alsrc
881: go to 260
882: 250 val=1.0d0/val
883: sens=-(value(lvnim1+node1)-value(lvnim1+node4))*
884: 1 (value(lvn +node1)-value(lvn +node4))/(val*val)
885: sensn=val*sens/100.0d0
886: write (6,101) alsrc,val,sens,sensn
887: c
888: c emitter resistance (re)
889: c
890: 260 val=value(locm+19)*area
891: if (val.ne.0.0d0) go to 270
892: write (6,186) alsre
893: go to 280
894: 270 val=1.0d0/val
895: sens=-(value(lvnim1+node3)-value(lvnim1+node6))*
896: 1 (value(lvn +node3)-value(lvn +node6))/(val*val)
897: sensn=val*sens/100.0d0
898: write (6,101) alsre,val,sens,sensn
899: c
900: c intrinsic parameters
901: c
902: 280 bf=value(locm+2)
903: br=value(locm+8)
904: csat=value(locm+1)*area
905: ova=value(locm+4)
906: ovb=value(locm+19)
907: oik=value(locm+5)/area
908: c2=value(locm+6)*area
909: xne=value(locm+7)
910: vte=xne*vt
911: oikr=value(locm+11)/area
912: c4=value(locm+12)*area
913: xnc=value(locm+13)
914: vtc=xnc*vt
915: vbe=type*(value(lvnim1+node5)-value(lvnim1+node6))
916: vbc=type*(value(lvnim1+node5)-value(lvnim1+node4))
917: vabe=type*(value(lvn+node5)-value(lvn+node6))
918: vabc=type*(value(lvn+node5)-value(lvn+node4))
919: vace=vabe-vabc
920: if (vbe.le.-vt) go to 320
921: evbe=dexp(vbe/vt/value(locm+3))
922: cbe=csat*(evbe-1.0d0)
923: gbe=csat*evbe/vt/value(locm+3)
924: if (c2.ne.0.0d0) go to 310
925: cben=0.0d0
926: gben=0.0d0
927: go to 350
928: 310 evben=dexp(vbe/vte)
929: cben=c2 *(evben-1.0d0)
930: gben=c2 *evben/vte
931: go to 350
932: 320 gbe=-csat/vbe
933: cbe=gbe*vbe
934: gben=-c2/vbe
935: cben=gben*vbe
936: 350 if (vbc.le.-vt) go to 370
937: evbc=dexp(vbc/vt/value(locm+9))
938: cbc=csat*(evbc-1.0d0)
939: gbc=csat*evbc/vt/value(locm+9)
940: if (c4.ne.0.0d0) go to 360
941: cbcn=0.0d0
942: gbcn=0.0d0
943: go to 400
944: 360 evbcn=dexp(vbc/vtc)
945: cbcn=c4 *(evbcn-1.0d0)
946: gbcn=c4 *evbcn/vtc
947: go to 400
948: 370 gbc=-csat/vbc
949: cbc=gbc*vbc
950: gbcn=-c4/vbc
951: cbcn=gbcn*vbc
952: 400 q1=1.0d0/(1.0d0-ova*vbc-ovb*vbe)
953: q2=oik*cbe+oikr*cbc
954: sqarg=dsqrt(1.0d0+4.0d0*q2)
955: qb=q1*(1.0d0+sqarg)/2.0d0
956: dqb=(cbe-cbc)/(qb*qb)
957: sqarg=dsqrt(1.0d0+4.0d0*q2)
958: dq1=dqb*(1.0d0+sqarg)+(q1*q1)/2.0d0
959: dq2=q1*dqb/sqarg
960: c
961: c compute sensitivities
962: c
963: c... bf
964: sens=-vabe*cbe/bf/bf
965: sensn=bf*sens/100.0d0
966: write (6,101) alsbf,bf,sens,sensn
967: c... jle
968: if (c2.ne.0.0d0) go to 430
969: write (6,186) alsc2
970: go to 440
971: 430 sens=vabe*cben/c2
972: sensn=c2*sens/100.0d0
973: write (6,101) alsc2,c2,sens,sensn
974: c... br
975: 440 sens=-vabc*cbc/br/br
976: sensn=br*sens/100.0d0
977: write (6,101) alsbr,br,sens,sensn
978: c... jlc
979: if (c4.ne.0.0d0) go to 450
980: write (6,186) alsc4
981: go to 460
982: 450 sens=vabc*cbcn/c4
983: sensn=c4*sens/100.0d0
984: write (6,101) alsc4,c4,sens,sensn
985: c... is
986: 460 sens=(vabe*(cbe/bf)+vabc*(cbc/br)
987: 1 +vace*(dqb*qb-dq2*q2))/csat
988: sensn=csat*sens/100.0d0
989: write (6,101) alsjs,csat,sens,sensn
990: c... ne
991: sens=-vabe*gben*vbe/xne
992: sensn=xne*sens/100.0d0
993: write (6,101) alsne,xne,sens,sensn
994: c... nc
995: sens=-vabc*gbcn*vbc/xnc
996: sensn=xnc*sens/100.0d0
997: write (6,101) alsnc,xnc,sens,sensn
998: c... ik
999: if (oik.ne.0.0d0) go to 470
1000: write (6,186) alsik
1001: go to 480
1002: 470 val=1.0d0/oik
1003: sens=vace*dq2*cbe/(val*val)
1004: sensn=val*sens/100.0d0
1005: write (6,101) alsik,val,sens,sensn
1006: c... ikr
1007: 480 if (oikr.ne.0.0d0) go to 490
1008: write (6,186) alsikr
1009: go to 500
1010: 490 val=1.0d0/oikr
1011: sens=vace*dq2*cbc/(val*val)
1012: sensn=val*sens/100.0d0
1013: write (6,101) alsikr,val,sens,sensn
1014: c... va
1015: 500 if (ova.ne.0.0d0) go to 510
1016: write (6,186) alsva
1017: go to 520
1018: 510 va=1.0d0/ova
1019: sens=vace*dq1*vbc/(va*va)
1020: sensn=va*sens/100.0d0
1021: write (6,101) alsva,va,sens,sensn
1022: c... vb
1023: 520 if (ovb.ne.0.0d0) go to 530
1024: write (6,186) alsvb
1025: go to 540
1026: 530 vb=1.0d0/ovb
1027: sens=vace*dq1*vbe/(vb*vb)
1028: sensn=vb*sens/100.0d0
1029: write (6,101) alsvb,vb,sens,sensn
1030: c
1031: c
1032: 540 loc=nodplc(loc)
1033: go to 220
1034: c
1035: c finished
1036: c
1037: 1000 continue
1038: return
1039: end
1040: subroutine asol
1041: implicit double precision (a-h,o-z)
1042: c
1043: c this routine evaluates the adjoint circuit response by doing a
1044: c forward/backward substitution on the transpose of the coefficient
1045: c matrix.
1046: c
1047: common /tabinf/ ielmnt,isbckt,nsbckt,iunsat,nunsat,itemps,numtem,
1048: 1 isens,nsens,ifour,nfour,ifield,icode,idelim,icolum,insize,
1049: 2 junode,lsbkpt,numbkp,iorder,jmnode,iur,iuc,ilc,ilr,numoff,isr,
1050: 3 nmoffc,iseq,iseq1,neqn,nodevs,ndiag,iswap,iequa,macins,lvnim1,
1051: 4 lx0,lvn,lynl,lyu,lyl,lx1,lx2,lx3,lx4,lx5,lx6,lx7,ld0,ld1,ltd,
1052: 5 imynl,imvn,lcvn,loutpt,nsnod,nsmat,nsval,icnod,icmat,icval
1053: common /cirdat/ locate(50),jelcnt(50),nunods,ncnods,numnod,nstop,
1054: 1 nut,nlt,nxtrm,ndist,ntlin,ibr,numvs
1055: common /blank/ value(1000)
1056: integer nodplc(64)
1057: complex*16 cvalue(32)
1058: equivalence (value(1),nodplc(1),cvalue(1))
1059: c
1060: c forward substitution
1061: c
1062: do 20 i=2,nstop
1063: io=nodplc(iorder+i)
1064: value(lvn+io)=value(lvn+io)/value(lynl+io)
1065: jstart=nodplc(iur+i)
1066: jstop=nodplc(iur+i+1)-1
1067: if (jstart.gt.jstop) go to 20
1068: if (value(lvn+io).eq.0.0d0) go to 20
1069: do 10 j=jstart,jstop
1070: jo=nodplc(iuc+j)
1071: jo=nodplc(iorder+jo)
1072: value(lvn+jo)=value(lvn+jo)-value(lyu+j)*value(lvn+io)
1073: 10 continue
1074: 20 continue
1075: c
1076: c backward substitution
1077: c
1078: k=nstop+1
1079: do 40 i=2,nstop
1080: k=k-1
1081: io=nodplc(iorder+k)
1082: jstart=nodplc(ilc+k)
1083: jstop=nodplc(ilc+k+1)-1
1084: if (jstart.gt.jstop) go to 40
1085: do 30 j=jstart,jstop
1086: jo=nodplc(ilr+j)
1087: jo=nodplc(iorder+jo)
1088: value(lvn+io)=value(lvn+io)-value(lyl+j)*value(lvn+jo)
1089: 30 continue
1090: 40 continue
1091: c
1092: c reorder right-hand side
1093: c
1094: do 50 i=2,nstop
1095: j=nodplc(iswap+i)
1096: value(ndiag+i)=value(lvn+j)
1097: 50 continue
1098: call copy8(value(ndiag+1),value(lvn+1),nstop)
1099: c
1100: c finished
1101: c
1102: return
1103: end
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.