C ******************************************************************* C rtc_init.f nstx power supply real time control program C C Initialization of COMMON block C C REVISONS: C -------- C 07/08/05 c.neumeyer Changed switchover of p13 alpha to 0.9999 of vmax C 03/02/05 c.neumeyer Add axial force and tf joint calc inputs (v6) C 04/28/04 t.gibney Eliminate CHI-specific code. C 11/11/03 t.gibney Change simulation estimates for dtinput and dtcalc. C Change default dt for simulation to 200usec. C 07/28/03 t.gibney Change from double-precision to single-precision. C 07/10/01 t.gibney Skip 2nd half of nfault(*,2) when reading o.d. C 10/25/00 t gibney use parameters for NBRANCHES, NCOILS, etc. C 05/25/00 t.gibney Remove structure refs, compatible with gnu fortran. C 02/22/00 t.gibney remove init of disused "raw2phys" in pcs_common. C 01/10/00 t.gibney Support for c.d version 01: chi arc checking. C 10/25/99 t.gibney Call rtcc_init (rtc configuration init), shared mem. C 07/28/99 t.gibney Pass icr as argument to subroutine build(). C 05/14/99 t.gibney Move some cycle-specific inits (e.g., npulse, C npcon) to newCycle_init subroutine. C When running with pcs, override certain c.d values. C 01/12/99 t.gibney Change current directory before opening input files. C - Sky ftn OPEN has max filename length of 39 chars! C 01/06/99 t.gibney Change usage of nselectreal and nselectint. C 12/29/98 t.gibney Initialize nfs to 2**(nbit-1), not 2**nbit-1. C Do *not* open d.o until writing it after shot. C Close all files after reading them. C 12/02/98 t.gibney allow RTC_INPUT: user-specified path to input. C 12/01/98 c neumeyer added fixes during PTP-ECS-034 C 11/17/98 c neumeyer added modifications for collect data C 11/06/98 c neumeyer added TF & CHI link requirement C 10/21/98 c neumeyer various post FDR fixes C 10/04/98 c neumeyer various fixes in preparation for FDR C 09/16/98 b davis bug fix in build call. Calls to dtime. C 09/09/98 c neumeyer C ******************************************************************* SUBROUTINE rtc_init IMPLICIT NONE !Include and common: include 'rtcdef.inc' include 'pcs_common.inc' include 'iboc_override.inc' include 'nstx_savedata.inc' !Arguments: ! *none* ! Functions: integer nonblankLen integer system ! Local variables: real*4 ecl,fsmeas,gamma,icr real*4 lc real*4 miniloc,mini2tmax,miniclp,maxicln,sumchk integer nlmeas,npmeas integer nin,nout integer nreal,nrealmax,nint,nintmax,nsave,nsavemax DIMENSION ecl(NCOILS),fsmeas(NBRANCHES,2),gamma(NCOILS) DIMENSION icr(NCOILS),lc(NCOILS) DIMENSION nlmeas(NBRANCHES,2),npmeas(NBRANCHES,2) real*4 dtcalc,dtinput,dtoutput,lclr,rclr,vac integer n,nn,n1,n2,nbit,nfs integer HT !Horizontal tab char integer iostat integer iversion_cd !version number of c.d file character*80 dmy character*60 path !path to input files: RTC_INPUT character*60 userDir !user's working directory integer nfaultOffset !because o.d file lists only first half ! of nfault(NCOILS,2). integer npath !length of path data HT /9/ !Horizontal Tab 900 FORMAT(/' *ERR* no input files in directory ''',A,'''', 1 /' --> Check RTC_INPUT environment variable'/) 901 FORMAT(A16) 902 FORMAT(/'rtc_init: Loading control data from ''',A,'''') 903 FORMAT(A,' inputDir=',A) 904 FORMAT(/'rtc_init: *EXIT* c.d file is version',I2,/) 905 FORMAT(/'rtc_init: *EXIT* old c.d file was version',I2,/ > ' It has been rewritten, and the old c.d renamed',/ > ' --> re run the program',/) 906 FORMAT(A7,I2) !"Version n" line (several files) 15-Nov'04 999 FORMAT(A) C C open files /////////////////////////////////////////// C C c.d contains control data C t.d contains timing and waveform data C s.d contains simulation data C o.d file is collect data selection file C write(6,911) 911 format(/'==> Inside rtc_init') !Initialize circuitnames: informational only !-------------------------------------------------- circuitname(1)='OH' circuitname(2)='PF1AU' circuitname(3)='PF1AL' circuitname(4)='PF1B' circuitname(5)='PF2U' circuitname(6)='PF2L' circuitname(7)='PF3U' circuitname(8)='PF3L' circuitname(9)='PF5' circuitname(10)='CHI' circuitname(11)='TF' circuitname(12)='PF4' call getenv('RTC_TEST',dmy) if (dmy(1:1).eq.'1') then testMode = .TRUE. else circuitname(13)='RWM1' circuitname(14)='RWM2' circuitname(15)='RWM3' testMode = .FALSE. endif call getenv('RTC_INPUT',path) !check env variable npath = index(path,' ') - 1 !idx of last non-blank char if (npath.eq.0) then path = './' npath = 2 else if (path(npath:npath).ne.'/') then npath = npath+1 path(npath:npath) = '/' endif endif if (system('test -f '//path(1:npath)//'c.d') .ne. 0) then write(6,900) path(1:npath) call pscExit(0,'c.d file not found') endif write(6,902) path(1:npath) write(6,*) call getwd(userDir) OPEN(10,FILE=path(1:npath)//'c.d',STATUS='old') OPEN(20,FILE=path(1:npath)//'t.d',STATUS='old') OPEN(30,FILE=path(1:npath)//'s.d',STATUS='old') OPEN(40,FILE=path(1:npath)//'o.d',STATUS='old') C C read control data /////////////////////////////////////////// C C rescale all data as required to fundamental engineering units C READ(10,901)dmy READ(dmy,'(5x,I2)') n if (n .ne. 6) then !Check current version num write(6,904) n call pscExit(0,'*ERR* bad c.d version') endif READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)nmode READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)nsubmode DO 90, n = 1,NCOILS READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)nlock(n) if(nlock(n).EQ.1)write(6,*)circuitname(n),' OUT' if(nlock(n).EQ.0)write(6,*)circuitname(n),' IN' 90 CONTINUE DO 100, n = 1,NCOILS READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)npc(n) 100 CONTINUE DO 103, n = 1,NCOILS DO 102, n1 = 1,2 DO 101, n2 = 1,2 READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)gain(n,n2,n1) 101 CONTINUE 102 CONTINUE 103 CONTINUE DO 105, n = 1,NBRANCHES DO 104, n1 = 1,2 READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)nmeas(n,n1) 104 CONTINUE 105 CONTINUE DO 111, n = 1,NBRANCHES DO 110, n1 = 1,2 READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)npmeas(n,n1) 110 CONTINUE 111 CONTINUE DO 113, n = 1,NBRANCHES DO 112, n1 = 1,2 READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)nlmeas(n,n1) 112 CONTINUE 113 CONTINUE DO 115, n = 1,NBRANCHES DO 114, n1 = 1,2 READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)fsmeas(n,n1) fsmeas(n,n1) = fsmeas(n,n1)*1.e3 114 CONTINUE 115 CONTINUE DO 116, n = 1,NBRANCHES READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)ncb(n) 116 CONTINUE DO 117, n = 1,NBRANCHES READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)kpb(n) 117 CONTINUE DO 118, n = 1,NBRANCHES READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)dimax(n) dimax(n) = dimax(n)*1.e3 118 CONTINUE DO 120, n = 1,NCOILS READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)nplink(n) 120 CONTINUE DO 121, n = 1,NCOILS READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)np(n) 121 CONTINUE DO 122, n = 1,NCOILS READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)ns(n) 122 CONTINUE DO 123, n = 1,NCOILS READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)iclp(n) iclp(n) = iclp(n)*1.e3 123 CONTINUE DO 126, n = 1,NCOILS READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)icln(n) icln(n) = icln(n)*1.e3 126 CONTINUE READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)vac vac = vac*1.e3 READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)fac DO 130, n = 1,NCOILS READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)gamma(n) 130 CONTINUE DO 131, n = 1,NCOILS READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)lc(n) lc(n) = lc(n)*1.e-6 131 CONTINUE DO 132, n = 1,NCOILS READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)vdotmax(n) 132 CONTINUE DO 133, n = 1,NCOILS READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)nap(n) 133 CONTINUE DO 134, n = 1,NCOILS READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)napblock(n) 134 CONTINUE DO 135, n = 1,NCOILS READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)ithresh(n) ithresh(n) = ithresh(n)*1.e3 135 CONTINUE DO 141, n = 1,NBRANCHES DO 140, n1 = 1,2 READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)iboc(n,n1) iboc(n,n1) = iboc(n,n1)*1.e3 140 CONTINUE 141 CONTINUE if (niboc .ne. 0) then !Optional overrided of IBOC DO 143 n = 1,NBRANCHES iboc(n,1) = ibocStbyOverride * 1.e3 iboc(n,2) = ibocPulseOverride * 1.e3 143 CONTINUE write(6,912) 912 format(/'rtc_init: *NOTE* iboc is being overridden'/) endif DO 145, n = 1,NCOILS READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)iloc(n) iloc(n) = iloc(n)*1.e3 145 CONTINUE DO 150, n = 1,NCOILS READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)i2tmax(n) i2tmax(n) = i2tmax(n)*1.e6 150 CONTINUE DO 156, n = 1,NPOLES DO 155, n1 = 1,2 READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)npsum(n,n1) 155 CONTINUE 156 CONTINUE DO 160, n = 1,NPOLES READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)r20pole(n) r20pole(n) = r20pole(n)/1.e6 160 CONTINUE DO 161, n = 1,NPOLES READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)trpole(n) 161 CONTINUE DO 162, n = 1,NPOLES READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)tcpole(n) 162 CONTINUE DO 163, n = 1,NPOLES READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)nppole(n) 163 CONTINUE READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)t0pole READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)tmpole READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)alaunch alaunch = alaunch*1.e6 READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)itfmin itfmin = itfmin*1.e3 READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)ichi0 ichi0 = ichi0*1.e3 DO 165, n = 1,NCOILS READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)scd(n) scd(n) = scd(n)*1.e6 165 CONTINUE DO 166, n = 1,NCOILS READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)trcoil(n) trcoil(n) = trcoil(n)/1.e3 166 CONTINUE DO 167, n = 1,NCOILS READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)tccoil(n) tccoil(n) = tccoil(n)*1.e3 167 CONTINUE READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)coeff DO 168, n = 1,NCOILS READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)t0coil(n) 168 CONTINUE DO 169, n = 1,NCOILS READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)tmcoil(n) 169 CONTINUE DO 170, n = 1,NCOILS READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)rpss(n) rpss(n) = rpss(n)/1.e3 170 CONTINUE READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)finv DO 171, n = 1,NCOILS READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)ecr(n) ecr(n) = ecr(n)/1.e3 171 CONTINUE DO 172, n = 1,NCOILS READ(10,901)dmy READ(dmy(1:index(dmy,char(HT))-1),*)ecl(n) ecl(n) = ecl(n)/1.e6 172 CONTINUE !Up to here, we have been reading "version 0" ! of the c.d file. ! Set the version number variable, and then read ! variables which have been tacked on for that ! version of c.d. ! Repeat as required for each new version ... iversion_cd = 1 !variables from c.d version 1 READ(10,901,iostat=iostat,end=178) dmy READ(dmy(1:index(dmy,char(HT))-1),*) zchimin READ(10,901) dmy READ(dmy(1:index(dmy,char(HT))-1),*) nchidelay nchiarc = 0 !added 10-Jan-00 TRG READ(10,901) dmy READ(dmy(1:index(dmy,char(HT))-1),*) vp13spa C C Additions for axial force fz and tf joint C DO 175, n=1,NCOILS_V0 READ(10,901) dmy READ(dmy(1:index(dmy,char(HT))-1),*)fzmax(n) fzmax(n)=fzmax(n)*1e3 175 CONTINUE READ(10,901) dmy READ(dmy(1:index(dmy,char(HT))-1),*) ptfjmax ptfjmax=ptfjmax*1e3 READ(10,901) dmy READ(dmy(1:index(dmy,char(HT))-1),*) ftfjbmax ftfjbmax=ftfjbmax*1e3 178 close(10) if (iostat .NE. 0) then call rewrite_cd(iversion_cd,path(1:npath)) write(6,905) iversion_cd-1 call exit(0) endif C C read timing data /////////////////////////////////////////// C READ(20,906) dmy,n !"Version n" line added 15-Nov'04 if (n.NE.6) then call pscExit(0,'t.d wrong version') endif READ(20,*)tsop READ(20,*)tpcstart READ(20,*)tpcend DO 180, n = 1,NCOILS READ(20,*)tireset(n) 180 CONTINUE DO 181, n = 1,NCOILS READ(20,*)nbpoint(n) 181 CONTINUE DO 185, n = 1,NCOILS DO 184, n1 = 1,nbpoint(n) READ(20,*)dat(n,n1,1),dat(n,n1,2),dat(n,n1,3) dat(n,n1,2) = dat(n,n1,2)*1.e3 184 CONTINUE 185 CONTINUE close(20) IF (nsubmode.EQ.1) THEN C C read simulation data //////////////////////////////////////// C READ(30,*)nchisdso READ(30,*)nchisdsc READ(30,*)nhcstf READ(30,*)nhcsoh READ(30,*)nhcspfchi READ(30,*)fmg READ(30,*)vmg READ(30,*)npcready ELSE ! added by R. Hatcher to give a npcready npcready = 1 ! in non-simulation mode 12/29/98 ENDIF close(30) C C read cntmtx.d C call chdir(path(1:npath)) !cd to "path" directory call getwd(inputDir) OPEN(17,FILE='cntmtx.d',STATUS='old') call chdir(userDir) !cd to original userDir READ(17,906) dmy,n !"Version n" line added 15-Nov'04 if (n.NE.6) then call pscExit(0,'cntmtx.d wrong version') endif DO n = 1,NCOILS READ(17,'(1p15e12.4)')(mc(n,n1),n1=1,NCOILS) ENDDO C WRITE(6,*)'cntmtx inductance read completed' READ(17,'(1p15e12.4)')(icr(n),n=1,NCOILS) DO n=1,NCOILS r20c(n)=icr(n) r20coil(n)=icr(n) END DO C C read collect data selections/////////////////////////////////// C C OPEN(10,FILE='d.o.trailer',STATUS='unknown') write(10,903) 'Channel names:',inputDir(1:nonblanklen(inputDir)) nrealmax=NMAX_VARIABLES nintmax=NMAX_VARIABLES/2 nsavemax=NMAX_VARIABLES nreal=0 nint=0 nsave=0 READ(40,906) dmy,n !"Version n" line added 15-Nov'04 if (n.NE.6) then call pscExit(0,'o.d wrong version') endif C C skip next line of file (column headings) C READ(40,*) DO 190 n=1,NSAVE_REAL READ(40,999)dmy !!! READ(dmy(1:index(dmy,char(HT))-1),*)nselectreal(n) !!! if(nselectreal(n).eq.1)then !!! nreal=nreal+1 !!! nsave=nsave+1 !!! end if READ(dmy(1:index(dmy,char(HT))-1),*) nn if (nn .eq. 1) then nreal=nreal+1 nselectreal(nreal) = n nsave=nsave+1 write(10,999) dmy(1:nonblankLen(dmy)) endif 190 CONTINUE if (nreal.LT.NMAX_VARIABLES) nselectreal(nreal+1) = 0 nfaultOffset = 0 !offset to skip 2nd half of nfault(*,2) DO 195 n=1,NSAVE_INT-3 C C skip navailable1,2,3 at end of common block C READ(40,999,end=196)dmy !!! READ(dmy(1:index(dmy,char(HT))-1),*)nselectint(n) !!! if(nselectint(n).eq.1)then !!! nint=nint+1 !!! nsave=nsave+1 !!! end if READ(dmy(1:index(dmy,char(HT))-1),*) nn if (n.EQ.7*NCOILS+1) then !nchiblock immediately follows nfault(*,2) nfaultOffset = NCOILS endif if (nn .eq. 1) then nint=nint+1 nselectint(nint) = n+nfaultOffset !adjust for nfault(*,2) nsave=nsave+1 write(10,999) dmy(1:nonblankLen(dmy)) endif 195 CONTINUE goto 197 196 CONTINUE !here if EOF on o.d file ... write(6,998) 998 format(/'--> *NOTE* o.d file is short!'/) if (n.EQ.92 .AND. nsave.LT.NMAX_VARIABLES) then nint = nint+1 nselectint(nint) = n+nfaultOffset nsave = nsave+1 write(10,997) char(HT),char(HT),char(HT),char(HT) 997 format('1',A,'nclampm',A,'(all)',A,'1.0',A,'n.a.') endif 197 CONTINUE if (nint.LT.NMAX_VARIABLES) nselectint(nint+1) = 0 close(10) close(40) C C if(nreal.gt.nrealmax.or.nint.gt.nintmax.or.nsave.gt.nsavemax)then write(6,*)'Error: Collect data request exceeds limit', & ' nreal,nint,nsave=',nreal,nint,nsave stop end if C C initialize ////////////////////////////////////////////////// C nckt = NCOILS !!! if (nsubmode .EQ. 0) then !!! dt = 0.200e-3 !!! else !!! dt = 1.0e-3 !!! endif dt = 0.200e-3 CXXXC CXXXC if simulation submode, issue sop event CXXXC CXXX if(nsubmode.eq.1)nsopeop=1 C C input, calculation, and output latencies C dtinput=150.e-6 dtcalc=350.e-6 dtoutput=100.e-6 dtinput=7.e-6 !new value 11-Nov-2003 TRG dtcalc=76.e-6 !new value 11-Nov-2003 TRG dtoutput=55.e-6 !new value 11-Nov-2003 TRG C C rclr and lclr are res and ind of current limiting reactors C rclr = 1.e-3 lclr = 265.e-6 C C nfs is bit count at measurement full scale value C nbit = 16 nfs = 2**(nbit-1) C C aomeas is simple moving average offset of measurement in amps C pi = 3.141592654 C C adjust number of series sections in twin circuits C DO 205, n = 1,NCOILS IF (nlock(n).EQ.0) GOTO 205 C C circuit is locked out C C if twin circuit is active, adjust number of ps sections C C PF1aU and PF1aL C IF (n.EQ.2.AND.nlock(3).NE.1) ns(3) = ns(2)+ns(3) IF (n.EQ.3.AND.nlock(2).NE.1) ns(2) = ns(2)+ns(3) C C PF2U and PF2L C IF (n.EQ.5.AND.nlock(6).NE.1) ns(6) = ns(5)+ns(6) IF (n.EQ.6.AND.nlock(5).NE.1) ns(5) = ns(5)+ns(6) C C PF3U and PF3L C IF (n.EQ.7.AND.nlock(8).NE.1) ns(7) = ns(7)+ns(8) IF (n.EQ.8.AND.nlock(7).NE.1) ns(8) = ns(7)+ns(8) 205 CONTINUE C C if twin circuits are operated as one then set trip/clamp C settings to the minimum of the two twin values C IF (nlock(2)+nlock(3).eq.1) THEN C C PF1a twins are operating as one circuit C miniloc=min(iloc(2),iloc(3)) iloc(2)=miniloc iloc(3)=miniloc miniloc=min(iloc(2),iloc(3)) iloc(2)=miniloc iloc(3)=miniloc mini2tmax=min(i2tmax(2),i2tmax(3)) i2tmax(2)=mini2tmax i2tmax(3)=mini2tmax miniclp=min(iclp(2),iclp(3)) iclp(2)=miniclp iclp(3)=miniclp maxicln=max(icln(2),icln(3)) icln(2)=maxicln icln(3)=maxicln END IF IF (nlock(5)+nlock(6).eq.1) THEN C C PF2 twins are operating as one circuit C miniloc=min(iloc(5),iloc(6)) iloc(5)=miniloc iloc(6)=miniloc miniloc=min(iloc(5),iloc(6)) iloc(5)=miniloc iloc(6)=miniloc mini2tmax=min(i2tmax(5),i2tmax(6)) i2tmax(5)=mini2tmax i2tmax(6)=mini2tmax miniclp=min(iclp(5),iclp(6)) iclp(5)=miniclp iclp(6)=miniclp maxicln=max(icln(5),icln(6)) icln(5)=maxicln icln(6)=maxicln END IF IF (nlock(7)+nlock(8).eq.1) THEN C C PF3 twins are operating as one circuit C miniloc=min(iloc(7),iloc(8)) iloc(7)=miniloc iloc(8)=miniloc miniloc=min(iloc(7),iloc(8)) iloc(7)=miniloc iloc(8)=miniloc mini2tmax=min(i2tmax(7),i2tmax(8)) i2tmax(7)=mini2tmax i2tmax(8)=mini2tmax miniclp=min(iclp(7),iclp(8)) iclp(7)=miniclp iclp(8)=miniclp maxicln=max(icln(7),icln(8)) icln(7)=maxicln icln(8)=maxicln END IF DO 210, n = 1,NCOILS ndat(n) = 1 vprtc(n) = 0. nclamp(n) = 0. nfault(n,1) = 0 nfault(n,2)=0 iload(n) = 0. i2t(n) = 0. if (tcoil(n).EQ.0.0) tcoil(n) = t0coil(n) disd(n) = scd(n)*dt !recalc'd in rtc_indat() each time if (n .LE. NCOILS_V0) then !init for Original coil set (non-RWM/SPA) vmax(n) = ns(n)*3.*sqrt(2.)*vac*0.75/pi/13.8 dvmax(n) = vdotmax(n)*2.*pi*fac*dt !recalc'd inside rtc_indat() damax = 2.*pi*fac*dt IF (nap(n).EQ.1) THEN if(napblock(n).eq.0)then napon(n)=1 else napon(n)=0 end if astdby(n) = 90.*pi/180. ELSE napon(n) = 0 astdby(n) = 165.*pi/180. ENDIF cosgamma(n) = cos(gamma(n)*pi/180.) amaxmult(n) =4.*pi*fac*lc(n)*13.8e3/sqrt(2.)/vac/750.0/np(n) else !..else init for RWM/SPA coils C C temporary fudge to allow local closed loop current control at SPA C set full scale to 3333.333 volts (interpreted as amps by SPA) C vmax(n) = vp13spa !max SPA output voltage C dvmax(n) = vp13spa ! vmax(n) = 3333.333 dvmax(n) = vmax(n) astdby(n) = 0.0 !SPA stdby reference (ctrl) voltage endif !..end init for Original or RWM/SPA coils a1(n) = astdby(n) a2(n) = astdby(n) a1old(n) = a1(n) a2old(n) = a2(n) ncvt1(n) = 0 ncvt2(n) = 0 nbb(n) = 0 210 CONTINUE C C set up spa dc source power supply p13 C if(nlock(IDX_RWM1).eq.0 &.or.nlock(IDX_RWM2).eq.0 &.or.nlock(IDX_RWM3).eq.0)then C C at least one spa active, set up p13 C nlockp13=0 if(vp13spa.gt.sqrt(2.)*vac*750./13.8e3)then write(6,*)'ERROR: VP13SPA NOT REALIZABLE' stop end if ap13=120.*pi/180.-asin(vp13spa/(sqrt(2.)*vac*750./13.8e3)) C C if asking for full p13 voltage then set to alp = 0 C if(vp13spa.ge.0.9999*sqrt(2.)*vac*750./13.8e3)ap13=0. else C C no spas active, p13 not required C nlockp13=1 ap13=165.*pi/180. end if C C initialize p13 control settings, store in rwm1 command set C nbb not used by spa, use nbb(IDX_RWM1) for p13 C a2 and ncvt2 not used by spa, use a2(IDX_RWM1) and ncvt2(IDX_RWM1) C for p13 C a2(IDX_RWM1) = 165.*pi/180. ncvt2(IDX_RWM1) = 0 nbb(IDX_RWM1)=0 C C adjust gains if bus links are reversed C DO 212, n = 1,NCOILS C C skip if links are normal C IF (nplink(n).GE.0) GOTO 212 C C links are reversed C C reverse polarity of all feedback gains C gain(n,1,1) = -1.*gain(n,1,1) gain(n,1,2) = -1.*gain(n,1,2) gain(n,2,1) = -1.*gain(n,2,1) gain(n,2,2) = -1.*gain(n,2,2) C 212 CONTINUE C C initialize current measurement multiplier & zero out offset C DO 214, n1 = 1,NBRANCHES DO 213, n2 = 1,2 aomeas(n1,n2)=0. kmeas(n1,n2) = fsmeas(n1,n2)/(nfs+1) & *npmeas(n1,n2) C C condition polarity of measurement if DCCT on PS side and C links are reversed C IF (nlmeas(n1,n2).EQ.-1) THEN C C DCCT is on the power supply side of the links C kmeas(n1,n2) = kmeas(n1,n2)*nplink(ncb(n1)) ENDIF IF (kmeas(n1,n2).EQ.0.) THEN WRITE(6,*)'Error: kmeas=0, Branch,DCCT#=',n1,n2 WRITE(6,*)'Check fsmeas,npmeas,nplink .ne. 0' STOP ENDIF 213 CONTINUE 214 CONTINUE npsrtcready = 1 DO 215, n = 1,NPOLES if (tpole(n).EQ.0.0) tpole(n) = t0pole 215 CONTINUE IF (nmode.EQ.1) THEN C C test mode, no plasma C C set all plasma control enables to zero in test mode C DO 220, n = 1,NCOILS npc(n) = 0 220 CONTINUE ENDIF C C establish time window for control C CALL window(nbpoint,dat,tsop,tmin,tmax) C C check if tpcstart and tpcend within tmin and tmax C IF(tpcstart.lt.tmin.or.tpcstart.gt.tmax)THEN WRITE(6,*)'Error: Time to start plasma control outside', & ' window of sop and breakpoint data waveforms' STOP ENDIF IF(tpcend.lt.tmin.or.tpcend.gt.tmax)THEN WRITE(6,*)'Error: Time to end plasma control outside', & ' window of sop and breakpoint data waveforms' STOP ENDIF C C If running in conjunction with PCS C -- override c.d plasmacontrol settings C if (npcsFlag .or. npcsTest) then if (npcsFlag) then tpcstart = tmin !Override c.d value tpcend = tmax !Override c.d value !do n=1,NCOILS ! npc(n) = 1 !Override c.d value !enddo if (nctFlag) then npc(IDX_CHI) = 0 !..CHI psrtc control endif endif !!! pcs_common init now inside rtcc_init: 25-May-2000 else !Check for npc() NE 0 in c.d ... do n=1,NCOILS if (npc(n).NE.0) npcsFlag=.TRUE. !for rtcc_init enddo endif C C set up simulation /////////////////////////////////////////// C DO 240, n = 1,NCOILS IF (nlock(n).EQ.1) THEN C C circuit is locked out C C locked out circuits are assumed to be open C C locked out twins are lumped with enabled twins: C - their inductances are lumped with enabled twins C - their former circuit assignment is treated as open C ntype(n) = 2 ELSE C C active circuits are voltage specified C ntype(n) = 0 ENDIF C C set power supply and external circuit impedances C rps(n) = float(ns(n))/float(np(n))*rpss(n) C C do not add CLR impedance if RWM C if (n .GT. NCOILS_V0) GO TO 240 C C add clr resistance to external circuit resistance C ecr(n) = ecr(n)+rclr/np(n) C C add clr inductance to external circuit inductance C ecl(n) = ecl(n)+lclr/np(n) 240 CONTINUE C C modify mutual inductance matrix for twin circuits operating as C common pairs C if(nlock(2).eq.0.and.nlock(3).eq.1)then C C PF1aU will control PF1aU & PF1aL C nin=2 nout=3 CALL lump(nin,nout,mc,r20c,ecr,ecl) end if if(nlock(2).eq.1.and.nlock(3).eq.0)then C C PF1aL will control PF1aU & PF1aL C nin=3 nout=2 CALL lump(nin,nout,mc,r20c,ecr,ecl) end if if(nlock(5).eq.0.and.nlock(6).eq.1)then C C PF2U will control PF2U & PF2L C nin=5 nout=6 CALL lump(nin,nout,mc,r20c,ecr,ecl) end if if(nlock(5).eq.1.and.nlock(6).eq.0)then C C PF2L will control PF2U & PF2L C nin=6 nout=5 CALL lump(nin,nout,mc,r20c,ecr,ecl) end if if(nlock(7).eq.0.and.nlock(8).eq.1)then C C PF3U will control PF3U & PF3L C nin=7 nout=8 CALL lump(nin,nout,mc,r20c,ecr,ecl) end if if(nlock(7).eq.1.and.nlock(8).eq.0)then C C PF3L will control PF3U & PF3L C nin=8 nout=7 CALL lump(nin,nout,mc,r20c,ecr,ecl) end if C C set up matricies for solution of l-r circuit equations C CALL build(ecl,icr) C C extract 20C coil resistances C CX DO 245, n = 1,NCOILS CX r20coil(n) = rc(n)-ecr(n) CX r20c(n) = r20coil(n) C C set coil resistance C C note: this resistance is used in both simulation and operation C submodes to simulate coil temperature C CX rcoil(n)=r20coil(n)*(1.+coeff*(t0coil(n)-20.)) CX245 CONTINUE C C ndtsim is number of sim dt's per control system dt C ndtsim = 10 dtsim = dt/ndtsim C C ndtin is number of sim dt's between time of last command update C and sampling of next set of input data C ndtin = int((dt-dtinput-dtcalc-dtoutput)*ndtsim/dt) IF(ndtin.gt.ndtsim)THEN write(6,*)'Error: ndtin>ndtsim' stop ENDIF C C simulation setup complete C !Initialize rtc config structure (rtcc, shared mem) !..also initialize gains in pcs_common. !-------------------------------------------------- call rtcc_init C C set up indices for force variable calculation array fv C number of force variables is NFV set up in rtcconfig.inc C nfvc(1)=IDX_OH nfvc(2)=IDX_PF1AU nfvc(3)=IDX_PF1AL nfvc(4)=IDX_PF1B nfvc(5)=IDX_PF2U nfvc(6)=IDX_PF2L nfvc(7)=IDX_PF3U nfvc(8)=IDX_PF3L nfvc(9)=IDX_PF5 nfvc(10)=IDX_PF4 nfvc(11)=IDX_TF nfvc(12)=IDX_TF nfvc(13)=IDX_TF nfvc(14)=IDX_TF C C Initialize Influence Matrix and TF Joint Coefficients C DATA (infmtx(1,nn),nn=1,NCOILS_V0)/0,6.3,-6.45,-71.6, &53.9,-54.035,27.36,-27.5,-0.14,0,0,-0.14/ DATA (infmtx(2,nn),nn=1,NCOILS_V0)/-6.4,0,-0.005,-0.02, &10.95,-0.094,0.015,-0.35,-2.361,0,0,-1.781/ DATA (infmtx(3,nn),nn=1,NCOILS_V0)/6.3,-0.082,0,-87.7, &0.006,-11.04,0.26,-0.104,2.188,0,0,1.608/ DATA (infmtx(4,nn),nn=1,NCOILS_V0)/72.4,0.193,89.5,0, &0.53,-20.45,1.524,9.5,11.085,0,0,8.895/ DATA (infmtx(5,nn),nn=1,NCOILS_V0)/-55.1,-11.42,-0.37,-0.635, &0,-1.91,-98.58,-6.6,-57.468,0,0,-47.88/ DATA (infmtx(6,nn),nn=1,NCOILS_V0)/55.1,0.39,11.43,21.46, &1.93,0,6.62,98.6,57.515,0,0,47.93/ DATA (infmtx(7,nn),nn=1,NCOILS_V0)/-28.4,-0.53,-0.85,-1.86, &98.8,-6.83,0,-25.634,-280.506,0,0,-260.8/ DATA (infmtx(8,nn),nn=1,NCOILS_V0)/28.3,0.84,0.52,-8.9, &6.82,-98.87,25.624,0,280.49,0,0,260.775/ DATA (infmtx(9,nn),nn=1,NCOILS_V0)/6,0.4,-1.96,-7.9, &15.327,-41.58,65.146,-214.25,200.106,0,0,-431.47/ DATA (infmtx(10,nn),nn=1,NCOILS_V0)/7.4,1.01,-0.79,-6.285, & 9.27,-38.7,35.8,-228.18,633.72,0,0,52.3/ DATA (infmtx(11,nn),nn=1,NCOILS_V0)/0.468,0.001,0.044,0.378, &0.011,1.175,0.045,0.885,0.418,0,0,0.314/ DATA (infmtx(12,nn),nn=1,NCOILS_V0)/5.012,0.008,0.528,4.302, &0.151,15.62,0.613,12.21,5.745,0,0,4.318/ DATA (infmtx(13,nn),nn=1,NCOILS_V0)/-154.25,-0.015,-6.22,-8.596, &-0.207,-0.769,-0.708,-1.416,-3.517,0,0,-2.514/ DATA (infmtx(14,nn),nn=1,NCOILS_V0)/1.349,0.000,0.000,0.000, &0.000,0.000,0.000,0.000,0.000,0,0,0.000/ sumchk=0. DO 950 n=1,NFV C write(6,*)(infmtx(n,nn),nn=1,NCOILS_V0) DO 950 nn=1,NCOILS_V0 sumchk=sumchk+infmtx(n,nn) 950 CONTINUE C write(6,*)sumchk DATA kip/5120,0.307,40000,5.86,-1.17E-04/ DATA koop/5120,1.970,40000,5.86,-7.50E-04/ DATA kb/3.848,6.375,-38.781,-5.18,-9.0,10.125,0.412,3.0E-07/ return end ********************************************************************* SUBROUTINE window(nbpoint,dat,tsop,tmin,tmax) ********************************************************************* C C determine overall time window C fill in breakpoints as required to cover entire time window C IMPLICIT NONE include 'rtcconfig.inc' C Arguments: C---------- integer nbpoint real*4 dat,tsop,tmin,tmax DIMENSION nbpoint(NCOILS),dat(NCOILS,101,3) C Local variables: C---------------- integer n,n1,n2 real*4 tdat tmin = 0. tmax = 0. DO 100, n = 1,NCOILS DO 90, n1 = 1,nbpoint(n) tdat = dat(n,n1,1) IF (tdat.LT.tmin) tmin = tdat IF (tdat.GT.tmax) tmax = tdat 90 CONTINUE 100 CONTINUE IF (tsop.GT.tmin) THEN C C tsop not early enough C WRITE(6,*) & 'Error: Breakpoint Data SpecifiedPrior to SOP:Tsop,Tmin,T' & //'max=',tsop,tmin,tmax STOP ENDIF IF (tmin.lt.-0.75) THEN C C tmin starting earlier than allowed by acq C WRITE(6,*) & 'Error: Breakpoint Data Specified Prior to -0.75 sec' STOP ENDIF C C fill in breakpoints prior to those specified C DO 200, n = 1,NCOILS IF (dat(n,1,1).EQ.tmin) THEN C C skip this circuit, data already exists at start of window C GOTO 200 ENDIF IF (dat(n,1,2).NE.0.) THEN C C non-zero data at start C WRITE(6,*) & 'Warning: Non-zero value for 1st breakpoint,reference w' & //'ill ramp to 1st breakpoint from zero' ENDIF C C shift data points up by one place to make room for new 1st point C DO 150, n1 = nbpoint(n),1,-1 DO 145, n2 = 1,3 dat(n,n1+1,n2) = dat(n,n1,n2) 145 CONTINUE 150 CONTINUE nbpoint(n) = nbpoint(n)+1 C C add new breakpoint at start C time value = tmin C data value = 0 C psmode = same as 1st original breakpoint C dat(n,1,1) = tmin dat(n,1,2) = 0. dat(n,1,3) = dat(n,2,3) 200 CONTINUE C C add breakpoints following those specified C DO 300, n = 1,NCOILS IF (dat(n,nbpoint(n),1).EQ.tmax) THEN C C skip this circuit, data already exists at end of window C GOTO 300 ENDIF IF (dat(n,nbpoint(n),2).NE.0.) THEN C C non-zero data at last original breakpoint C WRITE(6,*) & 'Warning: Non-zero value for last breakpoint,reference ' & //'will ramp from last breakpoint to zero' ENDIF nbpoint(n) = nbpoint(n)+1 C C add new breakpoint at end C time value = tmax C data value = 0 C psmode = same as last breakpoint C dat(n,nbpoint(n),1) = tmax dat(n,nbpoint(n),2) = 0. dat(n,nbpoint(n),3) = dat(n,nbpoint(n)-1,3) 300 CONTINUE write(6,*) WRITE(6,*)'tsop=',tsop WRITE(6,*)'tmin=',tmin WRITE(6,*)'tmax=',tmax write(6,*) WRITE(6,*) '==> end of rtc_init' write(6,*) RETURN END C ************************************************************ SUBROUTINE build(ecl,icr) C ************************************************************ C C This version of BUILD is a modification of the one used with C the LRSIM codes. The read of cntmtx.d is placed outside of the C subroutine so that it can be contracted via the LUMP subroutine C C read icc data/cull out open circuits/invert inductance matrix C C nckt is total number of circuits C C ntype gives ckt types C C if ntype = 0, specification is v C if ntype = 1, specification is i C if ntype = 2, circuit is open C if ntype = 3, circuit is shorted C if ntype = 4, voltage feedback control C if ntype = 5, current feedback control C C nsim is number of ckts after removing open ckts C C mc and rc are total ind and res matricies inc'l all ckts C C m and r are ind and res matricies after culling open ckts C C rc is sum of internal ckt res (icr) and external ckt res (ecr) C mc is sum of internal ckt ind and external ckt ind (ecl) C C x is inverse of m matrix C C ncs contains total ckt number as function of sim circuit number C nsc contains sim circuit number as function of total circuit number C IMPLICIT NONE include 'rtcdef.inc' C Arguments: real*4 ecl real*4 icr DIMENSION ecl(NCOILS) DIMENSION icr(NCOILS) C Local variables integer n,nn,nopen integer nrow,ncol do 300 n=1,NCOILS ncs(n)=n nsc(n)=n rc(n)=icr(n)+ecr(n) r(n)=rc(n) mc(n,n)=mc(n,n)+ecl(n) do 200 nn=1,NCOILS m(n,nn)=mc(n,nn) 200 continue 300 continue C C cull out open ciruits C nopen = 0 nsim = 0 DO 800, n = 1,NCOILS IF (ntype(n).EQ.2) THEN nopen = nopen+1 C C eliminate row in m and r matricies C DO 500, nrow = nsim+1,NCOILS-nopen r(nrow) = r(nrow+1) DO 400, ncol = 1,NCOILS-nopen+1 m(nrow,ncol) = m(nrow+1,ncol) 400 CONTINUE 500 CONTINUE C C eliminate column in m matrix C DO 700, nrow = 1,NCOILS-nopen DO 600, ncol = nsim+1,NCOILS-nopen m(nrow,ncol) = m(nrow,ncol+1) 600 CONTINUE 700 CONTINUE GOTO 800 ENDIF C C circuit is to be included C C ncs contains total ckt number as function of sim circuit number C nsc contains sim circuit number as function of total circuit number C nsim = nsim+1 nsc(n) = nsim ncs(nsim) = n 800 CONTINUE C WRITE(6,*)'Number total, open, sim ckts =',nckt,nopen,nsim CALL invert(m,x,nsim) C WRITE(6,*)'matrix inversion completed' CLOSE(17) ccc CLOSE(2) ccc CLOSE(3) RETURN END C ************************************************************ SUBROUTINE invert(a,x,n) C ************************************************************ IMPLICIT NONE include 'rtcconfig.inc' C C compute the inverse of a matrix C C reference Network Analysis Theory & Computer Methods C Jenkins & Watkins C Prentice-Hall, 1974 C p. 86 C C a - input n x n matrix C x - inverse of a C n - order of a C eps - zero tolerance C C Arguments: real*4 a,x,c DIMENSION a(NCOILS,NCOILS),x(NCOILS,NCOILS),c(NCOILS,2*NCOILS) C Local variables: integer i,im,ip,ip1,ist,j,k,m,n,n1 real*4 cl C C initialize subroutine values C real*4 eps DATA eps/1.0e-12/ n1 = 2*n C C load a matrix and construct identity matrix C DO 10, i = 1,n DO 10, j = 1,n 10 c(i,j) = a(i,j) DO 30, i = 1,n DO 30, j = 1,n IF (i.EQ.j) GOTO 20 c(i,j+n) = 0. GOTO 30 20 c(i,j+n) = 1. 30 CONTINUE C C invert matrix C DO 130, ip = 1,n C C find pivot element in column ip C im = ip IF (ip.GE.n) GOTO 50 ist = ip+1 DO 40, i = ist,n IF (abs(c(im,ip)).GE.abs(c(i,ip))) GOTO 40 im = i 40 CONTINUE 50 IF (abs(c(im,ip)).GE.eps) GOTO 70 C C near diagonal element flag C WRITE(6,60)ip,c(im,ip) 60 format(1h0,17hdiagonal element ,i2,2h =,1pe13.5) IF (c(im,ip).EQ.0.) RETURN 70 IF (im.EQ.ip) GOTO 90 C C interchance rows to position pivot element C DO 80, j = ip,n1 cl = c(ip,j) c(ip,j) = c(im,j) 80 c(im,j) = cl C C find multiplication constant,set c(i,i)=1 C 90 cl = c(ip,ip) c(ip,ip) = 1. C C divide element in row by pivot element C DO 100, j = ist,n1 100 c(ip,j) = c(ip,j)/cl C C zero column of pivot element C DO 120, i = 1,n IF (i.EQ.ip) GOTO 120 ip1 = ip+1 cl = c(i,ip) c(i,ip) = 0. DO 110, j = ip1,n1 110 c(i,j) = c(i,j)-cl*c(ip,j) 120 CONTINUE 130 CONTINUE C C load inverse matrix into x C m = n+1 DO 140, i = 1,n DO 140, j = m,n1 k = j-n 140 x(i,k) = c(i,j) RETURN END C ************************************************************ SUBROUTINE lump(nin,nout,mc,r20c,ecr,ecl) C ************************************************************ implicit real*4 (a-m,o-z) include 'rtcconfig.inc' dimension mc(NCOILS,NCOILS),r20c(NCOILS),ecr(NCOILS),ecl(NCOILS) C C this subroutine lumps together the impedance of two ciruits which C will be operated in series as one C C in the mutual inductance matrix the row and column values from C the abandoned circuit are added to the row and column values of the C active circuit C C also the 20C coil resistances and the external resistances are added C C add row of circuit out to circuit in C do 100 ncol=1,NCOILS mc(nin,ncol)=mc(nin,ncol)+mc(nout,ncol) 100 continue C C add column of circuit out to circuit in C do 200 nrow=1,NCOILS mc(nrow,nin)=mc(nrow,nin)+mc(nrow,nout) 200 continue C C lump resistance of commonly controlled twin circuits C r20c(nin)=r20c(nin)+r20c(nout) ecr(nin)=ecr(nin)+ecr(nout) ecl(nin)=ecl(nin)+ecl(nout) RETURN END