code 33015 ; procedure RK3N(X, A, B, Y, YA, Z, ZA, FXYJ, J, E, D, FI, N); value B, FI, N; integer J, N; real X, A, B, FXYJ; boolean FI; array Y, YA, Z, ZA, E, D; begin integer JJ; real XL, H, HMIN, INT, HL, ABSH, FHM, DISCRY, DISCRZ, TOLY, TOLZ, MU, MU1, FHY, FHZ; boolean LAST, FIRST, REJECT; array YL, ZL, K0, K1, K2, K3, K4, K5[1:N], EE[1:4 * N]; if FI then begin D[3]:= A; for JJ:= 1 step 1 until N do begin D[JJ + 3]:= YA[JJ]; D[N + JJ + 3]:= ZA[JJ] end D[1]:= 0; XL:= D[3]; for JJ:= 1 step 1 until N do begin YL[JJ]:= D[JJ + 3]; ZL[JJ]:= D[N + JJ + 3] end; if FI then D[2]:= B - D[3]; ABSH:= H:= ABS(D[2]); if B - XL < 0 then H:= - H; INT:= ABS(B - XL); HMIN:= INT * E[1] + E[2]; for JJ:= 2 step 1 until 2 * N do begin HL:= INT * E[2 * JJ - 1] + E[2 * JJ]; if HL < HMIN then HMIN:= HL end; for JJ:= 1 step 1 until 4 * N do EE[JJ]:= E[JJ] / INT; FIRST:= REJECT:= true; if FI then begin LAST:= true; goto STEP end; TEST: ABSH:= ABS(H); if ABSH < HMIN then begin H:= if H > 0 then HMIN else - HMIN; ABSH:= HMIN end; if H >= B - XL eqv H >= 0 then begin D[2]:= H; LAST:= true; H:= B - XL; ABSH:= ABS(H) end else LAST:= false; STEP: if REJECT then begin X:= XL; for JJ:= 1 step 1 until N do Y[JJ]:= YL[JJ]; for J:= 1 step 1 until N do K0[J]:= FXYJ * H end else begin FHY:= H / HL; for JJ:= 1 step 1 until N do K0[JJ]:= K5[JJ] * FHY end; X:= XL + .27639 3202250021 * H; for JJ:= 1 step 1 until N do Y[JJ]:= YL[JJ] + (ZL[JJ] * .276393202250021 + K0[JJ] * .038196601125011) * H; for J:= 1 step 1 until N do K1[J]:= FXYJ * H; X:= XL + .723606797749979 * H; for JJ:= 1 step 1 until N do Y[JJ]:= YL[JJ] + (ZL[JJ] * .723606797749979 + K1[JJ] * .261803398874989) * H; for J:= 1 step 1 until N do K2[J]:= FXYJ * H; X:= XL + H * .5; for JJ:= 1 step 1 until N do Y[JJ]:= YL[JJ] + (ZL[JJ] * .5 + K0[JJ] * .046875 + K1[JJ] * .079824155839840 - K2[JJ] * .00169 9155839840) * H; for J:= 1 step 1 until N do K4[J]:= FXYJ * H; X:= if LAST then B else XL + H; for JJ:= 1 step 1 until N do Y[JJ]:= YL[JJ] + (ZL[JJ] + K0[JJ] * .309016994374947 + K2[JJ] * .190983005625053) * H; for J:= 1 step 1 until N do K3[J]:= FXYJ * H; for JJ:= 1 step 1 until N do Y[JJ]:= YL[JJ] + (ZL[JJ] + K0[JJ] * .083333333333333 + K1[JJ] * .30150 2832395825 + K2[JJ] * .115163834270842) * H; for J:= 1 step 1 until N do K5[J]:= FXYJ * H; REJECT:= false; FHM:= 0; for JJ:= 1 step 1 until N do begin DISCRY:= ABS(( - K0[JJ] * .5 + K1[JJ] * 1.809016994374947 + K2[JJ] * .690983005625053 - K4[JJ] * 2) * H); DISCRZ:= ABS((K0[JJ] - K3[JJ]) * 2 - (K1[JJ] + K2[JJ]) * 10 + K4[JJ] * 16 + K5[JJ] * 4); TOLY:= ABSH * (ABS(ZL[JJ]) * EE[2 * JJ - 1] + EE[2 * JJ]); TOLZ:= ABS(K0[JJ]) * EE[2 * (JJ + N) - 1] + ABSH * EE[2 * (JJ + N)]; REJECT:= DISCRY > TOLY or DISCRZ > TOLZ or REJECT; FHY:= DISCRY / TOLY; FHZ:= DISCRZ / TOLZ; if FHZ > FHY then FHY:= FHZ; if FHY > FHM then FHM:= FHY end; MU:= 1 / (1 + FHM) + .45; if REJECT then begin if ABSH <= HMIN then begin D[1]:= D[1] + 1; for JJ:= 1 step 1 until N do begin Y[JJ]:= YL[JJ]; Z[JJ]:= ZL[JJ] end; FIRST:= true; goto NEXT end; H:= MU * H; goto TEST end REJ; if FIRST then begin FIRST:= false; HL:= H; H:= MU * H; goto ACC end; FHY:= MU * H / HL + MU - MU1; HL:= H; H:= FHY * H; ACC: MU1:= MU; for JJ:= 1 step 1 until N do Z[JJ]:= ZL[JJ] + (K0[JJ] + K3[JJ]) * .083333333333333 + (K1[JJ] + K2[JJ]) * .416666666666667; NEXT: if B ^= X then begin XL:= X; for JJ:= 1 step 1 until N do begin YL[JJ]:= Y[JJ]; ZL[JJ]:= Z[JJ] end; goto TEST end; if not LAST then D[2]:= H; D[3]:= X; for JJ:= 1 step 1 until N do begin D[JJ + 3]:= Y[JJ]; D[N + JJ + 3]:= Z[JJ] end end RK3N; eop