summaryrefslogtreecommitdiffstats
path: root/benchmarks/dc68.dat
diff options
context:
space:
mode:
authorAngelo Rossi <angelo.rossi.homelab@gmail.com>2023-06-21 12:04:16 +0000
committerAngelo Rossi <angelo.rossi.homelab@gmail.com>2023-06-21 12:04:16 +0000
commitb18347ffc9db9641e215995edea1c04c363b2bdf (patch)
treef3908dc911399f1a21e17d950355ee56dc0919ee /benchmarks/dc68.dat
Initial commit.
Diffstat (limited to 'benchmarks/dc68.dat')
-rw-r--r--benchmarks/dc68.dat1362
1 files changed, 1362 insertions, 0 deletions
diff --git a/benchmarks/dc68.dat b/benchmarks/dc68.dat
new file mode 100644
index 0000000..f90371b
--- /dev/null
+++ b/benchmarks/dc68.dat
@@ -0,0 +1,1362 @@
+BEGIN NEW DATA CASE
+C BENCHMARK DC-68
+C Illustration of MODELS using six of Laurent Dube's simple test cases
+C 1st of 16 data subcases. This was taken from Laurent's file TC1.DAT.
+C July of 1995, this is converted to the new MODELS STAND ALONE by
+C adding STAND ALONE to the MODELS declaration. Also, all data
+C between ENDMODELS & batch-mode plot cards (just BLANK) was deleted.
+ 0.1 1.00
+ 1 1 1 1 1
+MODELS STAND ALONE
+STORAGE integer pages:0 real pages:0 ENDSTORAGE
+MODEL tc1
+ VAR a
+ INPUT i
+ EXEC
+ a:=t
+ ENDEXEC
+ENDMODEL
+USE tc1 AS tc1a
+ INPUT i:=0
+ENDUSE
+RECORD tc1a.a AS a
+ENDMODELS
+C Step Time MODELS
+C A
+C 0 0.0 0.0
+C 1 0.1 0.1
+C 2 0.2 0.2
+C 3 0.3 0.3
+C 4 0.4 0.4
+C 5 0.5 0.5
+C 6 0.6 0.6
+C 7 0.7 0.7
+C 8 0.8 0.8
+C 9 0.9 0.9
+C 10 1.0 1.0
+C Variable maxima : 1.0
+C Times of maxima : 1.0
+C Variable minima : 0.0
+C Times of minima : 0.0
+BLANK card terminating plotting
+BEGIN NEW DATA CASE
+C Illustration of MODELS using six of Laurent Dube's simple test cases
+C 2nd of 16 data subcases. This was taken from Laurent's file TC2.DAT.
+ 0.001 0.01
+ 1 1 1 1 1
+MODELS
+INPUT elec90 {v(NOD1)}
+ elec91 {i(NOD2)}
+ elec93 {switch(NOD2)}
+ elec94 {imssv(NOD1)}
+ elec95 {imssi(NOD2)}
+MODEL TC2
+ INPUT i[90..95]
+ VAR a[90..95]
+ EXEC a[90..95]:=i[90..95] ENDEXEC
+ENDMODEL
+USE tc2 AS tc2
+ INPUT i[90]:=elec90
+ i[91]:=elec91
+ i[92]:=0
+ i[93]:=elec93
+ i[94]:=elec94
+ i[95]:=elec95
+ENDUSE
+RECORD tc2.a[90] AS a90
+ tc2.a[91] AS a91
+ tc2.a[92] AS a92
+ tc2.a[93] AS a93
+ tc2.a[94] AS a94
+ tc2.a[95] AS a95
+ENDMODELS
+ NOD1 NOD2 1.00
+BLANK card after last electric network branch
+ NOD2 MEASURING
+BLANK card ends switches
+14NOD1 100.0 50.0 -30. -1.0
+BLANK card after last electric network source
+ NOD1
+C Step Time NOD1 MODELS MODELS MODELS MODELS MODELS MODELS
+C A90 A91 A92 A93 A94 A95
+C *** Phasor I(0) = 8.6602540E+01 Switch "NOD2 " to " " closed in the steady-state.
+C 0 0.0 86.6025404 86.6025404 86.6025404 0.0 1.0 -50. -50.
+C 1 .1E-2 97.8147601 97.8147601 97.8147601 0.0 1.0 0.0 0.0
+C 2 .002 99.4521895 99.4521895 99.4521895 0.0 1.0 0.0 0.0
+C 3 .003 91.3545458 91.3545458 91.3545458 0.0 1.0 0.0 0.0
+BLANK card ends selective node voltage outputs
+C 8 .008 -40.673664 -40.673664 -40.673664 0.0 1.0 0.0 0.0
+C 9 .009 -66.913061 -66.913061 -66.913061 0.0 1.0 0.0 0.0
+C 10 .01 -86.60254 -86.60254 -86.60254 0.0 1.0 0.0 0.0
+C Variable maxima : 99.4521895 99.4521895 99.4521895 0.0 1.0 0.0 0.0
+C Times of maxima : .002 .002 .002 0.0 0.0 .1E-2 .1E-2
+C Variable minima : -86.60254 -86.60254 -86.60254 0.0 1.0 -50. -50.
+C Times of minima : .01 .01 .01 0.0 0.0 0.0 0.0
+BLANK card terminating plotting
+BEGIN NEW DATA CASE
+PRINTED NUMBER WIDTH, 12, 2, { Keep dT loop precision the same, but 2 blank separators
+C Illustration of MODELS using six of Laurent Dube's simple test cases
+C 3rd of 16 data subcases. This was taken from Laurent's file TC9.DAT.
+ 0.001 0.01
+ 1 1 1 1 0
+MODELS
+OUTPUT tacs1, tacs2, arcsh, arcch, arcth
+MODEL tc9
+ INPUT i {dflt:0}
+ VAR tacs[1..2]
+ OUTPUT tacs[1..2]
+ EXEC
+ tacs[1]:=cos(2*pi*1000/6*t)
+ tacs[2]:=sin(2*pi*1000/6*(t+startstep))
+ ENDEXEC
+ENDMODEL
+USE tc9 AS tc9
+ OUTPUT tacs1:=tacs[1]
+ tacs2:=tacs[2]
+ENDUSE
+C 21 September 2002, add HYPER to demonstrate corrected hyperbolic functions.
+C Previously, the inverse sinh, cosh, and tanh were missing. Dube did not halt
+C execution, however. Instead, he returned impossible function valu 88888.88..
+C as first reported by Orlando Hevia of Santa Fe, Argentina. WSM corrects XPR2
+C Note the internal names such as arcsinhyperbolic are impractically long.
+C This is done on purpose to demonstrate Dube's flexibility in this regard.
+MODEL HYPER -- module that serves only to verify 3 inverse hyperbolic function
+ VAR arcsinhyperbolic, arccoshyperbolic, arctanhyperbolic, arg
+ OUTPUT arcsinhyperbolic, arccoshyperbolic, arctanhyperbolic
+EXEC
+ arg := 200 * T -- argument of 3 function usages to follow
+ arcsinhyperbolic := ASINH( arg ) -- test hyperbolic arc sine over [0, 2]
+ IF arg >= 1.0 -- if arg of arc cosh is possible & arc tanh is impossible
+ THEN
+ arccoshyperbolic := ACOSH( arg ) -- test hyperbolic arc cosine over ...
+ arctanhyperbolic := -9999.
+ ELSE -- Alternate, if arg of arc cosh is impossible & arc tanh is possible:
+ arccoshyperbolic := -9999.
+ arctanhyperbolic := ATANH( arg ) -- test hyperbolic arc tangent over ...
+ ENDIF
+ENDEXEC
+ENDMODEL
+USE HYPER AS HYPER
+ OUTPUT ARCSH:=arcsinhyperbolic -- ATP output name ARCSH is connected to internal model name
+ OUTPUT ARCCH:=arccoshyperbolic -- ATP output name ARCCH is connected ...
+ OUTPUT ARCTH:=arctanhyperbolic -- ATP output name ARCTH is connected ...
+ENDUSE
+RECORD tacs1 AS tacs1
+ tacs2 AS tacs2
+ arcsh AS arcsh
+ arcch AS arcch
+ arcth AS arcth
+ENDMODELS
+ NOD1 NOD2 1.00
+ TACS2 NOD2 1.00
+ NOD2 NOD3 1.00
+BLANK card after last electric network branch
+ NOD2 MEASURING
+13NOD3 TACS1 1
+BLANK card ends switches
+17TACS2
+14NOD1 100.0 50.0 -1.0
+60TACS2
+BLANK card after last electric network source
+ TACS2
+C First 1 output variables are electric-network voltage differences (upper voltage minus lower voltage);
+C Next 5 output variables belong to MODELS (with "MODELS" an internally-added upper name of pair).
+C Step Time TACS2 MODELS MODELS MODELS MODELS MODELS
+C TACS1 TACS2 ARCSH ARCCH ARCTH
+C *** Phasor I(0) = 1.0000000E+02 Switch "NOD2 " to " " closed in the steady-state.
+C 0 0.0 0.0 1.0 .866025404 0.0 -9999. 0.0
+C Switch "NOD3 " to " " closing after 1.00000000E-03 sec.
+C 1 .1E-2 .866025404 0.5 .866025404 .19869011 -9999. .202732554
+C Switch "NOD3 " to " " opening after 2.00000000E-03 sec.
+C 2 .002 .866025404 -.5 .56655E-15 .39003532 -9999. .42364893
+C 3 .003 .56655E-15 -1. -.8660254 .568824899 -9999. .693147181
+BLANK card ends selective node voltage outputs
+C 4 .004 -.8660254 -.5 -.8660254 .732668256 -9999. 1.09861229
+C Switch "NOD3 " to " " closing after 5.00000000E-03 sec.
+C 5 .005 -.8660254 0.5 -.1133E-14 .881373587 0.0 -9999.
+C 6 .006 -.1133E-14 1.0 .866025404 1.01597313 .622362504 -9999.
+C 7 .007 .866025404 0.5 .866025404 1.13798205 .867014726 -9999.
+C Switch "NOD3 " to " " opening after 8.00000000E-03 sec.
+C 8 .008 .866025404 -.5 .36738E-15 1.24898333 1.04696792 -9999.
+C 9 .009 .36738E-15 -1. -.8660254 1.35044074 1.19291073 -9999.
+C 10 .01 -.8660254 -.5 -.8660254 1.44363548 1.3169579 -9999.
+BLANK card terminating plotting
+BEGIN NEW DATA CASE
+PRINTED NUMBER WIDTH, 11, 1, { Cancel preceding declaration (extra blank)
+C Illustration of MODELS using six of Laurent Dube's simple test cases
+C 4th of 16 data subcases. This was taken from Laurent's file TC10.DAT
+ 1.0 1.0
+ 1 1 1 1 1
+MODELS
+MODEL test
+COMMENT
+ Provides a quick check on the use of nested pointlist functions
+ for the representation of parametric curves.
+ Here the parameter is the variable 'r',
+ and the main functin is used as y=f(x,r)
+ using a=mainfun(x) at various values of r.
+ENDCOMMENT
+ FUNCTION fun1 POINTLIST (-inf,inf), (0,0), (inf,-inf)
+ FUNCTION fun2 POINTLIST (-inf,2*inf), (0,0), (inf,-2*inf)
+ FUNCTION fun3 POINTLIST (-inf,3*inf), (0,0), (inf,-3*inf)
+ FUNCTION fun4 POINTLIST (-inf,4*inf), (0,0), (inf,-4*inf)
+ FUNCTION mainfun POINTLIST
+ (0,0),(1,fun1(r)),(2,fun2(r)),(3,fun3(r)),(4,fun4(r))
+ FUNCTION ref(x,r):=-r*x --Same function in its analytic form
+ INPUT i
+ VAR a,r
+ INIT
+ r:=-1
+ WHILE r<5 DO
+ r:=r+1
+ write(' ')
+ write('r:',r)
+ FOR k:=0 to 3.5 by 0.5 DO
+ a:=mainfun(k)
+ write(' ',k,': ',a, ' ', ref(r,k))
+ ENDFOR
+ ENDWHILE
+ ENDINIT
+ EXEC
+ ENDEXEC
+ENDMODEL
+USE test AS test
+ INPUT i:=0
+ENDUSE
+ENDMODELS
+ NOD1 1.00
+BLANK card after last electric network branch
+BLANK card ends switches
+14NOD1 100.0 50.0 -1.0
+BLANK card after last electric network source
+ NOD1
+BLANK card ends selective node voltage outputs
+C Step Time NOD1
+C 0 0.0 100.
+C 1 1.0 100.
+C Variable maxima : 100.
+C Times of maxima : 0.0
+C Variable minima : 100.
+C Times of minima : 0.0
+BLANK card terminating plotting
+BEGIN NEW DATA CASE
+C Illustration of MODELS using six of Laurent Dube's simple test cases
+C 5th of 16 data subcases. This was taken from Laurent's file TC11.DAT
+ 0.1 1.00
+ 1 1 1 1 1
+MODELS
+MODEL tc11
+ VAR a
+ INPUT i
+ EXEC
+ a:=t
+ write ( 'Writing a line...' )
+ ENDEXEC
+ENDMODEL
+USE tc11 AS tc11
+ INPUT i:=0
+ENDUSE
+RECORD tc11.a AS a
+ENDMODELS
+ NOD1 1.00
+BLANK card after last electric network branch
+BLANK card ends switches
+14NOD1 100.0 50.0 -1.0
+BLANK card after last electric network source
+ NOD1
+C Step Time NOD1 TACS
+C A
+C 0 0.0 100. 0.0
+C 1 0.1 100. 0.1
+C 2 0.2 100. 0.2
+C 3 0.3 100. 0.3
+BLANK card ends selective node voltage outputs
+C 9 0.9 100. 0.9
+C 10 1.0 100. 1.0
+C Variable maxima : 100. 1.0
+C Times of maxima : 0.0 1.0
+C Variable minima : 100. 0.0
+C Times of minima : 0.0 0.0
+BLANK card terminating plotting
+BEGIN NEW DATA CASE
+C Illustration of MODELS using six of Laurent Dube's simple test cases
+C 6th of 16 data subcases. This was taken from Laurent's file TC13.DAT
+ 0.1 1.00
+ 1 1 1 1 1
+MODELS
+MODEL test
+ INPUT i
+ VAR a, b
+ HISTORY b
+ EXEC
+ a:=b {MAX:0.7}
+ b:=t
+ ENDEXEC
+ENDMODEL
+USE test AS test
+ INPUT i:=t
+ HISTORY b:=t
+ENDUSE
+RECORD test.a AS a
+ENDMODELS
+ NOD1 1.00
+BLANK card after last electric network branch
+BLANK card ends switches
+14NOD1 100.0 50.0 0.0
+BLANK card after last electric network source
+ NOD1
+C Step Time NOD1 TACS
+C A
+C 0 0.0 0.0 -.1
+C 1 0.1 100. 0.0
+C 2 0.2 100. 0.1
+C 3 0.3 100. 0.2
+C 4 0.4 100. 0.3
+BLANK card ends selective node voltage outputs
+C 9 0.9 100. 0.7
+C 10 1.0 100. 0.7
+C Variable maxima : 100. 0.7
+C Times of maxima : 0.1 0.9
+C Variable minima : 0.0 -.1
+C Times of minima : 0.0 0.0
+BLANK card terminating plotting
+BEGIN NEW DATA CASE
+C 7th of 16 subcases was added 22 August 1998 to illustrate a supplemental
+C variable that extends through column 80. Furthermore, it is unusual in
+C that it has no blank byte to the right of the equal sign in column 11.
+C The result is easily checked since, by inspection, it has value 9 + T.
+C For background text, see the July, 1998, newsletter. Dr. Yuan Bin called
+C our attention to this problem. But that was not the end. It was about a
+C year later that Carlos Mata pointed to the weakness of the correction: it
+C failed if column 80 contained a closed parenthesis. The correction was
+C fortified March 18, 1998, and new variable CARLOS confirms correct
+C operation. Note that the values of I1FG50 and CARLOS are equal.
+ABSOLUTE TACS DIMENSIONS
+ 10 90 100 20 30 400 350 60
+ .02 .10
+ 1 1 0 0
+TACS STAND ALONE
+99I1FG00 = 1.
+99I1FG02 = TIMEX
+99I1FG04 = 1.
+99I1FG06 = 1.
+99I1FG08 = 1.
+99I1FG10 = 1.
+99I1FG12 = 1.
+99I1FG14 = 1.
+99I1FG16 = 1.
+99I1FG18 = 1.
+98YUAN =I1FG00+I1FG02+I1FG04+I1FG06+I1FG08+I1FG10+I1FG12+I1FG14+I1FG16+I1FG18
+98CARLOS =I1FG00+I1FG02+I1FG04+I1FG06+I1FG08+I1FG10+I1FG12+ 1.5*(I1FG14+I1FG16)
+33YUAN CARLOS
+BLANK card ending all TACS data cards
+BLANK card ending plot cards
+BEGIN NEW DATA CASE
+C 8th of 16 subcases was added 18 August 1998 to illustrate new request words
+C that are related to repeatability of random numbers, or lack of it. Here
+C we illustrate repeatability, but the non-repeatable case really is the more
+C important. If the REPEATABLE card below is replaced by the comment card
+C that immediately precedes it, the new initialization of OVER8 (see the
+C October newsletter) is demonstrated. Previously, ATP GNU answers differed
+C from Salford and Watcom answers. Also illustrated is no output variable,
+C with MODELS doing all output. Of course, with only 2 energizations and
+C a dummy electric network, this is just a demonstrative shadow of the real
+C data case of interest (Gabor Furst's BACKFL.DAT). The dominant structural
+C feature is use of MODELS for all output of a SYSTEMATIC data case. There
+C _is_ a base case, although the engineering significance is not known.
+TRULY RANDOM NUMBERS { Alternative to following card removes repeatability
+REPEATABLE RANDOM NUMBERS { English-language override of NSEED 5 lines below
+C Note about preceding: Since the TRULY declaration is overridden by the
+C REPEATABLE declaration, in fact the TRULY declaration does nothing. But
+C it does illustrate acceptance by ATP.
+PRINTED NUMBER WIDTH, 13, 2, { Request maximum precision (for 8 output columns)
+ .060 .300 { Time step is artificially large (5 steps to completion)
+ 1 1 1 2 1 -1 -2
+C ISW ITEST IDIST IMAX IDICE KSTOUT NSEED
+ 0 1 0 0 1 1
+ 5 5 20 20 { Printout frequency changes for base case
+MODELS HYBRID
+STORAGE integer pages:0 real pages:0 ENDSTORAGE
+MODEL tc1
+ VAR angle
+ INPUT i
+ INIT angle:=0
+ ENDINIT
+ EXEC
+ angle:= random()
+C write(' T =', t, ' MODELS random() yields =', angle )
+--- Only program output within dT loop is preceding output of variable ANGLE ---
+C To preceding output, add the energization number KNT to demonstrate
+C that a 4-byte integer is handled properly by XPR1 for WRITE(. GNU ATP
+C uses different code of IEIGHT for this than for ISTEP. WSM, 28 Aug 2010
+ write(' T =', t, 'KNT =', atp(knt), 'MODELS random() ==>', angle )
+ ENDEXEC
+ENDMODEL
+USE tc1 AS tc1a
+ INPUT i:=0
+ENDUSE
+C RECORD tc1a.angle AS angle { If activates, variable angle becomes ATP output
+ENDMODELS
+ GENA A1 1.
+ A1 1.
+ ASW1 1.
+BLANK card ending branch cards
+ A1 ASW1 4.E-3 .2E-3 2 SYSTEMATIC
+BLANK card ending switch cards
+14GENA 303. 60. 0.0 -1.
+BLANK card ending source cards
+BLANK card ending output variable requests (node voltages only, here)
+ PRINTER PLOT
+BLANK card ending plot cards (for base case only)
+BLANK card ending statistical tabulate (none possible, since no output variable)
+BEGIN NEW DATA CASE
+C 9th of 16 subcases is a meaningful case from Trainee Rod Price of BPA, who
+C assembled and tested the data during the summer and fall of 1989. The
+C study involves the 500-kV Coulee-Raver overhead line, and the model
+C includes series capacitors at Columbia. These capacitors are protected
+C by nonlinear resistors (ZnO modeling) that are monitored and bypassed
+C using MODELS. To see more interesting simulation, the end-time TMAX
+C should be increased to 100 msec. Reformulated by Laurent Dube 8 Aug 91
+C following his correction of a MODELS bug for Bob Wilson at Univ. of Idaho
+C in Moscow. It was observed that the answer of this 7-th subcase changed,
+C too (i.e., before it was wrong, although perhaps not by much)!
+C 14 Sept 1993, WSM adds $BLANK DATA below to protect the user who has
+C NOBLAN of STARTUP equal to unity (equivalent to $BLANK HALT to halt if
+C a true blank card is found. This MODELS data does use blank lines as
+C comments, and MODELS ignores the blanks. The important thing is to
+C prevent EMTP data input from halting before MODELS receives the blanks!
+C 14 October 2000, the addition of one IF-statement to OVER16 (see 8th
+C subcase of DCNEW-16 and January, 2001, newsletter) changes the answers
+C of this subcase slightly. To document this, output variables are have
+C been made more relevant. Differences follow the fault to ground. From
+C the .DBG file, the new logic is used 3 times:
+C OVER16 bypasses compensation. J, T = 43 2.000000000000E-02
+C OVER16 bypasses compensation. J, T = 43 3.150000000000E-02
+C OVER16 bypasses compensation. J, T = 43 4.990000000000E-02
+C It is not surprising that the difference shows up after step 400, which
+C corresponds to that switch closure at T = 20 msec (1st of 3 uses).
+C The following comparison shows how extrema are affected only slightly:
+C Next 3 output variables are branch currents (flowing from the upper node to the lower node);
+C Next 5 output variables belong to MODELS (with "MODELS" an internally-added upper name of pair).
+C Step Time CR25A CR20A CRZ1A MODELS MODELS MODELS MODELS MODELS
+C CR30A TERRA CR20A VCAPA IZNA PZNA EZNA GAPA
+C Old maxima : 27237.7816 15421.0962 .445513E-4 39549.6402 .445515E-4 .139454E10 .500879 E7 1.0
+C New maxima : 26931.1785 15312.52 .445513E-4 39549.6402 .445515E-4 .135464E10 .4929451E7 1.0
+C Times of both : .03175 .0375 .0125 .0125 .0125 .028 .05 .0315
+C Old minima : -7552.2441 -28033.65 -7291.149 -191594.73 -7291.149 .432745E-9 0.0 0.0
+C New minima : -7466.2202 -28099.511 -7350.3614 -184341.43 -7350.3614 .432745E-9 0.0 0.0
+C Times of old : .0378 .02865 .02815 .02815 .02815 0.0 0.0 0.0
+C Times of new : .0378 .02865 .02805 .02805 .02805 0.0 0.0 0.0
+C About switching times, both have: Gap "CR25A " to "CR30A " closing
+C after 3.15000000E-02 sec. But openings differ by one time step:
+C Old: Gap "CR25A " to "CR30A " opening after 4.98500000E-02 sec.
+C New: Gap "CR25A " to "CR30A " opening after 4.99000000E-02 sec.
+PRINTED NUMBER WIDTH, 12, 2, { Cancel wider output of preceding subcase
+$BLANK DATA { Request old treatment of true blanks (not
+ .000050 .050 60.
+ 1 9 1 3 1 -1
+ 5 5 20 20 100 100
+C
+MODELS HYBRID
+
+INPUT -- declare names and function of all MODELS inputs from EMTP
+
+ vc3a {v(CR30A)} -- vc3a = voltage at node CR30A
+ vc2a {v(CR20A)} -- vc2a = voltage at node CR20A
+ vc3b {v(CR30B)} -- vc3b = voltage at node CR30B
+ vc2b {v(CR20B)} -- vc2b = voltage at node CR20B
+ vc3c {v(CR30C)} -- vc3c = voltage at node CR30C
+ vc2c {v(CR20C)} -- vc2c = voltage at node CR20C
+ iznoa {i(CRZ2A)} -- iznoa = current through switch at node CRZ2A
+ iznob {i(CRZ2B)} -- iznob = current through switch at node CRZ2B
+ iznoc {i(CRZ2C)} -- iznoc = current through switch at node CRZ2C
+
+OUTPUT -- declare names of MODELS outputs to EMTP
+
+ gapa, gapb, gapc -- switch control signals
+
+MODEL flash
+
+comment ****************************************************************
+* *
+* Function: set or cancel the gap firing control signal *
+* Inputs : voltage and current across ZnO resistor *
+* Output : the firing signal to the electrical ZnO component *
+* *
+************************************************************* endcomment
+
+ INPUT vcap -- voltage difference across series caps [Volts]
+ iczn -- ZNO current [Amps]
+
+ DATA Pset -- power setting [Megajoules/msec]
+ Eset -- energy setting [Megajoules]
+ fdel -- firing delay [msec]
+ fdur -- firing duration [msec]
+
+ VAR power -- power into ZnO resistor [Watts]
+ trip -- gap firing control signal [0 or 1]
+ energy -- energy into ZnO resistor [Joules]
+ tfire -- time at which the gap was last fired [sec]
+
+ OUTPUT trip
+
+ HISTORY INTEGRAL(power) {DFLT:0}
+
+ INIT trip:=0
+ tfire:=0
+ ENDINIT
+
+ EXEC
+ ------------------------------------------------------------------
+ power:=vcap*iczn
+ energy:=INTEGRAL(power)
+ ------------------------------------------------------------------
+ IF trip>0 -- is already firing
+ AND t-tfire>fdur*1.e-3 -- has exceeded firing duration
+ THEN
+ trip:=0 -- cancel the firing signal
+ tfire:=0 -- null the firing time
+ ENDIF
+ ------------------------------------------------------------------
+ IF trip=0 -- is not signaling to fire
+ AND tfire=0 -- firing condition not yet detected
+ AND ( power >= Pset * 1.e9 -- power setting exceeded
+ OR energy >= Eset * 1.e6 ) -- energy setting exceeded
+ THEN
+ tfire:=t -- set the firing detection time
+ ENDIF
+ ------------------------------------------------------------------
+ IF trip=0 -- is not signaling to fire
+ AND tfire>0 -- firing condition has been detected
+ AND t-tfire>=fdel*1.e-3 -- firing delay exceeded
+ THEN
+ trip:=1 -- set the firing signal
+ ENDIF
+ ------------------------------------------------------------------
+ ENDEXEC
+
+ENDMODEL
+
+USE flash AS flasha ----------------------------------- PHASE A ------
+
+ TIMESTEP MIN: 0.0005
+
+ INPUT vcap:=vc3a-vc2a -- voltage across caps
+ iczn:=iznoa -- current through ZnO resistor
+
+ DATA Pset:= 1 -- power setting [MJ/mS]
+ Eset:= 9 -- energy setting [MJ]
+ fdel:= 4 -- firing delay [msec]
+ fdur:=20 -- firing duration [msec]
+
+ OUTPUT gapa:=trip -- passes value of 'trip' to MODELS variable 'GAPA'
+
+ENDUSE -----------------------------------------------------------------
+
+USE flash AS flashb ----------------------------------- PHASE B ------
+
+ TIMESTEP MIN: 0.0005
+
+ INPUT vcap:=vc3b-vc2b -- voltage across caps
+ iczn:=iznob -- current through ZnO resistor
+
+ DATA Pset:= 1 -- power setting [MJ/mS]
+ Eset:= 9 -- energy setting [MJ]
+ fdel:= 4 -- firing delay [msec]
+ fdur:=20 -- firing duration [msec]
+
+ OUTPUT gapb:=trip -- passes value of 'trip' to MODELS variable 'GAPB'
+
+ENDUSE -----------------------------------------------------------------
+
+USE flash AS flashc ----------------------------------- PHASE C ------
+
+ TIMESTEP MIN: 0.0005
+
+ INPUT vcap:=vc3c-vc2c -- voltage across caps
+ iczn:=iznoc -- current through ZnO resistor
+
+ DATA Pset:= 1 -- power setting [MJ/mS]
+ Eset:= 9 -- energy setting [MJ]
+ fdel:= 4 -- firing delay [msec]
+ fdur:=20 -- firing duration [msec]
+
+ OUTPUT gapc:=trip -- passes value of 'trip' to MODELS variable 'GAPC'
+
+ENDUSE -----------------------------------------------------------------
+
+RECORD flasha.vcap AS VCAPA
+ flasha.iczn AS IZNA
+ flasha.power AS PZNA
+ flasha.energy AS EZNA
+ flasha.trip AS GAPA
+
+ENDMODELS
+C --------------------- End MODELS data, begin electric network:
+C **** Series Capacitors, Xc = 28.4 ohms
+C
+ CR20A CR30A 93.40
+ CR20B CR30B 93.40
+ CR20C CR30C 93.40
+C
+ CR20D CR30D 93.40
+ CR20E CR30E 93.40
+ CR20F CR30F 93.40
+C
+C **** Jumpers to help EMTP convergence with ZNO
+C
+ CRZ1A CRZ2A 0.01
+ CRZ1B CRZ2B 0.01
+ CRZ1C CRZ2C 0.01
+C
+ CRZ1D CRZ2D 0.01
+ CRZ1E CRZ2E 0.01
+ CRZ1F CRZ2F 0.01
+C
+C **** ZNO protection for Caps, I=1.8KA, 2.2 PU, alpha = 40, Vref =147.5
+C
+92CRZ1A CR20A 5555. 1
+ 147500. -1.
+ 1.0 40. .80
+ 9999.
+92CRZ1B CR20B CRZ1A CR20A 5555.
+92CRZ1C CR20C CRZ1A CR20A 5555.
+C
+92CRZ1D CR20D CRZ1A CR20A 5555.
+92CRZ1E CR20E CRZ1A CR20A 5555.
+92CRZ1F CR20F CRZ1A CR20A 5555.
+C
+C **** Source Impedances
+C
+51SRC1A RAVBA 0.1 20.0
+52SRC1B RAVBB 0.1 18.5
+53SRC1C RAVBC
+C
+51SRC1A RAVBA 300.0
+52SRC1B RAVBB 150.0
+53SRC1C RAVBC
+C
+C
+51SRC2A GRCBA 0.1 7.00
+52SRC2B GRCBB 0.1 10.7
+53SRC2C GRCBC
+C
+51SRC2A GRCBA 350.0
+52SRC2B GRCBB 150.0
+53SRC2C GRCBC
+C
+C *** LINE IMPEDANCES Values have been 1/4 to get high fault MVA ****
+C
+51CR50A CR30A 0.5 14.0
+52CR50B CR30B 8. 52.0
+53CR50C CR30C
+C
+51CR50D CR30D 0.5 14.0
+52CR50E CR30E 8. 52.0
+53CR50F CR30F
+C
+51CR01A CR20A 0.4 9.0
+52CR01B CR20B 5. 36.0
+53CR01C CR20C
+C
+51CR01D CR20D 0.4 9.0
+52CR01E CR20E 5. 36.
+53CR01F CR20F
+C
+C ***** BYPASS IMPEDANCE ****
+ CR20A CR25A 5. 0.23
+ CR20B CR25B 5. 0.23
+ CR20C CR25C 5. 0.23
+ CR20D CR25D 5. 0.23
+ CR20E CR25E 5. 0.23
+ CR20F CR25F 5. 0.23
+C **** PARALLEL RESISTOR
+ CR20A CR25A 200.
+ CR20B CR25B 200.
+ CR20C CR25C 200.
+ CR20D CR25D 200.
+ CR20E CR25E 200.
+ CR20F CR25F 200.
+BLANK card ending program branch cards.
+C **** TACS CONTROLLED CAP BYPASS SWITCHES *******
+C ** SWITCH 1
+12CR25A CR30A GAPA 11
+12CR25B CR30B GAPB 1
+12CR25C CR30C GAPC 1
+C 12CR25D CR30D GAPD 1
+C 12CR25E CR30E GAPE
+C 12CR25F CR30F GAPF
+ CR25D CR30D -1.0 10.
+ CR25E CR30E -1. 10.
+ CR25F CR30F -1. 10.
+C
+C **** Main Breaker Switches
+C
+ RAVBA CR50A -1.0 10.0
+ RAVBB CR50B -1.0 10.0
+ RAVBC CR50C -1.0 10.0
+C
+ RAVBA CR50D -1.0 10.0
+ RAVBB CR50E -1.0 10.0
+ RAVBC CR50F -1.0 10.0
+C
+ GRCBA CR01A -.006 10.0
+ GRCBB CR01B -.006 10.0
+ GRCBC CR01C -.006 10.0
+C
+ GRCBA CR01D -1.0 10.0
+ GRCBB CR01E -1.0 10.0
+ GRCBC CR01F -1.0 10.0
+C
+C *** ZNO current measuring switches ***
+C
+ CR30A CRZ2A -1. 1.
+ CR30B CRZ2B -1. 1.
+ CR30C CRZ2C -1. 1.
+ CR30D CRZ2D -1. 1.
+ CR30E CRZ2E -1. 1.
+ CR30F CRZ2F -1. 1.
+ CR20A 0.01998 10. { Fault switch, phase "a" to ground } 1
+BLANK card terminating program switch cards
+14SRC1A 440000. 60. -20.0 -1. 10.
+14SRC1B 440000. 60. 100.0 -1. 10.
+14SRC1C 440000. 60. 220.0 -1. 10.
+14SRC2A 440000. 60. 0.0 -1. 10.
+14SRC2B 440000. 60. 120.0 -1. 10.
+14SRC2C 440000. 60. 240.0 -1. 10.
+C --------------+------------------------------
+C From bus name | Names of all adjacent busses.
+C --------------+------------------------------
+C CR20A |TERRA *CR30A *CRZ1A *CR01A *CR25A *CR25A *
+C CR30A |CR20A *CRZ2A *CR50A *CR25A *
+C CR20B |CR30B *CRZ1B *CR01B *CR25B *CR25B *
+C CR30B |CR20B *CRZ2B *CR50B *CR25B *
+C CR20C |CR30C *CRZ1C *CR01C *CR25C *CR25C *
+C CR30C |CR20C *CRZ2C *CR50C *CR25C *
+C CR20D |CR30D *CRZ1D *CR01D *CR25D *CR25D *
+C CR30D |CR20D *CRZ2D *CR50D *CR25D *
+C Total network loss P-loss by summing injections = 8.848099600964E+07
+C Output for steady-state phasor switch currents.
+C Node-K Node-M I-real I-imag I-magn
+C CR25A CR30A Open Open Open
+C CR25B CR30B Open Open Open
+C CR25C CR30C Open Open Open
+C CR25D CR30D 9.29172435E+02 -1.77481409E+02 9.45970963E+02
+C CR25E CR30E -3.10882808E+02 8.93427637E+02 9.45970963E+02
+C CR25F CR30F -6.18289626E+02 -7.15946228E+02 9.45970963E+02
+BLANK card terminating program source cards.
+C End inject: SRC2C -220000. 440000. -1191.301327068 2345.3572211208
+C End inject: -381051.1776652 -120.0000000 -2020.272665456 -120.5267278
+C Step Time CR25A CR25A TACS TACS TACS TACS
+C CR30A CR30A VCAPA IZNA PZNA EZNA
+C *** Phasor I(0) = 9.2917243E+02 Switch "CR25D " to "CR30D " closed
+C *** Phasor I(0) = -3.1088281E+02 Switch "CR25E " to "CR30E " closed
+C *** Phasor I(0) = -6.1828963E+02 Switch "CR25F " to "CR30F " closed
+C < < Etc. (many switch currents) > >
+C 0 0.0 -143.44999 0.0 143.449993 .28697E-11 .411659E-9 0.0
+C 1 .5E-4 602.099261 0.0 -602.09926 -.68002E-6 .411659E-9 0.0
+C 2 .1E-3 1347.43459 0.0 -1347.4346 -.15198E-5 .411659E-9 0.0
+C 3 .15E-3 2092.2912 0.0 -2092.2912 -.23592E-5 .411659E-9 0.0
+BLANK card ending output variables requests (none here, since all column 80)
+C 300 .015 -23363.319 0.0 23363.3193 .263161E-4 .614832009 .014309206
+C *** Close switch "CR20A " to " " after 2.00000000E-02 sec.
+C 400 .02 37570.0647 0.0 -37570.065 -.42322E-4 1.59004745 .016928315
+C 500 .025 56444.8022 0.0 -56444.802 -.63583E-4 3.58891722 .022232887
+C 600 .03 183226.734 0.0 -183226.73 -5786.1028 .106017E10 .3650637E7
+C Gap "CR25A " to "CR30A " closing after 3.15000000E-02 sec.
+C 700 .035 0.0 -4433.8226 22742.1728 .256182E-4 .582612436 .4929451E7
+C 800 .04 0.0 -5450.0316 25576.0355 .288104E-4 .736855415 .4929451E7
+C 900 .045 0.0 4627.82167 -23165.41 -.26095E-4 .604499203 .4929451E7
+C Gap "CR25A " to "CR30A " opening after 4.99000000E-02 sec.
+C 1000 .05 -816.5929 0.0 3143.46664 .354099E-5 .011130982 .4929451E7
+C Variable max:184341.425 26931.1785 39549.6402 .445516E-4 .135464E10 .4929451E7
+C Times of max: .02805 .03175 .0125 .0125 .028 .05
+C Variable min: -39549.64 -7466.2202 -184341.43 -7350.3614 .411659E-9 0.0
+C Times of min: .0125 .0378 .02805 .02805 0.0 0.0
+ PRINTER PLOT
+ 194 5. 0.0 40. CR25A CR30A { Axis limits: (-0.747, 2.220)
+BLANK card ending plot cards
+BEGIN NEW DATA CASE
+PRINTED NUMBER WIDTH, 12, 2, { Keep dT loop precision the same, but 2 blank separators
+C 10th of 16 subcases was added 10 April 2011 to illustrate very simple use
+C of POINTLIST within MODELS. Data comes from Orlando Hevia of Santa Fe,
+C Argentina. It establishes the standard of comparison for his alternative
+C SEEK function (see 11th subcase, which produces the same answer). Here
+C the (X,Y) points of POINTLIST number only 7. For much larger numbers of
+C points (e.g., 20200 for Orlando's TFGNMOD.DAT), Dube's POINTLIST is much
+C too slow. For much larger numbers of points, use Orlando's new SEEK as
+C will be illustrated by the following subcase. WSM.
+ 1.0 10.00 { Take ten steps, from 1 second to 10 seconds. This is "T"
+ 1 -1 { Print every step; do not bother saving points for plotting
+MODELS STAND ALONE
+MODEL SAMPLE
+VAR
+ X,Y,Z
+FUNCTION R POINTLIST
+(0, 69)
+(2, 50)
+(4, 142)
+(5, 142)
+(6, 188)
+(8, 73)
+(10, 100)
+INIT
+ X:=0
+ENDINIT
+EXEC
+X:=R(T)
+Y:=SQRT(X)*0.3
+Z:=R(SQRT(X)*0.3)
+ENDEXEC
+ENDMODEL
+USE SAMPLE AS SAMPLE
+ENDUSE
+RECORD
+ SAMPLE.X AS X
+ SAMPLE.Y AS Y
+ SAMPLE.Z AS Z
+ENDRECORD
+ENDMODELS
+C Column headings for the 3 EMTP output variables follow. These are divided among the 5 possible classes as follows ....
+C Next 3 output variables belong to MODELS (with "MODELS" an internally-added upper name of pair).
+C Step Time MODELS MODELS MODELS
+C X Y Z
+C 0 0.0 69. 2.49198716 72.6314093
+C 1 1.0 59.5 2.31408729 64.4480155
+C 2 2.0 50. 2.12132034 55.5807358
+C 3 3.0 96. 2.93938769 93.2118338
+C 4 4.0 142. 3.57491259 122.445979
+C 5 5.0 142. 3.57491259 122.445979
+C 6 6.0 188. 4.11339276 144.608033
+C 7 7.0 130.5 3.4270979 115.646503
+C 8 8.0 73. 2.56320112 75.9072517
+C 9 9.0 86.5 2.79016129 86.3474191
+C 10 10. 100. 3.0 96.
+BLANK CARD ENDING PLOT CARDS
+BEGIN NEW DATA CASE
+$DEPOSIT, JSEEDR=0 { Reinitialize random # generator as if no preceding use
+C Note: the preceding $DEPOSIT is required because random # use here is
+C not the 1st of data case. Preceding subcases roll dice. Memory
+C of this is erased by setting JSEEDR to zero. WSM. 29 July 2011
+C This is for addition of GAUSS. Variable NORMAL below should be
+C identical in value to what is produced by TACS within subcase 2
+C of BENCHMARK DC-18. This works well for this subcase, but later
+C subcases then would have different random numbers. In order not
+C to change the random numbers of following subcases, we want to
+C restore at the end the JSEEDR value that exists at the beginning.
+C New logic of DEPOSI makes this possible. See 2nd $DEPOSIT below.
+$PREFIX, [] { $INCLUDE files are located in same place as this main data file
+C Preceding line allows the SEEK file to be remote if this main data file is
+PRINTED NUMBER WIDTH, 12, 2, { Keep dT loop precision the same, but 2 blank separators
+C 11th of 16 subcases was added 10 April 2011 to illustrate that Orlando
+C Hevia's new SEEK function of MODELS will produce the same answer as
+C the preceding 10th subcase, which uses POINTLIST. Data is too small for
+C there to be a reason to use SEEK instead of POINTLIST, of course. For
+C a more realistic illustration of the advantage of SEEK, see explanation
+C of the EEUG list server dated 11 June 2011. This includes the following
+C paragraph:
+C The inadequacy of MODELS POINTLIST for large data (many points of the
+C piecewise-linear X-Y function) was documented in detail in a paper that
+C was presented at the year-2006 EEUG conference. See "Determination of
+C an optimized energy storage size for a wind farm based on wind forecasts"
+C by Steve Völler, Christian Müller, and Johannes Verstege of the
+C University of Wuppertal in Germany. In this case, POINTLIST functioned
+C correctly, but was too slow to be practical. Those German users give new
+C meaning to the term industrial-strength. Would a reader believe 105K
+C pairs of (X, Y)? Yes, letter K as in one thousand. Mr. Dube obviously
+C did not imagine such use when he programmed his original logic. But
+C Orlando Hevia's new SEEK now provides an efficient alternative without
+C the need for compilation and linking of user-supplied source code.
+C This data too comes from Orlando Hevia of Santa Fe, Argentina. WSM.
+ 1.0 10.00 { Take ten steps, from 1 second to 10 seconds. This is "T"
+ 1 -1 { Print every step; do not bother saving points for plotting
+MODELS STAND ALONE
+C 3456789012345678901234567890123456789012345678 -- Ruler for following optional
+C < File Name> IPRMDL KOMPAC miscellaneous data card :
+MODELS MISC. DATA dc68seek. 0 0
+C Note about preceding. Special value KOMPAC = 8765 is required
+C to bypass the check on (X,Y) points of the function. Logic will
+C ensure that X never decreases (a fatal error). It also discards
+C all but the end points of each string of three or more identical
+C Y values. That is the 3rd of 3 params shown. As for the first,
+C see the substantial note below. Finally, the middle parameter,
+C IPRMDL, is just IPRSUP for the new MODELS code. I.e., it
+C provides diagnostic printout control. All 3 values will remain
+C in force unless redefined. So, for stacked subcases, there is
+C no need for repeated definition unless values are to be changed.
+MODEL SAMPLE
+VAR X,Y,Z, NORMAL
+FUNCTION seek FOREIGN seek {ixarg: 2}
+FUNCTION gauss FOREIGN gauss {ixarg: 2}
+EXEC
+C Of the three following EXEC lines, the first and the third use Orlando's
+C new SEEK function which is defined by the (X, Y) pairs of DC68SE37.DAT
+C where the root name "DC68SEEK." was read from the MODELS misc. data card.
+C The "37" that overwrites the "EK" of DC68SEEK corresponds to the second
+C argument of SEEK. If another function were involved, it would correspond
+C to a different subscript (any # < 100 million is legal). An arbitrary # of
+C such functions is allowed, and each can involve an arbitrary # of (X, Y)
+C as long as List 20 has adequate unused space. As for .DAT being the file
+C type, this follows from DATTYP of STARTUP because the user's file name
+C ended with a period. If there is no MODELS miscellaneous data, the default
+C disk file name SEEKDATA.DAT will be assumed by ATP. WSM.
+X:=SEEK(T, 37 )
+Y:=SQRT(X)*0.3
+Z:=SEEK(SQRT(X)*0.3, 37 )
+C In order that all data be seen, let's document the content of DC68SE37.DAT
+C which contains just the following 7 formatted lines:
+C 0, 69
+C 2, 50
+C 4, 142
+C 5, 142
+C 6, 188
+C 8, 73
+C 10, 100
+C This is just usual free-format data. Each line begins in column 1, and
+C the final line is known as ATP encounters an end-of-file mark. Neither
+C comment lines nor $-cards nor inline comments are allowed. This is not
+C an oversight. Rather, it is a way to speed execution for large files.
+C If the user wants to terminate his function prematurely, however, he
+C is allowed to add a line such as "-9999., 0" where the -9999 is
+C interpreted as a software EOF. For an illustration, see subcase # 15.
+C For some reason, Dube does not accept an in-line comment on the
+C following function. Perhaps not on any function? So, strip it off:
+C NORMAL:=gauss( T ) { Argument "T" in fact is unused, so is arbitrary
+NORMAL:=gauss( T )
+ENDEXEC
+ENDMODEL
+USE SAMPLE AS SAMPLE
+ENDUSE
+RECORD
+SAMPLE.X AS X
+SAMPLE.Y AS Y
+SAMPLE.Z AS Z
+SAMPLE.NORMAL AS NORMAL
+ENDRECORD
+ENDMODELS
+C MODELS "SEEK" function # 37 involves 7 pairs of (X,Y). These are stored in List-20 SCONST cells 3 through 16.
+C Column headings for the 4 EMTP output variables follow. These are divided among the 5 possible classes as follows ....
+C Next 4 output variables belong to MODELS (with "MODELS" an internally-added upper name of pair).
+C Step Time MODELS MODELS MODELS MODELS
+C X Y Z NORMAL
+C 0 0.0 69. 2.49198716 72.6314093 .890624441
+C 1 1.0 59.5 2.31408729 64.4480155 -.26382016
+C 2 2.0 50. 2.12132034 55.5807358 -.41782888
+C 3 3.0 96. 2.93938769 93.2118338 .66418476
+C 4 4.0 142. 3.57491259 122.445979 .242122394
+C 5 5.0 142. 3.57491259 122.445979 -.38464212
+C 6 6.0 188. 4.11339276 142. -.9889959
+C 7 7.0 130.5 3.4270979 115.646503 -.26589058
+C 8 8.0 73. 2.56320112 75.9072517 -.73063403
+C 9 9.0 86.5 2.79016129 86.3474191 -.75107385
+C 10 10. 100. 3.0 96. .288485151
+C Size 11-20: 0 3 -9999 -9999 1744 0 0 0 23 15
+C The preceding line from 80-column case-summary statistics confirms that
+C 15 cells of List 20 storage are being used by SEEK. There are 6 points
+C for a total of 12 values. But for each function, there are two extra
+C cells. Finally, the last function must be bounded by a blank, so the
+C total is 2 * 6 + 2 + 1 = 15 as shown.
+C That was before adding GAUSS use. After the addition, this becomes:
+C Size 11-20: 0 4 -9999 -9999 1838 0 0 0 23 17
+$DEPOSIT, JSEEDR= -98789 { Flag for DEPOSI to restore value before previous use
+C Note about preceding line: As this use begins on 31 July 2011, service is
+C limited to a single variable (in this case, JSEEDR). For use here, JSEEDR
+C is a random variable. The user wants to restore whatever value JSEEDR had
+C at the instant it was zeroed by the preceding $DEPOSIT. The user can not do
+C this manually (old logic) because he does not know what the value was. But
+C ATP remembers the value within DEPOSI, and can do it. WSM. 31 July 2011
+BLANK CARD ENDING PLOT CARDS
+BEGIN NEW DATA CASE
+C 12th of 16 subcases was added 13 April 2011 to illustrate that Orlando
+C Hevia's new EGAUSS function of MODELS will perform Gaussian elimination.
+C This was the 1st of 4 stacked subcases within Orlando's file TFGNMOD.DAT
+C Solve the following 3 linear equations in 3 unknowns:
+C X1 + 4*X2 + X3 = 7
+C X1 + 6*X2 - X3 = 13
+C 2*X1 - X2 + 2*X3 = 5
+C THE INITIAL SOLUTION IS X1= 5, X2= 1, X3= -2, x4=18
+C THE ELEMENT N+1 OF SOLUTION IS THE DETERMINANT OF the coefficient matrix
+ 0.2 1.0 { Take 5 time steps, changing the matrix at each new step
+ 1 -1 { Suppress accumulation of plot points, which are unused
+MODELS STAND ALONE
+C 3456789012345678901234567890123456789012345678 -- Ruler for following optional
+C < File Name> IPRMDL KOMPAC miscellaneous data card :
+MODELS MISC. DATA dc68seek. 0 0
+C THE FUNCTION CAN BE USED AS FOLLOW:
+C
+C SOLUTION[1..N+1]:= EGAUSS(N,MATRIX[1..N*N+N])
+C
+C EGAUSS FROM GAUSSIAN ELIMINATION
+C
+MODEL solve
+FUNCTION egauss FOREIGN egauss {ixarg:100}
+ VAR arg[1..12]
+ VAR x[1..4]
+EXEC
+if t=0 then
+ arg[1]:=1.0
+ arg[2]:=4.0
+ arg[3]:=1.0
+ arg[4]:=1.0
+ arg[5]:=6.0
+ arg[6]:=-1.0
+ arg[7]:=2.0
+ arg[8]:=-1.0
+ arg[9]:=2.0
+ arg[10]:=7.0
+ arg[11]:=13.0
+ arg[12]:=5.0
+endif
+x[1..4]:=egauss(3,arg[1..12])
+write( 'Solution vector X =', x[1], x[2], x[3], ' Determinant =', x[4] )
+arg[1]:= arg[1]*2.0
+arg[12]:= arg[12]+1.0
+C Warning about the following use of RANDOM, which is the MODELS connection
+C to ATP's random number generator. For KOMPAR = 0 or 1 (set in STARTUP),
+C answers will be different every time the data is simulated, and different
+C program versions can not be compared (e.g., GNU vs. Salford). For purposes
+C of a test case, we want repeatability and comparability, so set KOMPAR = 3.
+C Also, position within the data case (here, the 12th subcase) may be critical
+C since future random numbers depend on past random numbers.
+arg[2]:=random()
+C arg[2]:= arg[1] * arg[12] --- a deterministic alternative to preceding line
+write('arg(1,12,2) =', arg[1], arg[12], arg[2] )
+ENDEXEC
+ENDMODEL
+USE solve AS solve
+ENDUSE
+ENDMODELS
+C Because there is no electric network, there is no usual output of the
+C time-step loop. Once the dT loop begins, the only output is from MODELS:
+C Blank card ends electric sources. KCONST = 1 |BLANK card ending dummy source that ATP added to TACS or MODELS STAND ALONE
+C X vector = 5.0 1.0 -2. Determinant = 18.
+C arg(1,12,2) = 2.0 6.0 .2328306E-9
+C Blank card ending requests for output variables. |BLANK card ending output requests (none) that ATP added to this STAND ALONE
+C Column headings for the 0 EMTP output variables follow. These are divided among the 5 possible classes as follows ....
+C X vector = 3.111111111 1.777777778 .7777777784 Determinant = 8.999999999
+C arg(1,12,2) = 4.0 7.0 .1608161E-4
+C X vector = .7096617728 2.741941174 4.161308814 Determinant = 30.99993567
+C arg(1,12,2) = 8.0 8.0 .1107408979
+C X vector = .1640999973 3.031236365 5.351518185 Determinant = 74.55703641
+C arg(1,12,2) = 16. 9.0 .763080108
+C X vector = -.104457932 3.219802884 6.214359374 Determinant = 159.9476796
+C arg(1,12,2) = 32. 10. .1799780347
+C X vector = -.007328865 3.275392315 6.645025023 Determinant = 338.2800879
+C arg(1,12,2) = 64. 11. .902878134
+BLANK card ending plot requests
+BEGIN NEW DATA CASE
+C 13th of 16 subcases has the same solution as the preceding. It was added
+C 30 April 2011 to illustrate the LINSOL alternative to EGAUSS. The data
+C is compatible except that no value for the determinant of the coefficient
+C matrix is available. About history, SUBROUTINE LINSOL was provided by
+C MODELS author Laurent Dube, paid by BPA, as part of the original user-
+C supplied source code of MODELS. This dates to 1995 or 1996 (not part of
+C the original MODELS). Comment cards document the source and use:
+C Solves [a]*x=b using Crout's method with partial pivoting
+C Reference: Numerical Recipes chapter 2.
+C Modifies [a] and b, and returns b=x
+C Returns ierr=1 if [a] is singular, else returns ierr=0
+C Uses divzro as the smallest possible divisor value
+ 0.2 1.0 { Take 5 time steps, changing the matrix at each new step
+ 1 -1 { Suppress accumulation of plot points, which are unused
+MODELS STAND ALONE
+MODEL solve
+FUNCTION linsol FOREIGN linsol {ixarg:100}
+ VAR arg[1..12]
+ VAR x[1..4]
+EXEC
+if t=0 then
+C As written, LINSOL has been connected using the interface of Orlando's
+C EGAUSS. There is a complication, however. EGAUSS stored the matrix by
+C rows whereas LINSOL assumes storage by columns. Note the difference:
+C 1 2 3 4 5 6 7 8 9 -- Index of MODELS
+C 1,1 1,2 1,3 2,1 2,2, 2,3 3,1 3,2 3,3 -- Content for EGAUSS
+C 1,1 2,1 3,1 1,2 2,2, 3,2 1,3 2,3 3,3 -- Content for LINSOL
+C I.e., one is the transpose of the other. Of the 9 elements that must
+C be assigned, only the 3 diagonal elements are unchanged. One at a time:
+ arg[1]:=1.0 -- Element (1,1) of EGAUSS is a diagonal, so remains unchanged
+C arg[2]:=4.0 { Element (1,2) of EGAUSS had index 2, so index 4 of LINSOL:
+ arg[4]:=4.0
+C arg[3]:=1.0 { Element (1,3) of EGAUSS had index 3, so index 7 of LINSOL
+ arg[7]:=1.0
+C arg[4]:=1.0 { Element (2,1) of EGAUSS had index 4, so index 2 of LINSOL
+ arg[2]:=1.0
+ arg[5]:=6.0 -- Element (2,2) of EGAUSS is a diagonal, so remains unchanged
+C arg[6]:=-1.0 { Element (2,3) of EGAUSS had index 6, so index 8 of LINSOL
+ arg[8]:=-1.0
+C arg[7]:=2.0 { Element (3,1) of EGAUSS had index 7, so index 3 of LINSOL
+ arg[3]:=2.0
+C arg[8]:=-1.0 { Element (2,3) of EGAUSS had index 8, so index 6 of LINSOL
+ arg[6]:=-1.0
+ arg[9]:=2.0 -- Element (3,3) of EGAUSS is a diagonal, so remains unchanged
+C That was for the 3x3 square matrix. The right-hand-side vector is unchanged:
+ arg[10]:=7.0
+ arg[11]:=13.0
+ arg[12]:=5.0
+endif
+x[1..4]:=linsol(3,arg[1..12])
+write( 'Solution vector X =', x[1], x[2], x[3] )
+arg[1]:= arg[1]*2.0
+arg[12]:= arg[12]+1.0
+C arg[2]:=random() { Element (1,2) of EGAUSS had index 2, so index 4 of LINSOL
+ arg[4]:=random()
+C write('arg(1,12,2) =', arg[1], arg[12], arg[2] )
+write( 'arg(1,12,4) =', arg[1], arg[12], arg[4] )
+ENDEXEC
+ENDMODEL
+USE solve AS solve
+ENDUSE
+ENDMODELS
+C Because there is no electric network, there is no usual output of the
+C time-step loop. Once the dT loop begins, the only output is from MODELS:
+C Blank card ends electric sources. KCONST = 1 |BLANK card ending dummy source that ATP added to TACS or MODELS STAND ALONE
+C X vector = 5.0 1.0 -2.
+C arg(1,12,4) = 2.0 6.0 .2328306E-9
+C Blank card ending requests for output variables. |BLANK card ending output requests (none) that ATP added to this STAND ALONE
+C Column headings for the 0 EMTP output variables follow. These are divided among the 5 possible classes as follows ....
+C X vector = 3.111111111 1.777777778 .7777777784
+C arg(1,12,4) = 4.0 7.0 .1608161E-4
+C X vector = .7096617728 2.741941174 4.161308814
+C arg(1,12,4) = 8.0 8.0 .1107408979
+C X vector = .1640999973 3.031236365 5.351518185
+C arg(1,12,4) = 16. 9.0 .763080108
+C X vector = -.104457932 3.219802884 6.214359374
+C arg(1,12,4) = 32. 10. .1799780347
+C X vector = -.007328865 3.275392315 6.645025023
+C arg(1,12,4) = 64. 11. .902878134
+BLANK card ending plot requests
+BEGIN NEW DATA CASE
+C 14th of 16 subcases was added 17 April 2011 to illustrate that Orlando
+C Hevia's new INVERT function of MODELS will indeed invert a matrix.
+C Data comes from Orlando Hevia of Santa Fe, Argentina.
+ 1.0 5.0 { Take 5 time steps, changing the matrix at each new step
+ 1 -1 { Suppress accumulation of plot points, which are unused
+MODELS STAND ALONE
+C The input matrix is on the left, its inverse is on the right:
+C ORIGINAL matrix INVERSE CALCULATED WITH QPRO
+C 1 2 3 -0.33333 0.151515 0.121212
+C 4 2 5 0.666667 -0.48485 0.212121
+C 6 3 2 -2E-17 0.272727 -0.18182
+C
+C DEFINE ALL THE VARIABLES TO BE PASSED TO INVERT AS OUTPUT
+C Since the matrix to be inverted is 3x3, 9 scalar variables are used :
+VAR a1, a2, a3, a4, a5, a6, a7, a8, a9
+C Next come 4 scalar dimensions of the call to inversion routine INVERT:
+C IXDATA: NUMBER ELEMENTS IN AN ARRAY XDATA OR ARGUMENTS, FOR EXAMPLE, THE
+C NUMBER OF EQUATIONS, A KEY TO INVERT OR TO CALCULATE DETERMINANT
+C
+C IXIN: NUMBER OF ELEMENTS IN THE INPUT ARRAY XIN. SQUARE OF NUMBER OF
+C EQUATIONS IN THIS CASE
+C
+C IXOUT: NUMBER OF ELEMENTS IN THE OUTPUT ARRAY XOUT. SQUARE OF NUMBER OF
+C EQUATIONS IN THIS CASE, 1 FOR DETERMINANT (NOt IMPLEMENTED)
+C
+C IXVAR: NUMBER OF ELEMENTS IN HISTORY ARRAY, NOt USED IN THIS CASE
+C
+C The following is the call to invert the matrix. IXIN is # of elements (9)
+C of the input matrix, IXOUT is the # of elements (9) of the output matrix.
+C The 1st and the 4th parameters (IXDATA and IXVAR, respectively) do not
+C depend on the matrix order (3 in this case).
+MODEL invert FOREIGN invert {ixdata:1, ixin:9, ixout:9, ixvar:0}
+USE invert AS invert
+ DATA xdata[1]:=3 -- Order of the matrices involved (here, 3)
+C The following defines the 9 cells of the input matrix XIN. Storage is by
+C rows. Note the trailing "t" which adds simulation time to each element of
+C the matrix at each time step :
+ INPUT xin[1..9] := [ 1.0, 2.0, 3.0,
+ 4.0, 2.0, 5.0,
+ 6.0, 3.0, 2.0 ] + t
+C The following connects the 9 elements of the output matrix (the inverse)
+C with variable names A1 through A9. This is done to interface with the
+C ATP output vector (the preceding subcase had no such connection, note) :
+ OUTPUT a1:=xout[1],a2:=xout[2],a3:=xout[3],
+ a4:=xout[4],a5:=xout[5],a6:=xout[6],
+ a7:=xout[7],a8:=xout[8],a9:=xout[9]
+ENDUSE
+RECORD
+ a1 AS a11
+ a2 AS a12
+ a3 AS a13
+ a4 AS a21
+ a5 AS a22
+ a6 AS a23
+ a7 AS a31
+ a8 AS a32
+ a9 as a33
+ENDMODELS
+C Column headings for the 9 EMTP output variables follow. These are divided among the 5 possible classes as follows ....
+C Next 9 output variables belong to MODELS (with "MODELS" an internally-added upper name of pair).
+C Step Time MODELS MODELS MODELS MODELS MODELS MODELS MODELS MODELS MODELS
+C A11 A12 A13 A21 A22 A23 A31 A32 A33
+C 0 0.0 -.33333333 .151515152 .121212121 .666666667 -.48484848 .212121212 0.0 .272727273 -.18181818
+C 1 1.0 -.31914894 .14893617 .127659574 .574468085 -.46808511 .170212766 -.0212766 .276595745 -.19148936
+C 2 2.0 -.31147541 .147540984 .131147541 .524590164 -.45901639 .147540984 -.03278689 .278688525 -.19672131
+C 3 3.0 -.30666667 .146666667 .133333333 .493333333 -.45333333 .133333333 -.04 .28 -.2
+C 4 4.0 -.30337079 .146067416 .134831461 .471910112 -.4494382 .123595506 -.04494382 .280898876 -.20224719
+C 5 5.0 -.30097087 .145631068 .13592233 .45631068 -.44660194 .116504854 -.04854369 .281553398 -.2038835
+BLANK card ending plot requests
+BEGIN NEW DATA CASE
+$PREFIX, [] { $INCLUDE files are located in same place as this main data file
+C Preceding line allows the SEEK file to be remote if this main data file is
+PRINTED NUMBER WIDTH, 12, 2, { Keep dT loop precision the same, but 2 blank separators
+C 15th of 16 subcases was added 21 April 2011 to illustrate extensions to
+C Orlando Hevia's SEEK function that were not tested by the preceding 11th
+C subcase. Data is similar. That is, function # 82 began as a copy of # 32
+C but is more general in that two redundant (X, Y) points have been added
+C (to be omitted by ATP) and the software EOF is used, permitting arbitrary
+C comment information at the bottom of DC68SE82.DAT Erasure of a no longer
+C used SEEK function also is illustrated.
+ 1.0 10.00 { Take ten steps, from 1 second to 10 seconds. This is "T"
+ 1 -1 { Print every step; do not bother saving points for plotting
+MODELS STAND ALONE
+C 3456789012345678901234567890123456789012345678 -- Ruler for following optional
+C < File Name> IPRMDL KOMPAC miscellaneous data card :
+C MODELS MISC. DATA dc68seek. 9 0
+MODEL SAMPLE
+C Of the 6 following variables, the last is a dummy variable that is used for
+C SEEK function erasure only. The returned funct value must be put somewhere.
+VAR X,Y,Z, X37, X49, XDUM
+FUNCTION seek FOREIGN seek {ixarg: 2}
+EXEC
+C Before using function # 82, let's use a second function, the old # 37:
+ IF t<4.0 THEN -- Execute only at the first 3 time steps
+C About the following WRITE( statement, note lack of indentation and
+C use of apostrophes (not quotation marks) to delimit the text. These
+C two details are required to maintain any lower case within the text :
+write( ' Use SEEK model # 37 for T < 4 only.' )
+ X37:=SEEK(T, 37 )
+ ENDIF -- Terminate 4-line IF block
+ IF t=6.0 THEN -- execute only on the 6th time step
+write( ' Erase SEEK function # 37 at T =6.' )
+ XDUM:=SEEK(T, -37 )
+ ENDIF -- Terminate 4-line IF block
+C The following 3 lines use SEEK function # 82, which is the same as # 37:
+X :=SEEK(T, 82 )
+Y:=SQRT(X)*0.3
+Z:=SEEK(SQRT(X)*0.3, 82 )
+C The following nested IF loop results in action only for steps 8 and 9:
+ IF t>7.0 THEN -- Execute only for time steps 8 or later
+ IF t<10.0 THEN -- Execute only for time steps 9 or earlier
+write( ' Use SEEK model # 49 for T = 7 or 8 only.' )
+ X49:=SEEK(T, 49 )
+ ENDIF -- Terminate inner 4-line IF block
+ ENDIF -- Terminate outer 6-line IF block
+C The following IF loop disconnects SEEK function # 49 on step 10:
+ IF t=10.0 THEN -- execute only on the 10th time step
+write( ' Erase SEEK function # 49 at T =10.' )
+ XDUM:=SEEK(T, -49 )
+ ENDIF -- Terminate 4-line IF block
+C In order that all data be seen, let's document the content of DC68SE82.DAT
+C which contains the following formatted lines. Note the 4 consecutive
+C Y values 142 (points 3 through 6 inclusive). Thus points 4 and 5 will be
+C omitted. ATP will store only 7 (X,Y) pairs. Note also the software EOF
+C (value X = -9999.) which allows arbitrary comment information to follow.
+C 0, 69
+C 2, 50
+C 4, 142
+C 4.3, 142,
+C 4.7, 142,
+C 5, 142
+C 6, 188
+C 8, 73
+C 10, 100
+C -9999., 0,
+ENDEXEC
+ENDMODEL
+USE SAMPLE AS SAMPLE
+ENDUSE
+RECORD
+SAMPLE.X AS X
+SAMPLE.Y AS Y
+SAMPLE.Z AS Z
+SAMPLE.X37 AS X37
+SAMPLE.X49 AS X49
+ENDRECORD
+ENDMODELS
+C Use SEEK model # 37 for T < 4 only.
+C MODELS "SEEK" function # 37 involves 7 pairs of (X,Y). These are stored in List-20 SCONST cells 3 through 16.
+C 3 or more consecutive, identical Y has allowed the omission of 2 incoming (X,Y) points. The next line reflects this fact.
+C MODELS "SEEK" function # 82 involves 7 pairs of (X,Y). These are stored in List-20 SCONST cells 19 through 32.
+C Blank card ending requests for output variables. |BLANK card ending output requests (none) that ATP added to this STAND ALONE
+C Column headings for the 5 EMTP output variables follow. These are divided among the 5 possible classes as follows ....
+C Next 5 output variables belong to MODELS (with "MODELS" an internally-added upper name of pair).
+C Step Time MODELS MODELS MODELS MODELS MODELS
+C X Y Z X37 X49
+C 0 0.0 69. 2.49198716 72.6314093 69. 88888.8889
+C Use SEEK model # 37 for T < 4 only.
+C 1 1.0 59.5 2.31408729 64.4480155 59.5 88888.8889
+C Use SEEK model # 37 for T < 4 only.
+C 2 2.0 50. 2.12132034 55.5807358 50. 88888.8889
+C Use SEEK model # 37 for T < 4 only.
+C 3 3.0 96. 2.93938769 93.2118338 96. 88888.8889
+C 4 4.0 142. 3.57491259 122.445979 96. 88888.8889
+C 5 5.0 142. 3.57491259 122.445979 96. 88888.8889
+C Erase SEEK function # 37 at T =6.
+C MODELS erasure of SEEK function # 37 at time T = 6.00000E+00 sec has gained 16 cells of List 20. IFSEM = 16
+C 6 6.0 188. 4.11339276 142. 96. 88888.8889
+C 7 7.0 130.5 3.4270979 115.646503 96. 88888.8889
+C Use SEEK model # 49 for T = 7 or 8 only.
+C MODELS "SEEK" function # 49 involves 4 pairs of (X,Y). These are stored in List-20 SCONST cells 19 through 26.
+C 8 8.0 73. 2.56320112 75.9072517 96. 73.
+C Use SEEK model # 49 for T = 7 or 8 only.
+C 9 9.0 86.5 2.79016129 86.3474191 96. 86.5
+C Erase SEEK function # 49 at T =10.
+C MODELS erasure of SEEK function # 49 at time T = 1.00000E+01 sec has gained 10 cells of List 20. IFSEM = 16
+C 10 10. 100. 3.0 96. 96. 86.5
+BLANK CARD ENDING PLOT CARDS
+BEGIN NEW DATA CASE
+C 16th of 16 subcases was added 29 April 2011 to illustrate high-order
+C comparison between EGAUSS (see subcase 12) and LINSOL (see subcase 13).
+C The testing is ordered by using EGAUSS for order two (two equations in two
+C unknows). Normally the data would consist of 4 cells for the 2x2 matrix
+C [A] followed by 2 cells for the right hand side. Instead, cell 1 of the
+C 6 must be the special reserved key 23456. This is to be followed by five
+C parameters that dedicated code uses to generate the equation set. Matrix
+C order N is in cell 2. The matrix will be symmetrical, so transposition
+C (a difference between data of EGAUSS and data of LINSOL) is not an issue.
+C An offidiagonal in cell (I,J) will be given value -COEFF / SEPAR ** POWER
+C where SEPAR = IABS ( I - J ) is the distance from the diagonal. Each of
+C the diagonals A(I,I) will be the negative of the sum of all offidagonals
+C plus constant DIAG. The right hand side B(I) will equal RHS. So this
+C is the matrix of a set of resistors connecting each node with every other
+C node. Ground is a node, with DIAG being the admittance to ground. RHS
+C is the injected current at each node.
+C NEW LIST SIZES
+C 0 0 0 0 0 0 0 0 0 0
+C 0 0 0 0 0 0 0 0 0 100000
+C 0 0 0 0 0 10000 0 126000 0 0
+C 240000 742
+C About the preceding program dimensions, two are critical for high order use.
+C These are List 20 and List 28. Both are set to the limits of LISTSIZE.BPA
+C A third, List 26, must be twice the order of the equation set being solved.
+C In order to simplify things, the limiting value for this too has been used.
+C If order 31 below is changed to order 310 (this is arg[2]), uncomment the
+C preceding five NEW LIST SIZES data cards by removing "C " from the left.
+C Expect execution to take about 1000 times as long, of course. Note that 310
+C is close to the limit for List 20 equal to 100K since 3.1**2 is close to 10.
+ 0.2 1.0 { Take 5 time steps, changing the matrix at each new step
+ 1 -1 { Suppress accumulation of plot points, which are unused
+MODELS STAND ALONE
+MODEL solve
+FUNCTION egauss FOREIGN egauss {ixarg:100}
+ VAR arg[1..6]
+ VAR x[1..3]
+EXEC
+if t=0 then
+ arg[1]:=23456.0 -- Reserved key to request test of LINSOL and EGAUSS
+C arg[2]:=310.0 -- Matrix order of the equations being solved. This is N
+ arg[2]:= 31.0 -- Matrix order of the equations being solved. This is N
+ arg[3]:=2.0 -- Exponent to apply to separation. This is POWER
+ arg[4]:=1.0 -- Coefficient of matrix elements. This is COEFF
+ arg[5]:=1.0 -- Right hand side of each equation. This is RHS
+ arg[6]:=31.0 -- Constant term of each diagonal. This is DIAG
+endif
+x[1..2]:=egauss(2,arg[1..6]) -- The 1st arg, 2, is required for this test
+ENDEXEC
+ENDMODEL
+USE solve AS solve
+ENDUSE
+ENDMODELS
+BLANK card ending plot requests
+BEGIN NEW DATA CASE
+BLANK