diff --git a/tests/ALGOLXEM/YUSPAR.alg_m b/tests/ALGOLXEM/YUSPAR.alg_m new file mode 100644 index 0000000..d26b132 --- /dev/null +++ b/tests/ALGOLXEM/YUSPAR.alg_m @@ -0,0 +1,335 @@ +% FREE STANDING YU/SPARROW ROUTINE -- 3/5/70 -- P H KIMPEL 00337000 +BEGIN 00339000 +FILE IN YUSPAR (1,10); 00340000 +FILE OUT OUTPUT 18 "YU/SPAR" "RESULTS" (2,17); 00341000 +FORMAT FT ( "DATE: "A5" PROC TIME: "F8.2" I/O TIME: "F8.2" ELAPSED: " 00342000 + F8.2" SECS"), 00343000 + FMTYS (10(R20.12,/)); 00344000 +INTEGER 00350000 + STARTT, 00353000 + NPOINTS, 00356000 + NA, 00357000 + NB, 00358000 + NAB, 00359000 + N1 , 00360000 + N2, 00361000 + N12, 00362000 + KUSTAR, 00363000 + NUSTAR, 00364000 + NSUSTAR; 00365000 +REAL 00366000 + FN, 00367000 + GMDT, 00368000 + D, 00369000 + H, 00370000 + R, 00371000 + MU1, 00372000 + MU2, 00373000 + DPDZ, 00374000 + PI, 00375000 + MU12, 00376000 + ISUM, 00377000 + YSTRT, 00378000 + YSTEP, 00379000 + YEND, 00380000 + XSTEP, 00381000 + XEND, 00382000 + XX, 00383000 + YY, 00384000 + USUM, 00385000 + Q1Q1FULL, 00386000 + Q2Q2FULL, 00387000 + Q1STARQ1FULL, 00387100 + Q2STARQ2FULL, 00387200 + Q1Q1STAR, 00387300 + Q2Q2STAR, 00387400 + X, 00387500 + Y; 00387600 +REAL TAB21, TAB22, Q1Q2; 00388000 + 00389000 + % BEGIN EXECUTION IN OUTER BLOCK. 00390000 + NPOINTS ~ 100; NA ~ NB ~ 17; NAB ~ NA + NB; 00391000 + PI ~ 3.141592653589796; 00403000 + R ~ 1.0; 00405000 +FOR MU12 ~ 1.0,2.0,4.0,7.0,10.0,20.0,40.0,70.0,100.0,200.0,400.0, 00410000 + 700.0,1000.0,2000.0 DO 00411000 + FOR H ~ -0.8 STEP 0.2 UNTIL 0.8 DO 00412000 + BEGIN 00413000 + STARTT ~ TIME(1); 00414000 +%***********************************************************************00418000 +BEGIN % ***** INNER BLOCK ***** 00419000 + REAL ARRAY 00420000 + W[1:NPOINTS], 00421000 + Y[1:1, 1:NPOINTS+NAB], 00422000 + Z[1:NAB, 1:NPOINTS+NAB], 00423000 + X[1:NAB+1, 1:NPOINTS+NAB], 00423100 + DEV[1:1, 1:NPOINTS], 00423200 + COF[1:1, 1:NAB], 00423300 + STD[1:1], 00423400 + CV[1:NAB+1, 1:NAB], 00423500 + VCV[1:1, 1:NAB+1, 1:NAB], 00423600 + Q[1:1, 1:NAB+1], 00423700 + Q2, 00423800 + E, 00423900 + EP[1:1, 1:NAB], 00424000 + A[1:NAB, 1:NAB], 00424100 + GF[1:NAB+1], 00424200 + ENF[1:NAB]; 00424300 +LABEL DUMPIT; 00424500 +DUMP OUTPUT (D, MU1, MU2, H, MU12, N1, N2, GMDT, Q1Q1STAR, Q2Q2STAR, 00424600 + COF, STD, TAB21, TAB22, Q1Q2) DUMPIT:1; 00424700 + 00425000 +%***********************************************************************00426000 +REAL PROCEDURE G (I, X, Y); 00427000 + VALUE I, X, Y; REAL X, Y; INTEGER I; 00428000 + CASE (I - 1) OF 00429000 + BEGIN 00430000 + G ~ 1; % 100431000 + G ~ X; % 200432000 + G ~ Y; % 300433000 + G ~ X*2 - Y*2; % 400434000 + G ~ 2.0 | X | Y; % 500435000 + G ~ X*3 - 3.0 | X | Y*2; % 600436000 + G ~ 3.0 | X*2 | Y - Y*3; % 700437000 + G ~ X*4 + Y*4 - 6.0 | X*2 | Y*2; % 800438000 + G ~ 4.0 | X*3 | Y - 4.0 | X | Y*3; % 900439000 + G ~ X*5 - 10.0 | X*3 | Y*2 + 5.0 | X | Y*4; %1000440000 + G ~ Y*5 - 10.0 | Y*3 | X*2 + 5.0 | Y | X*4; %1100441000 + G ~ X*6 - 15.0 | X*4 | Y*2 + 15.0 | X*2 | Y*4 - Y*6; %1200442000 + G ~ 6.0 | X*5 | Y + 6.0 | X | Y*5 - 20.0 | X*3 | Y*3; %1300443000 + G ~ X*7 - 21.0 |X*5 | Y*2 + 35.0 | X*3 | Y*4 - 7.0 | X | Y*6; %1400444000 + G ~ 7.0 |X*6 | Y - 35.0 | X*4 | Y*3 + 21.0 | X*2 | Y*5 - Y*7; %1500445000 + G ~ X*8 + Y*8 - 28.0 |X*6 | Y*2 + 70.0 | X*4 | Y*4 - 28.0 00446000 + | X*2 | Y*6; %1600447000 + G ~ 8.0 | X*7 |Y - 56.0 | X*5 | Y*3 + 56.0 | X*3 | Y*5 00448000 + - 8.0 |X | Y*7; %1700449000 + END G; 00450000 + 00451000 +%***********************************************************************00452000 +REAL PROCEDURE DG (I, X, Y); 00453000 + VALUE I, X, Y; INTEGER I; REAL X, Y; 00454000 + CASE (I - 1) OF 00455000 + BEGIN 00456000 + DG ~ 0.0; % 100457000 + DG ~ 0.0; % 200458000 + DG ~ 1.0; % 300459000 + DG ~ -2.0 | Y; % 400460000 + DG ~ 2.0 | X; % 500461000 + DG ~ -6.0 | X |Y; % 600462000 + DG ~ 3.0 | X*2 - 3.0 | Y*2; % 700463000 + DG ~ 4.0 | Y*3 - 12.0 |X*2 | Y; % 800464000 + DG ~ 4.0 | X*3 - 12.0 | X | Y*2; % 900465000 + DG ~ -20.0 |X*3 | Y + 20.0 | X | Y*3; %1000466000 + DG ~ 5.0 | Y*4 - 30.0 | Y*2 |X*2 + 5.0 | X*4; %1100467000 + DG ~ -30.0 | X*4 |Y + 60.0 | X*2 | Y*3 - 6.0 | Y*5; %1200468000 + DG ~ 6.0 | X*5 + 30.0 | X | Y*4 - 60.0 | X*3 | Y*2; %1300469000 + DG ~ -42.0 | X*5 | Y + 140.0 | X*3 | Y*3 - 42.0 | X | Y*5; %1400470000 + DG ~ 7.0 | X*6 - 105.0 | X*4 | Y*2 + 105.0 | X*2 | Y*4 00471000 + - 7.0 | Y*6; %1500472000 + DG ~ 8.0 |Y*7 - 56.0 | X*6 | Y + 280.0 | X*4 | Y*3 00473000 + - 168.0 | X*2 | Y*5; %1600474000 + DG ~ 8.0 | X*7 - 168.0 | X*5 | Y*2 + 280.0 | X*3 | Y*4 00475000 + - 56.0 | X | Y*6; %1700476000 + END DG; 00477000 + 00478000 +%***********************************************************************00479000 +REAL PROCEDURE USTAR (I, X, Y); 00480000 + VALUE I, X, Y; INTEGER I; REAL X, Y; 00481000 + BEGIN 00482000 + NUSTAR ~ IF I = 1 THEN NA ELSE NA + NB; 00483000 + NSUSTAR ~ IF I = 1 THEN 1 ELSE NA + 1; 00484000 + USUM ~ 0.0; 00485000 + FOR KUSTAR ~ NSUSTAR STEP 1 UNTIL NUSTAR DO 00486000 + USUM ~ USUM + COF[1,KUSTAR] | G(KUSTAR-NSUSTAR+1, X, Y); 00487000 + USTAR ~ USUM; 00488000 + END USTAR; 00489000 + 00490000 +%***********************************************************************00491000 +REAL PROCEDURE INTEGRALUSTAR (I, S); 00492000 + VALUE I, S; INTEGER I; BOOLEAN S; 00493000 + BEGIN 00494000 + ISUM ~ 0.0; 00495000 + IF I = 1 THEN BEGIN YSTRT ~ H; YEND ~ 0.9999; END 00496000 + ELSE BEGIN YSTRT ~ -0.9999; YEND ~ H; END; 00497000 + YSTEP ~ (YEND - YSTRT)/30.00; 00498000 + FOR YY ~ YSTRT STEP YSTEP UNTIL YEND DO 00499000 + BEGIN XSTEP ~ (XX ~ SQRT(1.0 - YY*2))/15.0; 00500000 + FOR XX ~ -XX STEP XSTEP UNTIL 0.0 DO 00501000 + ISUM ~ ISUM + XSTEP | YSTEP | (USTAR (I, XX, YY) 00502000 + - (IF S THEN 0.25 | (XX*2 + YY*2) ELSE 0.0)); 00503000 + END; 00504000 + INTEGRALUSTAR ~ 2.0 | ISUM; 00505000 + END INTEGRALUSTAR; 00506000 + 00507000 +%***********************************************************************00508000 +REAL PROCEDURE ARCSIN (X); 00509000 + VALUE X; REAL X; 00510000 + ARCSIN ~ ARCTAN(X/SQRT(1.0 - X*2)); 00511000 + 00512000 +%***********************************************************************00513000 +REAL PROCEDURE EQUATION23 (A); 00514000 + VALUE A; REAL A; 00515000 + EQUATION23 ~ 0.5 - 2.0/3.0/PI | (A | (1.0 - A*2)*1.5 00516000 + + 1.5 | A |SQRT(1.0 - A*2) + 1.5 | ARCSIN(A)); 00517000 + 00518000 +%***********************************************************************00519000 + 00520000 +PROCEDURE ORTHO (W, Y, Z, N, M, P, COF); 00521000 +VALUE N, M, P; 00522000 +REAL ARRAY COF, W[1]; 00523000 +REAL ARRAY Y, Z [1, 1]; 00524000 +INTEGER N, M, P; 00525000 +% 00526000 +%COMMENT 00527000 +% -ORTHO- IS TAKEN FROM ACM ALGORITHM 127 [COMM. ACM, VOL.5, 00528000 +% OCTOBER 1962, P. 511, AUTHOR: PHILIP J. WALSH], 00529000 +% AND INCORPORATES THE SUGGESTIONS IN BARRODALE-S 00530000 +% CERTIFICATION [COMM. ACM, VOL.13, #2, FEB.1970, 122] 00531000 +% .... P H KIMPEL 3/25/70 U OF DELAWARE 00534000 +% 00535000 +BEGIN 00536000 +INTEGER NPP, M1, N2, M2, R1, RBAR, P2, BEI, I, J, K, THI, ALI, 00537000 + OMI, R; 00538000 +REAL ARRAY PK, XP [1:N+P], QK[1:M+1], X[1:M+1, 1:N+P]; 00539000 +REAL SUM, DK2, DK; 00540000 +LABEL BOX1, BOX2, BOX3, BOX4, BOX5, BOX6, BOX7, OM1, OM2, BOX8, BOX9, 00541000 + TH1, TH2, TH3, BOX10, AL1, AL2, BOX12, BE1, BE2, BOX14, FINAL; 00542000 +SWITCH BE ~ BE1, BE2; 00543000 +SWITCH TH ~ TH1, TH2, TH3; SWITCH AL ~ AL1, AL2; 00544000 +SWITCH OM ~ OM1, OM2; 00545000 + NPP ~ N + P; R ~ 1; M1 ~ M-1; N2 ~ N+1; M2 ~ M+1; 00546000 + R1 ~ 0; RBAR ~ R; P2 ~ P+1; BEI ~ 1; 00547000 +BOX1: FOR I ~ 1 STEP 1 UNTIL M DO BEGIN 00548000 + FOR J ~ 1 STEP 1 UNTIL N DO 00549000 + X[I,J] ~ Z[I,J] END; 00550000 +BOX2: FOR I ~ 1 STEP 1 UNTIL M DO BEGIN 00551000 + FOR J ~ N2 STEP 1 UNTIL NPP DO 00552000 + X[I,J] ~ 0.0; X[I,N+I] ~ 1.0 END; GO TO BOX3; 00553000 +BOX3: K ~ 1; 00554000 +BOX4: THI ~ 1; 00555000 +BOX5: ALI ~ OMI ~ 1; 00556000 + FOR J ~ 1 STEP 1 UNTIL P DO PK[N+J] ~ 0.0; 00557000 +BOX6: FOR I ~ 1 STEP 1 UNTIL N DO 00558000 + PK[I] ~ X[K,I] | W[I]; 00559000 +BOX7: GO TO OM[OMI]; 00560000 + OM1: FOR I ~ 1 STEP 1 UNTIL K DO BEGIN SUM ~ 0.0; 00561000 + FOR J ~ 1 STEP 1 UNTIL NPP DO 00562000 + SUM ~ SUM + PK[J] | X[I,J]; QK[I] ~ SUM END; 00563000 + GO TO BOX8; 00564000 + OM2: DK2 ~ 0.0; FOR I ~ 1 STEP 1 UNTIL NPP DO 00565000 + DK2 ~ DK2 + PK[I] | X[K,I]; 00566000 + DK ~ SQRT(DK2); 00567000 + FOR I ~ 1 STEP 1 UNTIL NPP DO 00568000 + X[K,I] ~ X[K,I]/DK;; 00569000 + OMI ~ 1; GO TO BOX6; 00570000 +BOX8: FOR I ~ 1 STEP 1 UNTIL K-1 DO 00571000 + QK[I] ~ -QK[I]; QK[K] ~ 1.0; 00572000 + FOR I ~ 1 STEP 1 UNTIL NPP DO BEGIN 00573000 + SUM ~ 0.0; FOR J ~ 1 STEP 1 UNTIL K DO 00574000 + SUM ~ SUM + X[J,I] | QK[J]; 00575000 + XP[I] ~ SUM END; GO TO BOX9; 00576000 +BOX9: GO TO TH[THI]; 00577000 + TH1: FOR I ~ 1 STEP 1 UNTIL NPP DO 00578000 + X[K,I] ~ XP[I]; GO TO BOX10; 00579000 + TH2: FOR I ~ 1 STEP 1 UNTIL P DO 00580000 + COF[I] ~ -XP[N+I]; THI ~ 3; GO TO TH1; 00581000 + TH3: GO TO BOX14; 00582000 +BOX10: GO TO AL[ALI]; 00583000 + AL1: OMI ~ ALI ~ 2; GO TO BOX6; 00584000 + AL2: IF K < M THEN BEGIN K ~ K + 1; GO TO BOX4; END 00585000 + ELSE GO TO BOX12; 00586000 +BOX12: GO TO BE[BEI]; 00587000 + BE1: BEI ~ THI ~ 2; K ~ K + 1; GO TO BOX14; 00588000 + BE2: GO TO BOX14; 00589000 +BOX14: IF RBAR = 0 THEN GO TO FINAL ELSE RBAR ~ RBAR - 1; 00590000 + R1 ~ R1 + 1; THI ~ 2; 00591000 + FOR I ~ 1 STEP 1 UNTIL N DO 00592000 + X[M2,I] ~ Y[R1,I]; 00593000 + FOR I ~ 1 STEP 1 UNTIL P DO 00594000 + X[M2,N+I] ~ 0.0; GO TO BOX5; 00595000 +FINAL: END ORTHO ; 00596000 +%***********************************************************************00597000 + 00598000 + 00598400 +%***********************************************************************00598500 +%***********************************************************************00599000 +% BEGIN EXECUTION (INNER BLOCK) %00600000 +%***********************************************************************00601000 +%***********************************************************************00601500 +BEGIN 00602000 + COMMENT INITALIZE & SET UP DATA FOR -ORTHO-; 00604000 + INTEGER I, J; 00605000 + REAL THETA1, ARCC1, ARCC2, ARCC12, DS1, DS2, DS12, XX, YY, TOTALARC; 00606000 +%% CALCULATE ARC LENGTHS OF EACH INTERFACE. 00607000 + ARCC12 ~ 2.0 |SQRT(1.0 - H*2); 00608000 + ARCC1 ~ 2.0 | (THETA1 ~ PI/2.0 - ARCTAN(2.0|H/ARCC12)); 00609000 + ARCC2 ~ 2.0 |PI - ARCC1; 00610000 +%% CALCULATE NO. OF POINTS ON EACH ARC. 00611000 + N1 ~ NPOINTS | ARCC1/(TOTALARC ~ ARCC1 + ARCC2 + 2|ARCC12); 00612000 + N2 ~ NPOINTS |ARCC2/TOTALARC; 00613000 + IF (N12 ~ NPOINTS - N1 - N2)MOD 2 = 0 THEN 00614000 + N12 ~ N12 DIV 2 00615000 + ELSE BEGIN 00616000 + N12 ~ N12 DIV 2 + 1; 00617000 + N2 ~ N2 - 1; 00618000 + END; 00619000 +%% CALCULATE THE LENGTH OF LINE SEGMENTS BETWEEN POINTS ON EACH 00620000 +%% INTERFACE. 00621000 + DS1 ~ ARCC1/(N1 - 1); 00622000 + DS2 ~ ARCC2/(N2 - 1); 00623000 + DS12 ~ ARCC12/(N12 - 1); 00624000 + 00625000 + FOR J ~ 1 STEP 1 UNTIL NPOINTS-N12 DO 00626000 + IF J { N1 THEN %% DATA FOR ARC C1. 00627000 + BEGIN 00628000 + W[J] ~ 1.0; Y[1,J] ~ 0.25; 00629000 + XX ~ SIN(YY ~ DS1|(J-1)-THETA1); YY ~ COS(YY); 00630000 + FOR I ~ 1 STEP 1 UNTIL NAB DO 00631000 + Z[I,J] ~ IF I { NA THEN G(I, XX, YY) ELSE 0.0; 00632000 + END ELSE 00633000 + IF J { N1 + N2 THEN %% DATA FOR ARC C2. 00634000 + BEGIN 00635000 + W[J] ~ DS2/DS1; Y[1,J] ~ 0.25; 00636000 + XX ~ SIN(YY ~ - THETA1 - DS2|(J-N1-1)); YY ~ COS(YY); 00637000 + FOR I ~ 1 STEP 1 UNTIL NAB DO 00638000 + Z[I,J] ~ IF I { NA THEN 0.0 ELSE G(I-NA, XX, YY); 00639000 + END ELSE %% DATA FOR ARC C12. 00640000 + BEGIN 00641000 + W[J] ~ W[J+N12] ~ DS12/DS1; 00642000 + Y[1,J] ~ 0.25 | ((XX ~ DS12|(J-N1-N2-1)-ARCC12/2)*2 + H*2)00643000 + | (1.0 - MU12); 00644000 + Y[1,J+N12] ~ 0.0; 00645000 + FOR I ~ 1 STEP 1 UNTIL NAB DO 00646000 + IF I { NA THEN 00647000 + BEGIN 00648000 + Z[I,J] ~ G(I, XX, H); 00649000 + Z[I,J+N12] ~ DG(I, XX, H); 00650000 + END ELSE 00651000 + BEGIN 00652000 + Z[I,J] ~ -G (I-NA, XX, H) | MU12; 00653000 + Z[I,J+N12] ~ -DG (I-NA, XX, H); 00654000 + END; 00655000 + END; 00656000 + END INITIALIZATION; 00657000 + 00658000 +%%% THE NUMBER SMASHING COMMENCES. 00659000 + ORTHO (W, Y, Z, NPOINTS, NAB, NAB, COF[1,*]); %######################00660000 + 00661000 +%%% UPON EXIT FROM ORTHO, THE YU/SPARROW COEFFICIENTS, A[J], B[J], 00662000 +%%%% ARE CONTAINED IN THE VECTOR COF[1,*]. 00663000 + 00664000 + Q1Q1FULL ~ 8.0/PI | INTEGRALUSTAR(1,TRUE); 00665000 + Q2Q2FULL ~ 8.0/PI | INTEGRALUSTAR(2,TRUE); 00666000 + Q1STARQ1FULL ~ EQUATION23(H); 00667000 + Q2STARQ2FULL ~ 1.0 - Q1STARQ1FULL; 00668000 + TAB21 ~ Q1Q1FULL | PI/8.0; TAB22 ~ Q2Q2FULL | PI/8.0; 00669000 + Q1Q2 ~ Q1Q1FULL/Q2Q2FULL/MU12; 00670000 + Q1Q1STAR ~ Q1Q1FULL/Q1STARQ1FULL; 00671000 + Q2Q2STAR ~ Q2Q2FULL/Q2STARQ2FULL; 00671500 + 00672000 + 00673000 +DUMPIT: 00680000 +WRITE(OUTPUT, FT, TIME(0), TIME(2)/60, TIME(3)/60, (TIME(1)-STARTT)/60);00681000 +END INNER BLOCK; 00682000 + END H FOR LOOP; 00684000 +END. 00686000 diff --git a/tests/E-mode/ORTHO-XEM.alg_m b/tests/E-mode/ORTHO-XEM.alg_m new file mode 100644 index 0000000..ee31999 --- /dev/null +++ b/tests/E-mode/ORTHO-XEM.alg_m @@ -0,0 +1,257 @@ +% -ORTHO- TEST, 3/3/70 00000100120429PK +BEGIN 00000200 +INTEGER STARTT; 00000300120429PK +FORMAT FT (//"DATE: "A6" PROC TIME: "F8.2" I/O TIME: "F8.2" ELAPSED: " 00000400120429PK + F8.2" SEC"); 00000500120429PK +INTEGER I, J, N, M, P, R, AI, AUI, ZEI, MUI; 00000600120429PK +REAL SUM, GMDT, FN; 00000700120429PK +FILE PARAMS (KIND=READER, DEPENDENTSPECS); 00000800120826PK +FILE FLUSH (KIND=PRINTER, MAXRECSIZE=132, FRAMESIZE=8, FILEUSE=OUT); 00000900120826PK + 00001000120429PK +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%00002000120429PK +PROCEDURE ORTHO (W, Y, Z, N, FN, M, P, R, AI, AUI, MUI, ZEI, X, DEV, 00003000 + COF, STD, CV, VCV, GMDT, Q, Q2, E, EP, A, GF, ENF); 00004000 +VALUE N, M, P, R, AI, AUI, MUI, ZEI; 00005000 +REAL FN, GMDT; 00006000 +REAL ARRAY STD, GF, W, ENF [1]; 00007000 +REAL ARRAY Y, Z, X, DEV, COF, CV, Q, Q2, E, EP, A [1, 1]; 00008000 +REAL ARRAY VCV [1, 1, 1]; 00009000 +INTEGER N, M, P, R, AI, AUI, ZEI, MUI; 00010000 + 00011000120429PK +COMMENT 00012000 + ORTHO IS TAKEN FROM ACM ALGORITHM 127 [COMM. ACM, VOL.5, 00013000 + OCTOBER 1962, P. 511, AUTHOR: PHILIP J. WALSH]; 00014000 + 00015000120429PK +BEGIN 00016000 +INTEGER NPP, NPM, M1, N2, M2, R1, RBAR, P2, BEI, RHI, I18, GAI, SII, I,00017000 + J, DEI, NUI, E1Z2, E1Z1, K, THI, ALI, OMI, NII; 00018000 +REAL ARRAY PK, XP [1:N+P], QK[1:M+1]; 00019000 +REAL DENOM, SUM, DK2, DK, FI, SS, SSQ; 00020000 +LABEL BOX1, AT1, AT2, BOX2, AU1, AU2, BOX3, BOX4, BOX5, BOX6, MU1, 00021000 + MU2, BOX7, OM1, OM2, BOX8, DE1, BOX8A, BOX8B, DE2, BOX8C, BOX8D, 00022000 + BOX9, TH1, TH2, TH3, BOX10, AL1, AL2, BOX11, NU1, NU2, BOX12, 00023000 + BE1, BE2, BOX13, GA1, GA2, BOX14, RH1, ZE1, ZE2, RH2, SI1, SI2, 00024000 + FINAL; 00025000 +SWITCH AT:= AT1, AT2; SWITCH ZE:= ZE1, ZE2; 00026000120826PK +SWITCH AU:= AU1, AU2; SWITCH MU:= MU1, MU2; 00027000120826PK +SWITCH BE:= BE1, BE2; SWITCH RH:= RH1, RH2; SWITCH GA:= GA1, GA2; 00028000120826PK +SWITCH SI:= SI1, SI2; SWITCH DE:= DE1, DE2; SWITCH NU:= NU1, NU2; 00029000120826PK +SWITCH TH:= TH1, TH2, TH3; SWITCH AL:= AL1, AL2; 00030000120826PK +SWITCH OM:= OM1, OM2; 00031000120826PK + NPP:= N+P; NPM:= N+M; M1:= M-1; N2:= N+1; M2:= M+1; 00032000120826PK + R1:= 0; RBAR:= R; P2:= P+1; DENOM:= IF N=M THEN 1.0 00033000120826PK + ELSE SQRT(N-M); BEI:= RHI:= I18:= 1; 00034000120826PK + IF (P ^= 0) THEN GAI:= SII:= 2 ELSE GAI:= SII:= 1; 00035000120826PK +BOX1: GO TO AT[AI]; 00036000 + AT1: FOR J:= 1 STEP 1 UNTIL N DO BEGIN 00037000120826PK + X[2,J]:= Z[1,J]; X[1,J]:= 1.0 END; 00038000120826PK + FOR I:= 2 STEP 1 UNTIL M1 DO BEGIN 00039000120826PK + FOR J:= 1 STEP 1 UNTIL N DO 00040000120826PK + X[I+1,J]:= X[I,J] * X[2,J] END; GO TO BOX2; 00041000120826PK + AT2: FOR I:= 1 STEP 1 UNTIL M DO BEGIN 00042000120826PK + FOR J:= 1 STEP 1 UNTIL N DO 00043000120826PK + X[I,J]:= Z[I,J] END; 00044000120826PK +BOX2: IF P = 0 THEN GO TO BOX3 ELSE GO TO AU[AUI]; 00045000 + AU1: FOR I:= 1 STEP 1 UNTIL M DO BEGIN 00046000120826PK + FOR J:= N2 STEP 1 UNTIL NPP DO 00047000120826PK + X[I,J]:= 0.0; X[I,N+I]:= 1.0 END; GO TO BOX3; 00048000120826PK + AU2: FOR I:= 1 STEP 1 UNTIL M DO BEGIN 00049000120826PK + FOR J:= N2 STEP 1 UNTIL NPP DO 00050000120826PK + X[I,J]:= Z[I,J] END; 00051000120826PK +BOX3: DEI:= NUI:= E1Z1:= E1Z2:= K:= 1; 00052000120826PK +BOX4: THI:= 1; 00053000120826PK +BOX5: ALI:= OMI:= 1; IF P = 0 THEN GO TO BOX6 ELSE 00054000120826PK + FOR J:= 1 STEP 1 UNTIL P DO PK[N+J]:= 0.0; 00055000120826PK +BOX6: GO TO MU[MUI]; 00056000 + MU1: FOR I:= 1 STEP 1 UNTIL N DO PK[I]:= X[K,I]; 00057000120826PK + GO TO BOX7; 00058000 + MU2: FOR I:= 1 STEP 1 UNTIL N DO 00059000120826PK + PK[I]:= X[K,I] * W[I]; GO TO BOX7; 00060000120826PK +BOX7: GO TO OM[OMI]; 00061000 + OM1: FOR I:= 1 STEP 1 UNTIL K DO BEGIN SUM:= 0.0; 00062000120826PK + FOR J:= 1 STEP 1 UNTIL NPP DO 00063000120826PK + SUM:= SUM + PK[J] * X[I,J]; QK[I]:= SUM END; 00064000120826PK + GO TO BOX8; 00065000 + OM2: DK2:= 0.0; FOR I:= 1 STEP 1 UNTIL NPP DO 00066000120826PK + DK2:= DK2 + PK[I] * X[K,I]; 00067000120826PK + DK:= SQRT(DK2); 00068000120826PK + GF[I18]:= DK; I18:= I18 + 1; 00069000120826PK + FOR I:= 1 STEP 1 UNTIL NPP DO 00070000120826PK + X[K,I]:= X[K,I]/DK;; 00071000120826PK + OMI:= 1; GO TO BOX6; 00072000120826PK +BOX8: GO TO DE[DEI]; 00073000 + DE1: E1Z1:= -E1Z1; IF E1Z1 < 0 THEN GO TO BOX8B ELSE 00074000120826PK + GO TO BOX8A; 00075000 +BOX8A: FOR I:= 1 STEP 1 UNTIL K-1 DO 00076000120826PK + QK[I]:= -QK[I]; QK[K]:= 1.0; 00077000120826PK + FOR I:= 1 STEP 1 UNTIL NPP DO BEGIN 00078000120826PK + SUM:= 0.0; FOR J:= 1 STEP 1 UNTIL K DO 00079000120826PK + SUM:= SUM + X[J,I] * QK[J]; 00080000120826PK + XP[I]:= SUM END; GO TO BOX9; 00081000120826PK +BOX8B: ENF[I18]:= SQRT(QK[K]); GO TO BOX8A; 00082000120826PK + DE2: E1Z2:= -E1Z2; IF E1Z2 < 0 THEN GO TO BOX8C ELSE 00083000120826PK + GO TO BOX8A; 00084000 +BOX8C: FOR I:= 1 STEP 1 UNTIL M DO BEGIN 00085000120826PK + Q[R1,I]:= QK[I]; Q2[R1,I]:= QK[I] * QK[I] END; 00086000120826PK + Q[R1,M2]:= QK[M2]; E[R1,1]:= Q[R1,M2] - Q2[R1,1]; 00087000120826PK + FOR J:= 2 STEP 1 UNTIL M DO 00088000120826PK + E[R1,J]:= E[R1,J-1] - Q2[R1,J]; 00089000120826PK + FI:= 1.0; 00090000120826PK + FOR I:= 1 STEP 1 UNTIL M DO BEGIN 00091000120826PK + IF (FN - FI) > 0.0 THEN BEGIN IF E[R1,I] < 0.0 THEN 00093000 + BEGIN EP[R1,I]:= -SQRT(ABS(E[R1,I])/(FN - FI)); 00094000120826PK + GO TO BOX8D; END 00095000 + ELSE EP[R1,I]:= SQRT(E[R1,I]/(FN - FI)); 00096000120826PK + GO TO BOX8D; END ELSE E[R1,I]:= -1.0; 00097000120826PK +BOX8D: FI:= FI + 1.0; END; GO TO BOX8A; 00098000120826PK +BOX9: GO TO TH[THI]; 00099000 + TH1: FOR I:= 1 STEP 1 UNTIL NPP DO 00100000120826PK + X[K,I]:= XP[I]; GO TO BOX10; 00101000120826PK + TH2: FOR I:= 1 STEP 1 UNTIL N DO 00102000120826PK + DEV[R1,I]:= XP[I]; 00103000120826PK + FOR I:= 1 STEP 1 UNTIL P DO 00104000120826PK + COF[R1,I]:= -XP[N+I]; THI:= 3; GO TO TH1; 00105000120826PK + TH3: GO TO BOX11; 00106000 +BOX10: GO TO AL[ALI]; 00107000 + AL1: OMI:= ALI:= 2; GO TO BOX6; 00108000120826PK + AL2: IF K < M THEN BEGIN K:= K + 1; GO TO BOX4; END 00109000120826PK + ELSE GO TO BOX12; 00110000 +BOX11: GO TO NU[NUI]; 00111000 + NU1: NUI:= 2; GO TO BOX14; 00112000120826PK + NU2: SS:= DK/DENOM; SSQ:= SS * SS; 00113000120826PK + STD[R1]:= SS; GO TO BOX14; 00114000120826PK +BOX12: GO TO BE[BEI]; 00115000 + BE1: FOR I:= 1 STEP 1 UNTIL M DO BEGIN 00116000120826PK + FOR J:= 1 STEP 1 UNTIL P DO 00117000120826PK + A[I,J]:= X[I,N+J] END; 00118000120826PK + GMDT:= 1.0; FOR I:= 1 STEP 1 UNTIL M DO 00119000120826PK + GMDT:= GMDT * (GF[I]/ENF[I]); 00120000120826PK + GMDT:= GMDT * GMDT; DEI:= BEI:= THI:= 2; 00121000120826PK + K:= K + 1; GO TO BOX13; 00122000120826PK + BE2: GO TO BOX11; 00123000 +BOX13: GO TO GA[GAI]; 00124000 + GA1: GO TO BOX11; 00125000 + GA2: FOR I:= 1 STEP 1 UNTIL P DO BEGIN 00126000120826PK + FOR J:= I STEP 1 UNTIL P DO BEGIN 00127000120826PK + SUM:= 0.0; 00128000120826PK + FOR NII:= 1 STEP 1 UNTIL M DO 00129000120826PK + SUM:= SUM + A[NII,I] * A[NII,J]; 00130000120826PK + CV[I,J]:= SUM END END; 00131000120826PK + FOR I:= 1 STEP 1 UNTIL P DO 00132000120826PK + CV[P2,I]:= SQRT(CV[I,I]); GAI:= 1; GO TO BOX11; 00133000120826PK +BOX14: GO TO RH[RHI]; 00134000 + RH1: IF RBAR = 0 THEN GO TO FINAL ELSE RBAR:= RBAR - 1; 00135000120826PK + R1:= R1 + 1; THI:= RHI:= 2; GO TO ZE[ZEI]; 00136000120826PK + ZE1: FOR I:= 1 STEP 1 UNTIL N DO 00137000120826PK + X[M2,I]:= Y[R1,I]; 00138000120826PK + FOR I:= 1 STEP 1 UNTIL P DO 00139000120826PK + X[M2,N+I]:= 0.0; GO TO BOX5; 00140000120826PK + ZE2: FOR I:= 1 STEP 1 UNTIL NPP DO 00141000120826PK + X[M2,I]:= Y[R1,I]; GO TO BOX5; 00142000120826PK + RH2: GO TO SI[SII]; 00143000 + SI1: GO TO RH1; 00144000 + SI2: FOR I:= 1 STEP 1 UNTIL P DO BEGIN 00145000120826PK + FOR J:= I STEP 1 UNTIL P DO 00146000120826PK + VCV[R1,I,J]:= SSQ * CV[I,J] END; 00147000120826PK + FOR I:= 1 STEP 1 UNTIL P DO 00148000120826PK + VCV[R1, P2, I]:= SS * CV[P2,I]; GO TO RH1; 00149000120826PK +FINAL: END ORTHO ; 00150000 +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%00151000120429PK + 00152000120429PK +PROCEDURE DUMP1 (CAPTION, A); 00800000120826PK + VALUE CAPTION; 00800100120826PK + ALPHA CAPTION; 00800200120826PK + ARRAY A[0]; 00800300120826PK + BEGIN 00800400120826PK + INTEGER 00800500120826PK + UB1, 00800600120826PK + X; 00800700120826PK + 00800800120826PK + UB1:= SIZE(A)-1; 00800900120826PK + WRITE (FLUSH, , CAPTION); 00801000120826PK + WRITE (FLUSH, <6E19.11>, 00801100120826PK + FOR X:= 0 STEP 1 UNTIL UB1 DO A[X]); 00801200120826PK + WRITE (FLUSH); 00801300120826PK + END DUMP1; 00801400120826PK + 00810000120826PK +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%00810100120826PK +PROCEDURE DUMP2 (CAPTION, A); 00810200120826PK + VALUE CAPTION; 00810300120826PK + ALPHA CAPTION; 00810400120826PK + ARRAY A[0,0]; 00810500120826PK + BEGIN 00810600120826PK + INTEGER 00810700120826PK + UB1, 00810800120826PK + UB2, 00810900120826PK + X, 00811000120826PK + Y; 00811100120826PK + 00811200120826PK + UB1:= SIZE(A)-1; 00811300120826PK + WRITE (FLUSH, , CAPTION); 00811400120826PK + FOR X:= 0 STEP 1 UNTIL UB1 DO 00811500120826PK + BEGIN 00811600120826PK + UB2:= SIZE(A[X,*])-1; 00811700120826PK + WRITE (FLUSH, <6E19.11>, 00811800120826PK + FOR Y:= 0 STEP 1 UNTIL UB2 DO A[X,Y]); 00811900120826PK + END FOR X; 00812000120826PK + WRITE (FLUSH); 00812100120826PK + END DUMP2; 00812200120826PK + 00900000120826PK +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%00900100120826PK +STARTT:= TIME(1); 00901000120826PK + %%READ (PARAMS, /, N, M, P); 00902000120826PK + %%CLOSE(PARAMS); 00903000120826PK + N:= 7; M:= 5; P:= 5; % PARAMETER VALUES ORIGINALLY READ IN 00903100120826PK + R:= 1; AI:= 2; AUI:= 1; 00904000120826PK + ZEI:= 1; MUI:= 2; FN:= N; 00905000120826PK + 00906000120429PK +BEGIN %INNER BLOCK 00907000120429PK +LABEL DUMPIT; 00908000120429PK +REAL ARRAY W[1:N], Y[1:R, 1:N+P], Z[1:M, 1:N+P], 00909000120429PK + X[1:M+1, 1:N+P], DEV[1:R, 1:N], COF[1:R, 1:P], STD[1:R], 00910000120429PK + CV[1:P+1, 1:P], VCV[1:R, 1:P+1, 1:P], Q[1:R, 1:M+1], 00911000120429PK + Q2, E, EP [1:R, 1:M], A[1:M, 1:P], GF[1:M+R], ENF[1:M]; 00912000120429PK +%%DUMP FLUSH (N, M, P, R, AI, AUI, ZEI, MUI, 00913000120826PK +%% GMDT, FN, 00914000120826PK +%% W, Y, Z, X, DEV, COF, STD, CV, VCV, Q, Q2, A, GF, ENF, E, 00915000120826PK +%% EP) DUMPIT:1; 00916000120826PK + FOR I:= 1 STEP 1 UNTIL N DO 00917000120826PK + BEGIN 00918000120429PK + W[I]:= 1.0; SUM:= 0.0; 00919000120826PK + FOR J:= 1 STEP 1 UNTIL M DO 00920000120826PK + SUM:= SUM + J*(Z[J,I]:= I*J); 00921000120826PK + Y[1,I]:= SUM; 00922000120826PK + END; 00922500120430PK + ORTHO (W, Y, Z, N, FN, M, P, R, AI, AUI, MUI, ZEI, X, DEV, 00923000120429PK + COF, STD, CV, VCV, GMDT, Q, Q2, E, EP, A, GF, ENF); 00924000120429PK + 00924500120430PK + WRITE(FLUSH, <"N = ",I13,/, 00924510120826PK + "M = ",I13,/, 00924520120826PK + "P = ",I13,/, 00924530120826PK + "R = ",I13,/, 00924540120826PK + "AI = ",I13,/, 00924550120826PK + "AUI = ",I13,/, 00924560120826PK + "ZEI = ",I13,/, 00924570120826PK + "MUI = ",I13,/, 00924580120826PK + "GMDT = ",E19.11,/, 00924590120826PK + "FN = ",E19.11,/>, 00924600120826PK + N, M, P, R, AI, AUI, ZEI, MUI, GMDT, FN); 00924610120826PK + DUMP1("W ", W); 00924620120826PK + DUMP2("Y ", Y); 00924630120826PK + DUMP2("Z ", Z); 00924640120826PK + DUMP2("X ", X); 00924650120826PK + DUMP2("DEV ", DEV); 00924660120826PK + DUMP2("COF ", COF); 00924670120826PK + DUMP1("STD ", STD); 00924680120826PK + DUMP2("CV ", CV); 00924690120826PK + DUMP2("VCV ", VCV[1,*,*]); 00924700120826PK + DUMP2("Q ", Q); 00924710120826PK + DUMP2("Q2 ", Q2); 00924720120826PK + DUMP2("A ", A); 00924730120826PK + DUMP1("GF ", GF); 00924740120826PK + DUMP1("ENF ", ENF); 00924750120826PK + DUMP2("E ", E); 00924760120826PK + DUMP2("EP ", EP); 00924770120826PK +DUMPIT: END INNER BLOCK; 00925000120429PK +WRITE(FLUSH, FT, TIME(0), TIME(2)/60, TIME(3)/60, (TIME(1)-STARTT)/60); 00926000120429PK +END. 00999000 diff --git a/tests/E-mode/YUSPAR-XEM.alg_m b/tests/E-mode/YUSPAR-XEM.alg_m new file mode 100644 index 0000000..f1b8202 --- /dev/null +++ b/tests/E-mode/YUSPAR-XEM.alg_m @@ -0,0 +1,338 @@ +% FREE STANDING YU/SPARROW ROUTINE -- 3/5/70 -- P H KIMPEL 00337000 +BEGIN 00339000 +FILE YUSPAR (KIND=READER, DEPENDENTSPECS); 00340000 +FILE OUTPUT (KIND=PRINTER, MAXRECSIZE=132, FRAMESIZE=8, FILEUSE=OUT); 00341000 +FORMAT FT ( "DATE: "A5" PROC TIME: "F8.2" I/O TIME: "F8.2" ELAPSED: " 00342000 + F8.2" SECS"), 00343000 + FMTYS (10(R20.12,/)); 00344000 +INTEGER 00350000 + STARTT, 00353000 + NPOINTS, 00356000 + NA, 00357000 + NB, 00358000 + NAB, 00359000 + N1 , 00360000 + N2, 00361000 + N12, 00362000 + KUSTAR, 00363000 + NUSTAR, 00364000 + NSUSTAR; 00365000 +REAL 00366000 + FN, 00367000 + GMDT, 00368000 + D, 00369000 + H, 00370000 + R, 00371000 + MU1, 00372000 + MU2, 00373000 + DPDZ, 00374000 + PI, 00375000 + MU12, 00376000 + ISUM, 00377000 + YSTRT, 00378000 + YSTEP, 00379000 + YEND, 00380000 + XSTEP, 00381000 + XEND, 00382000 + XX, 00383000 + YY, 00384000 + USUM, 00385000 + Q1Q1FULL, 00386000 + Q2Q2FULL, 00387000 + Q1STARQ1FULL, 00387100 + Q2STARQ2FULL, 00387200 + Q1Q1STAR, 00387300 + Q2Q2STAR, 00387400 + X, 00387500 + Y; 00387600 +REAL TAB21, TAB22, Q1Q2; 00388000 + 00389000 + % BEGIN EXECUTION IN OUTER BLOCK. 00390000 + NPOINTS:= 100; NA:= NB:= 17; NAB:= NA + NB; 00391000 + PI:= 3.141592653589796; 00403000 + R:= 1.0; 00405000 +FOR MU12:= 1.0,2.0,4.0,7.0,10.0,20.0,40.0,70.0,100.0,200.0,400.0, 00410000 + 700.0,1000.0,2000.0 DO 00411000 + FOR H:= -0.8 STEP 0.2 UNTIL 0.8 DO 00412000 + BEGIN 00413000 + STARTT:= TIME(1); 00414000 +%***********************************************************************00418000 +BEGIN % ***** INNER BLOCK ***** 00419000 + REAL ARRAY 00420000 + W[1:NPOINTS], 00421000 + Y[1:1, 1:NPOINTS+NAB], 00422000 + Z[1:NAB, 1:NPOINTS+NAB], 00423000 + X[1:NAB+1, 1:NPOINTS+NAB], 00423100 + DEV[1:1, 1:NPOINTS], 00423200 + COF[1:1, 1:NAB], 00423300 + STD[1:1], 00423400 + CV[1:NAB+1, 1:NAB], 00423500 + VCV[1:1, 1:NAB+1, 1:NAB], 00423600 + Q[1:1, 1:NAB+1], 00423700 + Q2, 00423800 + E, 00423900 + EP[1:1, 1:NAB], 00424000 + A[1:NAB, 1:NAB], 00424100 + GF[1:NAB+1], 00424200 + ENF[1:NAB]; 00424300 +ARRAY COF1[1] = COF[1,*]; 00424400 +LABEL DUMPIT; 00424500 +DUMP OUTPUT (D, MU1, MU2, H, MU12, N1, N2, GMDT, Q1Q1STAR, Q2Q2STAR, 00424600 + COF1, STD, TAB21, TAB22, Q1Q2) DUMPIT:1; 00424700 + 00425000 +%***********************************************************************00426000 +REAL PROCEDURE G (I, X, Y); 00427000 + VALUE I, X, Y; REAL X, Y; INTEGER I; 00428000 + CASE (I - 1) OF 00429000 + BEGIN 00430000 + G:= 1; % 100431000 + G:= X; % 200432000 + G:= Y; % 300433000 + G:= X**2 - Y**2; % 400434000 + G:= 2.0 * X * Y; % 500435000 + G:= X**3 - 3.0 * X * Y**2; % 600436000 + G:= 3.0 * X**2 * Y - Y**3; % 700437000 + G:= X**4 + Y**4 - 6.0 * X**2 * Y**2; % 800438000 + G:= 4.0 * X**3 * Y - 4.0 * X * Y**3; % 900439000 + G:= X**5 - 10.0 * X**3 * Y**2 + 5.0 * X * Y**4; %1000440000 + G:= Y**5 - 10.0 * Y**3 * X**2 + 5.0 * Y * X**4; %1100441000 + G:= X**6 - 15.0 * X**4 * Y**2 + 15.0 * X**2 * Y**4 - Y**6; %1200442000 + G:= 6.0 * X**5 * Y + 6.0 * X * Y**5 - 20.0 * X**3 * Y**3; %1300443000 + G:= X**7 - 21.0 *X**5 * Y**2 + 35.0 * X**3 * Y**4 - 7.0 * X * %1400444000 + Y**6; 00444100 + G:= 7.0 *X**6 * Y - 35.0 * X**4 * Y**3 + 21.0 * X**2 * Y**5 - %1500445000 + Y**7; 00445100 + G:= X**8 + Y**8 - 28.0 *X**6 * Y**2 + 70.0 * X**4 * Y**4 - 28.0 00446000 + * X**2 * Y**6; %1600447000 + G:= 8.0 * X**7 *Y - 56.0 * X**5 * Y**3 + 56.0 * X**3 * Y**5 00448000 + - 8.0 *X * Y**7; %1700449000 + END G; 00450000 + 00451000 +%***********************************************************************00452000 +REAL PROCEDURE DG (I, X, Y); 00453000 + VALUE I, X, Y; INTEGER I; REAL X, Y; 00454000 + CASE (I - 1) OF 00455000 + BEGIN 00456000 + DG:= 0.0; % 100457000 + DG:= 0.0; % 200458000 + DG:= 1.0; % 300459000 + DG:= -2.0 * Y; % 400460000 + DG:= 2.0 * X; % 500461000 + DG:= -6.0 * X *Y; % 600462000 + DG:= 3.0 * X**2 - 3.0 * Y**2; % 700463000 + DG:= 4.0 * Y**3 - 12.0 *X**2 * Y; % 800464000 + DG:= 4.0 * X**3 - 12.0 * X * Y**2; % 900465000 + DG:= -20.0 *X**3 * Y + 20.0 * X * Y**3; %1000466000 + DG:= 5.0 * Y**4 - 30.0 * Y**2 *X**2 + 5.0 * X**4; %1100467000 + DG:= -30.0 * X**4 *Y + 60.0 * X**2 * Y**3 - 6.0 * Y**5; %1200468000 + DG:= 6.0 * X**5 + 30.0 * X * Y**4 - 60.0 * X**3 * Y**2; %1300469000 + DG:= -42.0 * X**5 * Y + 140.0 * X**3 * Y**3 - 42.0 * X * Y**5;%1400470000 + DG:= 7.0 * X**6 - 105.0 * X**4 * Y**2 + 105.0 * X**2 * Y**4 00471000 + - 7.0 * Y**6; %1500472000 + DG:= 8.0 *Y**7 - 56.0 * X**6 * Y + 280.0 * X**4 * Y**3 00473000 + - 168.0 * X**2 * Y**5; %1600474000 + DG:= 8.0 * X**7 - 168.0 * X**5 * Y**2 + 280.0 * X**3 * Y**4 00475000 + - 56.0 * X * Y**6; %1700476000 + END DG; 00477000 + 00478000 +%***********************************************************************00479000 +REAL PROCEDURE USTAR (I, X, Y); 00480000 + VALUE I, X, Y; INTEGER I; REAL X, Y; 00481000 + BEGIN 00482000 + NUSTAR:= IF I = 1 THEN NA ELSE NA + NB; 00483000 + NSUSTAR:= IF I = 1 THEN 1 ELSE NA + 1; 00484000 + USUM:= 0.0; 00485000 + FOR KUSTAR:= NSUSTAR STEP 1 UNTIL NUSTAR DO 00486000 + USUM:= USUM + COF[1,KUSTAR] * G(KUSTAR-NSUSTAR+1, X, Y); 00487000 + USTAR:= USUM; 00488000 + END USTAR; 00489000 + 00490000 +%***********************************************************************00491000 +REAL PROCEDURE INTEGRALUSTAR (I, S); 00492000 + VALUE I, S; INTEGER I; BOOLEAN S; 00493000 + BEGIN 00494000 + ISUM:= 0.0; 00495000 + IF I = 1 THEN BEGIN YSTRT:= H; YEND:= 0.9999; END 00496000 + ELSE BEGIN YSTRT:= -0.9999; YEND:= H; END; 00497000 + YSTEP:= (YEND - YSTRT)/30.00; 00498000 + FOR YY:= YSTRT STEP YSTEP UNTIL YEND DO 00499000 + BEGIN XSTEP:= (XX:= SQRT(1.0 - YY**2))/15.0; 00500000 + FOR XX:= -XX STEP XSTEP UNTIL 0.0 DO 00501000 + ISUM:= ISUM + XSTEP * YSTEP * (USTAR (I, XX, YY) 00502000 + - (IF S THEN 0.25 * (XX**2 + YY**2) ELSE 0.0)); 00503000 + END; 00504000 + INTEGRALUSTAR:= 2.0 * ISUM; 00505000 + END INTEGRALUSTAR; 00506000 + 00507000 +%***********************************************************************00508000 +REAL PROCEDURE ARCSIN (X); 00509000 + VALUE X; REAL X; 00510000 + ARCSIN:= ARCTAN(X/SQRT(1.0 - X**2)); 00511000 + 00512000 +%***********************************************************************00513000 +REAL PROCEDURE EQUATION23 (A); 00514000 + VALUE A; REAL A; 00515000 + EQUATION23:= 0.5 - 2.0/3.0/PI * (A * (1.0 - A**2)**1.5 00516000 + + 1.5 * A *SQRT(1.0 - A**2) + 1.5 * ARCSIN(A)); 00517000 + 00518000 +%***********************************************************************00519000 + 00520000 +PROCEDURE ORTHO (W, Y, Z, N, M, P, COF); 00521000 +VALUE N, M, P; 00522000 +REAL ARRAY COF, W[1]; 00523000 +REAL ARRAY Y, Z [1, 1]; 00524000 +INTEGER N, M, P; 00525000 +% 00526000 +%COMMENT 00527000 +% -ORTHO- IS TAKEN FROM ACM ALGORITHM 127 [COMM. ACM, VOL.5, 00528000 +% OCTOBER 1962, P. 511, AUTHOR: PHILIP J. WALSH], 00529000 +% AND INCORPORATES THE SUGGESTIONS IN BARRODALE-S 00530000 +% CERTIFICATION [COMM. ACM, VOL.13, #2, FEB.1970, 122] 00531000 +% .... P H KIMPEL 3/25/70 U OF DELAWARE 00534000 +% 00535000 +BEGIN 00536000 +INTEGER NPP, M1, N2, M2, R1, RBAR, P2, BEI, I, J, K, THI, ALI, 00537000 + OMI, R; 00538000 +REAL ARRAY PK, XP [1:N+P], QK[1:M+1], X[1:M+1, 1:N+P]; 00539000 +REAL SUM, DK2, DK; 00540000 +LABEL BOX1, BOX2, BOX3, BOX4, BOX5, BOX6, BOX7, OM1, OM2, BOX8, BOX9, 00541000 + TH1, TH2, TH3, BOX10, AL1, AL2, BOX12, BE1, BE2, BOX14, FINAL; 00542000 +SWITCH BE:= BE1, BE2; 00543000 +SWITCH TH:= TH1, TH2, TH3; SWITCH AL:= AL1, AL2; 00544000 +SWITCH OM:= OM1, OM2; 00545000 + NPP:= N + P; R:= 1; M1:= M-1; N2:= N+1; M2:= M+1; 00546000 + R1:= 0; RBAR:= R; P2:= P+1; BEI:= 1; 00547000 +BOX1: FOR I:= 1 STEP 1 UNTIL M DO BEGIN 00548000 + FOR J:= 1 STEP 1 UNTIL N DO 00549000 + X[I,J]:= Z[I,J] END; 00550000 +BOX2: FOR I:= 1 STEP 1 UNTIL M DO BEGIN 00551000 + FOR J:= N2 STEP 1 UNTIL NPP DO 00552000 + X[I,J]:= 0.0; X[I,N+I]:= 1.0 END; GO TO BOX3; 00553000 +BOX3: K:= 1; 00554000 +BOX4: THI:= 1; 00555000 +BOX5: ALI:= OMI:= 1; 00556000 + FOR J:= 1 STEP 1 UNTIL P DO PK[N+J]:= 0.0; 00557000 +BOX6: FOR I:= 1 STEP 1 UNTIL N DO 00558000 + PK[I]:= X[K,I] * W[I]; 00559000 +BOX7: GO TO OM[OMI]; 00560000 + OM1: FOR I:= 1 STEP 1 UNTIL K DO BEGIN SUM:= 0.0; 00561000 + FOR J:= 1 STEP 1 UNTIL NPP DO 00562000 + SUM:= SUM + PK[J] * X[I,J]; QK[I]:= SUM END; 00563000 + GO TO BOX8; 00564000 + OM2: DK2:= 0.0; FOR I:= 1 STEP 1 UNTIL NPP DO 00565000 + DK2:= DK2 + PK[I] * X[K,I]; 00566000 + DK:= SQRT(DK2); 00567000 + FOR I:= 1 STEP 1 UNTIL NPP DO 00568000 + X[K,I]:= X[K,I]/DK;; 00569000 + OMI:= 1; GO TO BOX6; 00570000 +BOX8: FOR I:= 1 STEP 1 UNTIL K-1 DO 00571000 + QK[I]:= -QK[I]; QK[K]:= 1.0; 00572000 + FOR I:= 1 STEP 1 UNTIL NPP DO BEGIN 00573000 + SUM:= 0.0; FOR J:= 1 STEP 1 UNTIL K DO 00574000 + SUM:= SUM + X[J,I] * QK[J]; 00575000 + XP[I]:= SUM END; GO TO BOX9; 00576000 +BOX9: GO TO TH[THI]; 00577000 + TH1: FOR I:= 1 STEP 1 UNTIL NPP DO 00578000 + X[K,I]:= XP[I]; GO TO BOX10; 00579000 + TH2: FOR I:= 1 STEP 1 UNTIL P DO 00580000 + COF[I]:= -XP[N+I]; THI:= 3; GO TO TH1; 00581000 + TH3: GO TO BOX14; 00582000 +BOX10: GO TO AL[ALI]; 00583000 + AL1: OMI:= ALI:= 2; GO TO BOX6; 00584000 + AL2: IF K < M THEN BEGIN K:= K + 1; GO TO BOX4; END 00585000 + ELSE GO TO BOX12; 00586000 +BOX12: GO TO BE[BEI]; 00587000 + BE1: BEI:= THI:= 2; K:= K + 1; GO TO BOX14; 00588000 + BE2: GO TO BOX14; 00589000 +BOX14: IF RBAR = 0 THEN GO TO FINAL ELSE RBAR:= RBAR - 1; 00590000 + R1:= R1 + 1; THI:= 2; 00591000 + FOR I:= 1 STEP 1 UNTIL N DO 00592000 + X[M2,I]:= Y[R1,I]; 00593000 + FOR I:= 1 STEP 1 UNTIL P DO 00594000 + X[M2,N+I]:= 0.0; GO TO BOX5; 00595000 +FINAL: END ORTHO ; 00596000 +%***********************************************************************00597000 + 00598000 + 00598400 +%***********************************************************************00598500 +%***********************************************************************00599000 +% BEGIN EXECUTION (INNER BLOCK) %00600000 +%***********************************************************************00601000 +%***********************************************************************00601500 +BEGIN 00602000 + COMMENT INITALIZE & SET UP DATA FOR -ORTHO-; 00604000 + INTEGER I, J; 00605000 + REAL THETA1, ARCC1, ARCC2, ARCC12, DS1, DS2, DS12, XX, YY, TOTALARC; 00606000 +%% CALCULATE ARC LENGTHS OF EACH INTERFACE. 00607000 + ARCC12:= 2.0 *SQRT(1.0 - H**2); 00608000 + ARCC1:= 2.0 * (THETA1:= PI/2.0 - ARCTAN(2.0*H/ARCC12)); 00609000 + ARCC2:= 2.0 *PI - ARCC1; 00610000 +%% CALCULATE NO. OF POINTS ON EACH ARC. 00611000 + N1:= NPOINTS * ARCC1/(TOTALARC:= ARCC1 + ARCC2 + 2*ARCC12); 00612000 + N2:= NPOINTS *ARCC2/TOTALARC; 00613000 + IF (N12:= NPOINTS - N1 - N2)MOD 2 = 0 THEN 00614000 + N12:= N12 DIV 2 00615000 + ELSE BEGIN 00616000 + N12:= N12 DIV 2 + 1; 00617000 + N2:= N2 - 1; 00618000 + END; 00619000 +%% CALCULATE THE LENGTH OF LINE SEGMENTS BETWEEN POINTS ON EACH 00620000 +%% INTERFACE. 00621000 + DS1:= ARCC1/(N1 - 1); 00622000 + DS2:= ARCC2/(N2 - 1); 00623000 + DS12:= ARCC12/(N12 - 1); 00624000 + 00625000 + FOR J:= 1 STEP 1 UNTIL NPOINTS-N12 DO 00626000 + IF J <= N1 THEN %% DATA FOR ARC C1. 00627000 + BEGIN 00628000 + W[J]:= 1.0; Y[1,J]:= 0.25; 00629000 + XX:= SIN(YY:= DS1*(J-1)-THETA1); YY:= COS(YY); 00630000 + FOR I:= 1 STEP 1 UNTIL NAB DO 00631000 + Z[I,J]:= IF I <= NA THEN G(I, XX, YY) ELSE 0.0; 00632000 + END ELSE 00633000 + IF J <= N1 + N2 THEN %% DATA FOR ARC C2. 00634000 + BEGIN 00635000 + W[J]:= DS2/DS1; Y[1,J]:= 0.25; 00636000 + XX:= SIN(YY:= - THETA1 - DS2*(J-N1-1)); YY:= COS(YY); 00637000 + FOR I:= 1 STEP 1 UNTIL NAB DO 00638000 + Z[I,J]:= IF I <= NA THEN 0.0 ELSE G(I-NA, XX, YY); 00639000 + END ELSE %% DATA FOR ARC C12. 00640000 + BEGIN 00641000 + W[J]:= W[J+N12]:= DS12/DS1; 00642000 + Y[1,J]:= 0.25*((XX:= DS12*(J-N1-N2-1)-ARCC12/2)**2 + H**2)00643000 + * (1.0 - MU12); 00644000 + Y[1,J+N12]:= 0.0; 00645000 + FOR I:= 1 STEP 1 UNTIL NAB DO 00646000 + IF I <= NA THEN 00647000 + BEGIN 00648000 + Z[I,J]:= G(I, XX, H); 00649000 + Z[I,J+N12]:= DG(I, XX, H); 00650000 + END ELSE 00651000 + BEGIN 00652000 + Z[I,J]:= -G (I-NA, XX, H) * MU12; 00653000 + Z[I,J+N12]:= -DG (I-NA, XX, H); 00654000 + END; 00655000 + END; 00656000 + END INITIALIZATION; 00657000 + 00658000 +%%% THE NUMBER SMASHING COMMENCES. 00659000 + ORTHO (W, Y, Z, NPOINTS, NAB, NAB, COF[1,*]); %######################00660000 + 00661000 +%%% UPON EXIT FROM ORTHO, THE YU/SPARROW COEFFICIENTS, A[J], B[J], 00662000 +%%%% ARE CONTAINED IN THE VECTOR COF[1,*]. 00663000 + 00664000 + Q1Q1FULL:= 8.0/PI * INTEGRALUSTAR(1,TRUE); 00665000 + Q2Q2FULL:= 8.0/PI * INTEGRALUSTAR(2,TRUE); 00666000 + Q1STARQ1FULL:= EQUATION23(H); 00667000 + Q2STARQ2FULL:= 1.0 - Q1STARQ1FULL; 00668000 + TAB21:= Q1Q1FULL * PI/8.0; TAB22:= Q2Q2FULL * PI/8.0; 00669000 + Q1Q2:= Q1Q1FULL/Q2Q2FULL/MU12; 00670000 + Q1Q1STAR:= Q1Q1FULL/Q1STARQ1FULL; 00671000 + Q2Q2STAR:= Q2Q2FULL/Q2STARQ2FULL; 00671500 + 00672000 + 00673000 +DUMPIT: 00680000 +WRITE(OUTPUT, FT, TIME(0), TIME(2)/60, TIME(3)/60, (TIME(1)-STARTT)/60);00681000 +END INNER BLOCK; 00682000 + END H FOR LOOP; 00684000 +END. 00686000 diff --git a/tests/E-mode/YUSPARA-65.alg_m b/tests/E-mode/YUSPARA-65.alg_m new file mode 100644 index 0000000..2bc2fc9 --- /dev/null +++ b/tests/E-mode/YUSPARA-65.alg_m @@ -0,0 +1,344 @@ +%%%%%%%%%%%%%%%%%%%%FILTERED MANUALLY BY P.KIMPEL... 2012-08-25 00000100120825PK +% YUSPARA 00336000 +% FREE STANDING YU/SPARROW ROUTINE -- 3/5/70 -- P H KIMPEL 00337000 +% TEST OF ALGEBRAIC "G" AND "DG" FUNCTIONS -- 8/8/70 -- PHK 00338000 +BEGIN 00339000 +FILE YUSPAR (KIND=READER, DEPENDENTSPECS); 00340000120825PK +FILE OUTPUT (KIND=PRINTER, MAXRECSIZE=22, FRAMESIZE=48, FILEUSE=OUT); 00341000120825PK +FORMAT 00342000 + F1 (//X10"H/R = "F4.1,X15"MU1/MU2 = "I5), 00343000 + F2 (2("A[J]"X16"B[J]"X20)/8(2(2(I3,E14.5,X3)X5)X5,2A6"="E13.5/), 00344000 + 2(I3,E14.5,X3)X55"Q1/Q2"X7"="E13.5), 00345000 + F3 (X19"YU/SPARROW ANALYSIS -- ALGEBRAIC """G""" FUNCTIONS"X8,I2, 00346000 + 2("/"I2)X11"PAGE:"I4), 00347000120825PK + FE (//"RUN TOTALS:"), 00348000 + FT (/X9"PROCESSOR TIME: "R8.2,X15"ELAPSED TIME: "R8.2" SEC"); 00349000 +INTEGER 00350000 + I, J, % INDICIES. 00351000 + PT, % PROCESSOR TIME [1/64-S SEC] 00352000 + ET, % ELAPSED TIME [SAME]. 00353000 + PAG, % SET NUMBER. 00354000 + MO, % MONTH-OF-YEAR OF RUN. 00355000 + NPOINTS, % TOTAL # POINTS FOR APPROXIMATION. 00356000 + NA, % NUMBER OF A-PHASE COEFFICIENTS. 00357000 + NB, % NUMBER OF B-PHASE COEFFICIENTS. 00358000 + NAB, % NA + NB. 00359000 + N1 , % NUMBER OF POINTS ON ARC C1. 00360000 + N2, % NUMBER OF POINTS ON ARC C2. 00361000 + N12, % NUMBER OF POINTS ON ARC C12. 00362000 + DA, YR, % DAY- AND YEAR-OF-RUN. 00362100120825PK + KUSTAR, % COUNTER FOR "USTAR" 00363000 + NUSTAR, % COUNT LIMITS FOR USTAR. 00364000 + NSUSTAR; % DITTO. 00365000 +REAL 00366000 + H, % INTERFACE DISTANCE FROM PIPE CENTER. 00370000 + R, % PIPE RADIUS (= 1.0) 00371000 + MU1, % VISCOSITY OF THE MOST VISCOUS PHASE. 00372000 + MU2, % VISCOSITY OF THE LESS VISCOUS PHASE. 00373000 + PI, % USUAL MEANING: 3.14159... 00375000 + MU12, % MU1/MU2. 00376000 + ISUM, % USED BY "INTEGRALUSTAR". 00377000 + YSTRT, % DITTO 00378000 + YSTEP, % DITTO 00379000 + YEND, % DITTO 00380000 + XSTEP, % DITTO 00381000 + XEND, % DITTO 00382000 + XX, % COORDINATE WORKING VARIABLES. 00383000 + YY, % DITTO 00384000 + USUM, % WORKING VARIABLE VOR "INTEGRALUSTAR" 00385000 + Q1Q2, % Q1/Q2. 00386000 + X, Y, % WORKING VARIABLES. 00387000 + THETA1, % INTERFACE ANGLE FROM VERTICAL DIAMETER. 00388000 + ARCC1, % LENGTH OF ARC C1. 00388100 + ARCC2, % LENGTH OF ARC C2. 00388200 + ARCC12, % LENGTH OF ARC C12. 00388300 + DS1, % LENGTH OF ARC BETWEEN POINTS ON C1. 00388400 + DS2, % LENGTH OF ARC BETWEEN POINTS ON C2. 00388500 + DS12, % LENGTH OF ARC BETWEEN POINTS ON C12. 00388600 + TOTALARC; % ARCC1 + ARCC2 + ACRR12. 00388700 +LABEL JMP; 00389000 +INTEGER ARRAY NDX1, NDX2, NAME1, NAME2 [0:7]; 00391000 +ARRAY CALC[0:7]; 00392000 +DEFINE Q1Q1FULL = CALC[0]#, Q2Q2FULL = CALC[1]#, Q1STARQ1FULL=CALC[2]#, 00393000 + Q2STARQ2FULL = CALC[3]#, TAB22 = CALC[4]#, QMU12 = CALC[5]#, 00394000 + Q1Q1STAR = CALC[6]#, Q2Q2STAR = CALC[7]#; 00395000 + % BEGIN EXECUTION IN OUTER BLOCK. 00396000 +FILL NDX1[*] WITH 1, 3, 4, 7, 8, 11, 12, 15; 00397000 +FILL NDX2[*] WITH 2, 5, 6, 9, 10, 13, 14, 17; 00398000 +FILL NAME1[*] WITH "Q1/Q1F", "Q2/Q2F", "Q1*/Q1", "Q2*/Q2", "TABLE ", 00399000 + "QMU12 ", "Q1/Q1*", "Q2/Q2*"; 00400000 +FILL NAME2[*] WITH "ULL ", "ULL ", "FULL ", "FULL ", "2.2 ", 00401000 + " ", " ", " "; 00402000 + PI := 3.141592653589793; 00403000120825PK + NPOINTS := 100; NA := NB := 17; NAB := NA + NB; 00404000120825PK + R := 1; PAG := 2; 00405000120825PK + DA := TIME(0); YR := DA.[29:6]*10 + DA.[23:6]; 00406000120825PK + DA := DA.[17:6]*100 + DA.[11:6]*10 + DA.[5:6]; 00407000120825PK + FOR I := 31, IF YR MOD 4 = 0 THEN 29 ELSE 28, 31, 00408000120825PK + 30, 31, 30, 31, 31, 30, 31, 30, 31 DO 00409000 + BEGIN 00410000 + MO := MO + 1; 00411000120825PK + IF DA > I THEN DA := DA-I ELSE GO TO JMP; 00412000120825PK + END; 00413000 +JMP: 00414000120825PK +FOR H := -0.8, -0.6, -0.4, -0.2, 0.0, 0.2, 0.4, 0.6, 0.8 DO 00416000120825PK + FOR MU12 := 1, 2, 4, 10, 25, 100, 250, 1000 DO 00417000120825PK + 00418000 +BEGIN % ***** INNER BLOCK ***** 00419000 + REAL ARRAY 00420000 + W[1:NPOINTS], 00421000 + Y[1:1, 1:NPOINTS+NAB], 00422000 + Z[1:NAB, 1:NPOINTS+NAB], 00423000 + COF[1:1, 1:NAB]; 00424000 + 00425000 +%***********************************************************************00426000 +REAL PROCEDURE G (I, X, Y); 00427000 + VALUE I, X, Y; REAL X, Y; INTEGER I; 00428000 + CASE (I - 1) OF 00429000 + BEGIN 00430000 + G := 1; % 100431000120825PK + G := X; % 200432000120825PK + G := Y; % 300433000120825PK + G := X**2 - Y**2; % 400434000120825PK + G := 2.0 * X * Y; % 500435000120825PK + G := X**3 - 3.0 * X * Y**2; % 600436000120825PK + G := 3.0 * X**2 * Y - Y**3; % 700437000120825PK + G := X**4 + Y**4 - 6.0 * X**2 * Y**2; % 800438000120825PK + G := 4.0 * X**3 * Y - 4.0 * X * Y**3; % 900439000120825PK + G := X**5 - 10.0 * X**3 * Y**2 + 5.0 * X * Y**4; %1000440000120825PK + G := Y**5 - 10.0 * Y**3 * X**2 + 5.0 * Y * X**4; %1100441000120825PK + G := X**6 - 15.0 * X**4 * Y**2 + 15.0 * X**2 * Y**4 - Y**6; %1200442000120825PK + G := 6.0 * X**5 * Y + 6.0 * X * Y**5 - 20.0 * X**3 * Y**3; %1300443000120825PK + G := X**7 - 21.0 *X**5 * Y**2 + 35.0 * X**3 * Y**4 - 7.0 * X *%1400444000120825PK + Y**6; 00444100120825PK + G := 7.0 *X**6 * Y - 35.0 * X**4 * Y**3 + 21.0 * X**2 * Y**5 - %1500445000120825PK + Y**7; 00445100120825PK + G := X**8 + Y**8 - 28.0 *X**6 * Y**2 + 70.0 * X**4 * Y**4 - 28.0 00446000120825PK + * X**2 * Y**6; %1600447000120825PK + G := 8.0 * X**7 *Y - 56.0 * X**5 * Y**3 + 56.0 * X**3 * Y**5 00448000120825PK + - 8.0 *X * Y**7; %1700449000120825PK + END G; 00450000 +%***********************************************************************00452000 +REAL PROCEDURE DG (I, X, Y); 00453000 + VALUE I, X, Y; INTEGER I; REAL X, Y; 00454000 + CASE (I - 1) OF 00455000 + BEGIN 00456000 + DG := 0.0; % 100457000120825PK + DG := 0.0; % 200458000120825PK + DG := 1.0; % 300459000120825PK + DG := -2.0 * Y; % 400460000120825PK + DG := 2.0 * X; % 500461000120825PK + DG := -6.0 * X *Y; % 600462000120825PK + DG := 3.0 * X**2 - 3.0 * Y**2; % 700463000120825PK + DG := 4.0 * Y**3 - 12.0 *X**2 * Y; % 800464000120825PK + DG := 4.0 * X**3 - 12.0 * X * Y**2; % 900465000120825PK + DG := -20.0 *X**3 * Y + 20.0 * X * Y**3; %1000466000120825PK + DG := 5.0 * Y**4 - 30.0 * Y**2 *X**2 + 5.0 * X**4; %1100467000120825PK + DG := -30.0 * X**4 *Y + 60.0 * X**2 * Y**3 - 6.0 * Y**5; %1200468000120825PK + DG := 6.0 * X**5 + 30.0 * X * Y**4 - 60.0 * X**3 * Y**2; %1300469000120825PK + DG := -42.0 * X**5 * Y + 140.0 * X**3 * Y**3 - 42.0 * X * Y**5%1400470000120825PK + ; 00470100120825PK + DG := 7.0 * X**6 - 105.0 * X**4 * Y**2 + 105.0 * X**2 * Y**4 00471000120825PK + - 7.0 * Y**6; %1500472000120825PK + DG := 8.0 *Y**7 - 56.0 * X**6 * Y + 280.0 * X**4 * Y**3 00473000120825PK + - 168.0 * X**2 * Y**5; %1600474000120825PK + DG := 8.0 * X**7 - 168.0 * X**5 * Y**2 + 280.0 * X**3 * Y**4 00475000120825PK + - 56.0 * X * Y**6; %1700476000120825PK + END DG; 00477000 +%***********************************************************************00479000 +REAL PROCEDURE USTAR (I, X, Y); 00480000 + VALUE I, X, Y; INTEGER I; REAL X, Y; 00481000 + BEGIN 00482000 + NUSTAR := IF I = 1 THEN NA ELSE NA + NB; 00483000120825PK + NSUSTAR := IF I = 1 THEN 1 ELSE NA + 1; 00484000120825PK + USUM := 0.0; 00485000120825PK + FOR KUSTAR := NSUSTAR STEP 1 UNTIL NUSTAR DO 00486000120825PK + USUM := USUM + COF[1,KUSTAR] * G(KUSTAR-NSUSTAR+1, X, Y); 00487000120825PK + USTAR := USUM; 00488000120825PK + END USTAR; 00489000 +%***********************************************************************00491000 +REAL PROCEDURE INTEGRALUSTAR (I); 00492000 + VALUE I; INTEGER I; 00493000 + BEGIN 00494000 + ISUM := 0.0; 00495000120825PK + IF I = 1 THEN BEGIN YSTRT := H; YEND := 0.9999; END 00496000120825PK + ELSE BEGIN YSTRT := -0.9999; YEND := H; END; 00497000120825PK + YSTEP := (YEND - YSTRT)/30.00; 00498000120825PK + FOR YY := YSTRT STEP YSTEP UNTIL YEND DO 00499000120825PK + BEGIN XSTEP := (XX := SQRT(1.0 - YY**2))/15.0; 00500000120825PK + FOR XX := -XX STEP XSTEP UNTIL 0.0 DO 00501000120825PK + ISUM := ISUM + XSTEP * YSTEP * (USTAR (I, XX, YY) 00502000120825PK + - (XX**2 + YY**2)/4); 00503000120825PK + END; 00504000 + INTEGRALUSTAR := 2.0 * ISUM; 00505000120825PK + END INTEGRALUSTAR; 00506000 +%***********************************************************************00508000 +REAL PROCEDURE ARCSIN (X); 00509000 + VALUE X; REAL X; 00510000 + ARCSIN := ARCTAN(X/SQRT(1.0 - X**2)); 00511000120825PK +%***********************************************************************00513000 +REAL PROCEDURE EQUATION23 (A); 00514000 + VALUE A; REAL A; 00515000 + EQUATION23 := 0.5 - 2.0/3.0/PI * (A * (1.0 - A**2)**1.5 00516000120825PK + + 1.5 * A *SQRT(1.0 - A**2) + 1.5 * ARCSIN(A)); 00517000120825PK +%***********************************************************************00519000 +PROCEDURE ORTHO (W, Y, Z, N, M, P, COF); 00521000 +VALUE N, M, P; 00522000 +REAL ARRAY COF, W[1]; 00523000 +REAL ARRAY Y, Z [1, 1]; 00524000 +INTEGER N, M, P; 00525000 +% 00526000 +%COMMENT 00527000 +% -ORTHO- IS TAKEN FROM ACM ALGORITHM 127 [COMM. ACM, VOL.5, 00528000 +% OCTOBER 1962, P. 511, AUTHOR: PHILIP J. WALSH], 00529000 +% AND INCORPORATES THE SUGGESTIONS IN BARRODALE-S 00530000 +% CERTIFICATION [COMM. ACM, VOL.13, #2, FEB.1970, 122] 00531000 +% CODE REDUCED TO ONLY THAT NEEDED TO PERFORM LEAST-SQUARES 00532000 +% DETERMINATION OF COEFFICIENTS [4/70]. 00533000 +% .... P H KIMPEL 3/25/70 U OF DELAWARE 00534000 +% 00535000 +BEGIN 00536000 +INTEGER NPP, M1, N2, M2, R1, RBAR, P2, BEI, I, J, K, THI, ALI, 00537000 + OMI, R; 00538000 +REAL ARRAY PK, XP [1:N+P], QK[1:M+1], X[1:M+1, 1:N+P]; 00539000 +REAL SUM, DK2, DK; 00540000 +LABEL BOX1, BOX2, BOX3, BOX4, BOX5, BOX6, BOX7, OM1, OM2, BOX8, BOX9, 00541000 + TH1, TH2, TH3, BOX10, AL1, AL2, BOX12, BE1, BE2, BOX14, FINAL; 00542000 +SWITCH BE := BE1, BE2; 00543000120825PK +SWITCH TH := TH1, TH2, TH3; SWITCH AL := AL1, AL2; 00544000120825PK +SWITCH OM := OM1, OM2; 00545000120825PK + NPP := N + P; R := 1; M1 := M-1; N2 := N+1; M2 := M+1; 00546000120825PK + R1 := 0; RBAR := R; P2 := P+1; BEI := 1; 00547000120825PK +BOX1: FOR I := 1 STEP 1 UNTIL M DO BEGIN 00548000120825PK + FOR J := 1 STEP 1 UNTIL N DO 00549000120825PK + X[I,J] := Z[I,J] END; 00550000120825PK +BOX2: FOR I := 1 STEP 1 UNTIL M DO BEGIN 00551000120825PK + FOR J := N2 STEP 1 UNTIL NPP DO 00552000120825PK + X[I,J] := 0.0; X[I,N+I] := 1.0 END; GO TO BOX3; 00553000120825PK +BOX3: K := 1; 00554000120825PK +BOX4: THI := 1; 00555000120825PK +BOX5: ALI := OMI := 1; 00556000120825PK + FOR J := 1 STEP 1 UNTIL P DO PK[N+J] := 0.0; 00557000120825PK +BOX6: FOR I := 1 STEP 1 UNTIL N DO 00558000120825PK + PK[I] := X[K,I] * W[I]; 00559000120825PK +BOX7: GO TO OM[OMI]; 00560000 + OM1: FOR I := 1 STEP 1 UNTIL K DO BEGIN SUM := 0.0; 00561000120825PK + FOR J := 1 STEP 1 UNTIL NPP DO 00562000120825PK + SUM := SUM + PK[J] * X[I,J]; QK[I] := SUM END; 00563000120825PK + GO TO BOX8; 00564000 + OM2: DK2 := 0.0; FOR I := 1 STEP 1 UNTIL NPP DO 00565000120825PK + DK2 := DK2 + PK[I] * X[K,I]; 00566000120825PK + DK := SQRT(DK2); 00567000120825PK + FOR I := 1 STEP 1 UNTIL NPP DO 00568000120825PK + X[K,I] := X[K,I]/DK;; 00569000120825PK + OMI := 1; GO TO BOX6; 00570000120825PK +BOX8: FOR I := 1 STEP 1 UNTIL K-1 DO 00571000120825PK + QK[I] := -QK[I]; QK[K] := 1.0; 00572000120825PK + FOR I := 1 STEP 1 UNTIL NPP DO BEGIN 00573000120825PK + SUM := 0.0; FOR J := 1 STEP 1 UNTIL K DO 00574000120825PK + SUM := SUM + X[J,I] * QK[J]; 00575000120825PK + XP[I] := SUM END; GO TO BOX9; 00576000120825PK +BOX9: GO TO TH[THI]; 00577000 + TH1: FOR I := 1 STEP 1 UNTIL NPP DO 00578000120825PK + X[K,I] := XP[I]; GO TO BOX10; 00579000120825PK + TH2: FOR I := 1 STEP 1 UNTIL P DO 00580000120825PK + COF[I] := -XP[N+I]; THI := 3; GO TO TH1; 00581000120825PK + TH3: GO TO BOX14; 00582000 +BOX10: GO TO AL[ALI]; 00583000 + AL1: OMI := ALI := 2; GO TO BOX6; 00584000120825PK + AL2: IF K < M THEN BEGIN K := K + 1; GO TO BOX4; END 00585000120825PK + ELSE GO TO BOX12; 00586000 +BOX12: GO TO BE[BEI]; 00587000 + BE1: BEI := THI := 2; K := K + 1; GO TO BOX14; 00588000120825PK + BE2: GO TO BOX14; 00589000 +BOX14: IF RBAR = 0 THEN GO TO FINAL ELSE RBAR := RBAR - 1; 00590000120825PK + R1 := R1 + 1; THI := 2; 00591000120825PK + FOR I := 1 STEP 1 UNTIL N DO 00592000120825PK + X[M2,I] := Y[R1,I]; 00593000120825PK + FOR I := 1 STEP 1 UNTIL P DO 00594000120825PK + X[M2,N+I] := 0.0; GO TO BOX5; 00595000120825PK +FINAL: END ORTHO ; 00596000 +%***********************************************************************00597000 +%***********************************************************************00599000 +% BEGIN EXECUTION (INNER BLOCK) %00600000 +%***********************************************************************00601000 + ET := TIME(1); PT := TIME(2); 00602000120825PK + COMMENT INITALIZE & SET UP DATA FOR -ORTHO-; 00604000 +%% CALCULATE ARC LENGTHS OF EACH INTERFACE. 00607000 + ARCC12 := 2.0 *SQRT(1.0 - H**2); 00608000120825PK + ARCC1 := 2.0 * (THETA1 := PI/2.0 - ARCTAN(2.0*H/ARCC12)); 00609000120825PK + ARCC2 := 2.0 *PI - ARCC1; 00610000120825PK +%% CALCULATE NO. OF POINTS ON EACH ARC. 00611000 + N1 := NPOINTS * ARCC1/(TOTALARC := ARCC1 + ARCC2 + 2*ARCC12); 00612000120825PK + N2 := NPOINTS * ARCC2 / TOTALARC; 00613000120825PK + IF (N12 := NPOINTS - N1 - N2) MOD 2 = 0 THEN 00614000120825PK + N12 := N12 DIV 2 00615000120825PK + ELSE BEGIN 00616000 + N12 := N12 DIV 2 + 1; 00617000120825PK + N2 := N2 - 1; 00618000120825PK + END; 00619000 +%% CALCULATE THE LENGTH OF LINE SEGMENTS BETWEEN POINTS ON EACH 00620000 +%% INTERFACE. 00621000 + DS1 := ARCC1/(N1 - 1); 00622000120825PK + DS2 := ARCC2/(N2 - 1); 00623000120825PK + DS12 := ARCC12/(N12 - 1); 00624000120825PK + FOR I := 1 STEP 1 UNTIL NAB DO 00626000120825PK + BEGIN 00627000 + FOR J := 1 STEP 1 UNTIL N1 DO %% DATA FOR ARC C1. 00628000120825PK + BEGIN 00629000 + W[J] := 1; Y[1,J] := 0.25; 00630000120825PK + XX := SIN(YY := DS1*(J-1)-THETA1); YY := COS(YY); 00631000120825PK + Z[I,J] := IF I LEQ NA THEN G(I, XX, YY) ELSE 0; 00632000120825PK + END; 00633000 + FOR J := N1+1 STEP 1 UNTIL N1+N2 DO %% DATA FOR ARC C2. 00634000120825PK + BEGIN 00635000 + W[J] := DS2/DS1; Y[1,J] := 0.25; 00636000120825PK + XX := SIN(YY := -THETA1 - DS2*(J-N1-1)); YY := COS(YY); 00637000120825PK + Z[I,J] := IF I LEQ NA THEN 0 ELSE G(I-NA, XX, YY); 00638000120825PK + END; 00639000 + FOR J := N1+N2+1 STEP 1 UNTIL NPOINTS-N12 DO %% DATA FOR ARC C12.00640000120825PK + BEGIN 00641000 + W[J] := W[J+N12] := DS12/DS1; 00642000120825PK + Y[1,J] := 0.25 * ((XX := DS12*(J-N1-N2-1)-ARCC12/2)**2 + H**2)00643000120825PK + * (1 - MU12); 00644000120825PK + Y[1,J+N12] := 0; 00645000120825PK + IF I LEQ NA THEN 00646000120825PK + BEGIN 00647000 + Z[I,J] := G(I, XX, H); 00648000120825PK + Z[I,J+N12] := DG(I, XX, H); 00649000120825PK + END ELSE 00650000 + BEGIN 00651000 + Z[I,J] := -G(I-NA, XX, H) * MU12; 00652000120825PK + Z[I,J+N12] := -DG(I-NA, XX, H); 00653000120825PK + END; 00654000 + END; 00655000 + END INITIALIZATION; 00656000 + 00658000 +%%% THE NUMBER SMASHING COMMENCES. 00659000 + ORTHO (W, Y, Z, NPOINTS, NAB, NAB, COF[1,*]); %######################00660000 +%%% UPON EXIT FROM ORTHO, THE YU/SPARROW COEFFICIENTS, A[J], B[J], 00662000 +%%%% ARE CONTAINED IN THE VECTOR COF[1,*]. 00663000 + Q1Q1FULL := 8/PI * INTEGRALUSTAR(1); 00665000120825PK + Q2Q2FULL := 8/PI * INTEGRALUSTAR(2); 00666000120825PK + Q1STARQ1FULL := EQUATION23(H); 00667000120825PK + Q2STARQ2FULL := 1 - Q1STARQ1FULL; 00668000120825PK + TAB22 := Q2Q2FULL * PI/8; 00669000120825PK + Q1Q2 := Q1Q1FULL/Q2Q2FULL/MU12; 00670000120825PK + QMU12 := Q1Q2 * MU12; 00670500120825PK + Q1Q1STAR := Q1Q1FULL/Q1STARQ1FULL; Q2Q2STAR := Q2Q2FULL/Q2STARQ2FULL;00671000120825PK + IF (PAG := PAG+1) MOD 3 = 0 THEN 00674000120825PK + BEGIN 00674800 + WRITE (OUTPUT[SKIP 1]); 00674900120825PK + WRITE (OUTPUT, F3, MO, DA, YR, PAG DIV 3); 00675000 + END; 00675100 + WRITE (OUTPUT[SPACE 2], F1, H, MU12); 00676000120825PK + WRITE (OUTPUT, F2, FOR I := 0 STEP 1 UNTIL 7 DO 00677000120825PK + [NDX1[I], COF[1,NDX1[I]], NDX1[I], COF[1,NDX1[I]+NA], NDX2[I], 00678000 + COF[1,NDX2[I]], NDX2[I], COF[1,NDX2[I]+NA], NAME1[I], 00679000 + NAME2[I], CALC[I]], 16, COF[1,16], 16, COF[1,16], Q1Q2); 00680000 + WRITE (OUTPUT, FT, (TIME(2) - PT)/60, (TIME(1) - ET)/60); 00681000 +END INNER BLOCK; 00682000 +WRITE (OUTPUT[SPACE 2], FE); 00683000120825PK +WRITE (OUTPUT[SPACE 2], *, "TOTALS:"); 00684000120825PK +WRITE (OUTPUT, FT, TIME(2)/60, TIME(1)/60); 00685000 +END. 00686000