code 34432; procedure PRAXIS(N, X, FUNCT, IN, OUT); value N; integer N; array X, IN, OUT; real procedure FUNCT; begin commentTHIS PROCEDURE MINIMIZES FUNCT(N,X),WITH THE PRINCIPAL AXIS METHOD (SEE BRENT,R.P, 1973, ALGORITHMS FOR MINIMIZATION WITHOUT DERIVATIVES,CH.7); procedure SORT; begin integer I, J, K; real S; for I:= 1 step 1 until N - 1 do begin K:= I; S:= D[I]; for J:= I+1 step 1 until N do if D[J]>S then begin K:= J; S:= D[J] end; if K>I then begin D[K]:= D[I]; D[I]:= S; for J:= 1 step 1 until N do begin S:=V[J,I]; V[J,I]:= V[J,K]; V[J,K]:= S end end end end SORT procedure MIN(J, NITS, D2, X1, F1, FK); value J, NITS, FK; integer J, NITS; real D2, X1, F1; boolean FK; begin real procedure FLIN(L); value L; real L; begin integer I; array T[1:N]; if J > 0 then begin for I:= 1 step 1 until N do T[I]:= X[I] + L * V[I,J] end else begin comment SEARCH ALONG PARABOLIC SPACE CURVE; QA:= L * (L - QD1) / (QD0 * (QD0 + QD1)); QB:= (L + QD0) * (QD1 - L) /(QD0 * QD1); QC:= L * (L + QD0) / (QD1 * (QD0 + QD1)); for I:= 1 step 1 until N do T[I]:= QA * Q0[I] +QB * X[I] + QC * Q1[I] end; NF:= NF + 1; FLIN:= FUNCT(N, T) end FLIN; integer K; boolean DZ; real X2, XM, F0, F2, FM, D1, T2, S, SF1, SX1; SF1:= F1; SX1:= X1; K:= 0; XM:= 0; F0:= FM:= FX; DZ:= D2 < RELTOL; S:= SQRT(VECVEC(1,N,0,X,X)); T2:= M4 * SQRT(ABS(FX) / (if DZ then DMIN else D2) + S * LDT) + M2 * LDT; S:= S * M4 + ABSTOL; if DZ and T2 > S then T2:= S; if T2<SMALL"THEN"T2:= SMALL; if T2>0.01*H then T2:= 0.01*H; if FK"AND"F1<=FM then begin XM:=X1; FM:= F1 end; if ^ FK"OR"ABS(X1)<T2"THEN" begin X1:=if X1>0 then T2"ELSE"-T2; F1:= FLIN(X1) end; if F1<= FM"THEN" begin XM:= X1; FM:= F1 end; L0: if DZ then begin commentEVALUATE FLIN AT ANOTHER POINT AND ESTIMATE THE SECOND DERIVATIVE; X2:= if F0 < F1 then -X1 else X1 * 2; F2:= FLIN(X2); if F2 <= FM then begin XM:= X2; FM:= F2 end; D2:=(X2*(F1-F0)-X1*(F2-F0))/(X1*X2*(X1-X2)) end; commentESTIMATE FIRST DERIVATIVE AT 0; D1:=(F1-F0)/X1-X1*D2; DZ:=true; X2:= if D2<=SMALL"THEN" (if D1<0then H"ELSE"-H) else -0.5*D1/D2; if ABS(X2)>H"THEN"X2:=if X2>0then H"ELSE"-H; L1: F2:=FLIN(X2); if K<NITS"AND"F2>F0"THEN" begin K:=K+1; if F0<F1"AND"X1*X2>0then goto L0; X2:= 0.5*X2; goto L1 end; NL:= NL+1; if F2>FM"THEN"X2:=XM"ELSE"FM:=F2; D2:=if ABS(X2*(X2-X1))>SMALL"THEN" (X2*(F1-F0)-X1*(FM-F0))/(X1*X2*(X1-X2)) else if K>0then 0else D2; if D2<=SMALL"THEN"D2:=SMALL; X1:=X2; FX:=FM; if SF1<FX"THEN" begin FX:=SF1; X1:=SX1 end; if J>0then ELMVECCOL(1,N,J,X,V,X1) end MIN; procedure QUAD; begin integer I; real L, S; S:= FX; FX:= QF1; QF1:= S; QD1:= 0; for I:= 1 step 1 until N do begin S:=X[I]; X[I]:= L:= Q1[I]; Q1[I]:= S; QD1:= QD1 + (S - L) ** 2 end; L:=QD1:=SQRT(QD1); S:= 0; if (QD0*QD1>DWARF)and NL>=3*N*N"THEN" begin MIN(0,2,S,L,QF1,true ); QA:= L*(L-QD1)/(QD0*(QD0+QD1)); QB:=(L+QD0)*(QD1-L)/(QD0*QD1); QC:= L*(L+QD0)/(QD1*(QD0+QD1)) end else begin FX:= QF1; QA:= QB:= 0; QC:= 1 end; QD0:= QD1;for I:= 1step 1until N"DO" begin S:=Q0[I]; Q0[I]:=X[I]; X[I]:= QA*S + QB*X[I]+QC*Q1[I] end end QUAD; boolean ILLC; integer I, J, K, K2, NL, MAXF, NF, KL, KT, KTM; real S, SL, DN, DMIN, FX, F1, LDS, LDT, SF, DF, QF1, QD0, QD1, QA, QB, QC, M2, M4, SMALL, VSMALL, LARGE, VLARGE, SCBD, LDFAC,T2, MACHEPS, RELTOL, ABSTOL, H; array V[1:N,1:N], D, Y, Z, Q0, Q1[1:N]; MACHEPS:= IN[0]; RELTOL:= IN[1]; ABSTOL:= IN[2]; MAXF:= IN[5]; H:= IN[6]; SCBD:= IN[7]; KTM:= IN[8]; ILLC:= IN[9] < 0; SMALL:= MACHEPS ** 2; VSMALL:= SMALL ** 2; LARGE:= 1/SMALL; VLARGE:= 1/VSMALL; M2:= RELTOL; M4:= SQRT(M2); SETRANDOM(0.5); LDFAC:= if ILLC then 0.1 else 0.01; KT:=NL:=0; NF:=1; OUT[3]:= QF1:=FX:=FUNCT(N,X); ABSTOL:=T2:= SMALL+ABS(ABSTOL); DMIN:= SMALL; if H<ABSTOL*100then H:=ABSTOL*100; LDT:=H; INIMAT(1,N,1,N,V,0); for I:=1step 1until N"DO"V[I,I]:= 1; D[1]:= QD0:= 0; DUPVEC(1,N,0,Q1,X); INIVEC(1,N,Q0,0); commentMAIN LOOP; L0: SF:=D[1]; D[1]:= S:= 0; MIN(1,2,D[1],S,FX,false ); if S <= 0 then MULCOL(1, N, 1, 1, V, V, -1); if SF <= 0.9 * D[1] or 0.9 * SF >= D[1] then INIVEC(2,N,D,0); for K:= 2step 1until N"DO" begin DUPVEC(1,N,0,Y,X); SF:=FX; ILLC:= ILLC or KT>0; L1: KL:=K; DF:= 0; if ILLC then begin commentRANDOM STOP TO GET OFF RESULTION VALLEY; for I:= 1 step 1until N"DO" begin S:=Z[I]:=(0.1*LDT+T2*10**KT) *(RANDOM-0.5); ELMVECCOL(1,N,I,X,V,S) end; FX:= FUNCT(N,X); NF:= NF+1 end; for K2:= K step 1 until N do begin SL:=FX; S:= 0; MIN (K2, 2, D[K2], S, FX, false ); S:=if ILLC then D[K2] * (S + Z[K2]) ** 2 else SL-FX;if DF<S"THEN" begin DF:=S;KL:= K2"END"; end; if ^ILLC and DF < ABS(100 * MACHEPS * FX) then begin ILLC:= true; goto L1 end; for K2:= 1step 1until K-1do begin S:= 0; MIN(K2, 2, D[K2], S, FX, false ) end; F1:= FX; FX:= SF; LDS:= 0; for I:= 1 step 1 until N do begin SL:= X[I]; X[I]:= Y[I]; SL:= Y[I]:= SL - Y[I]; LDS:= LDS + SL * SL end; LDS:= SQRT(LDS); if LDS > SMALL then begin for I:= KL - 1 step -1 until K do begin for J:= 1 step 1 until N do V[J, I + 1]:= V[J,I]; D[I + 1]:= D[I] end; D[K]:= 0; DUPCOLVEC(1, N, K, V, Y); MULCOL(1, N, K, K, V, V, 1 / LDS); MIN(K, 4, D[K], LDS, F1, true ); if LDS <= 0 then begin LDS:= LDS; MULCOL(1, N, K, K, V, V, -1) end end; LDT:= LDFAC * LDT; if LDT < LDS then LDT:= LDS; T2:= M2 * SQRT(VECVEC(1, N, 0, X, X)) + ABSTOL; KT:= if LDT > 0.5 * T2 then 0 else KT + 1; if KT > KTM then begin OUT[1]:= 0; goto L2 end end; QUAD; DN:= 0;for I:= 1step 1until N"DO" begin D[I]:= 1/SQRT(D[I]); if DN<D[I]then DN:=D[I] end; for J:= 1step 1until N"DO" begin S:= D[J]/DN; MULCOL(1,N,J,J,V,V,S)end; if SCBD>1then begin S:=VLARGE; for I:=1 step 1until N"DO" begin SL:= Z[I]:= SQRT(MATTAM(1, N, I, I, V, V)); if SL<M4"THEN"Z[I]:= M4; if S>SL then S:= SL end; for I:=1step 1until N"DO" begin SL:=S/Z[I];Z[I]:= 1/SL; if Z[I]>SCBD"THEN" begin SL:=1/SCBD; Z[I]:= SCBD"END"; MULROW(1, N, I, I, V, V, SL) end end; for I:= 1 step 1 until N do ICHROWCOL(I + 1, N, I, I, V); begin array A[1:N,1:N], EM[0:7]; EM[0]:= EM[2]:= MACHEPS; EM[4]:= 10 * N; EM[6]:= VSMALL; DUPMAT(1, N, 1, N, A, V); if QRISNGVALDEC(A, N, N, D, V, EM) ^= 0 then begin OUT[1]:= 2; goto L2 end; end; if SCBD>1then begin for I:=1step 1until N"DO" MULROW(1,N,I,I,V,V,Z[I]); for I:= 1step 1until N"DO" begin S:= SQRT(TAMMAT(1,N,I,I,V,V)); D[I]:= S*D[I]; S:= 1/S; MULCOL(1,N,I,I,V,V,S) end end; for I:= 1 step 1 until N do begin S:= DN * D[I]; D[I]:= if S > LARGE then VSMALL else if S < SMALL then VLARGE else S ** (-2) end; SORT; DMIN:= D[N]; if DMIN < SMALL then DMIN:= SMALL; ILLC:= (M2 * D[1]) > DMIN; if NF < MAXF then goto L0 else OUT[1]:= 1; L2: OUT[2]:= FX; OUT[4]:= NF; OUT[5]:= NL; OUT[6]:= LDT end PRAXIS; eop