code 34441; procedure GSSNEWTON(M, N, PAR, RV, JJINV, FUNCT, JACOBIAN, IN, OUT); value M, N; integer M, N; array PAR, RV, JJINV, IN, OUT; boolean procedure FUNCT; procedure JACOBIAN; begin integer I, J, INR, MIT, TEXT, IT, ITMAX, INRMAX, TIM, FEVAL, FEVALMAX; real RHO, RES1, RES2, RN, RELTOLPAR, ABSTOLPAR, ABSTOLRES, STAP, NORMX; boolean CONV, TESTTHF, DAMPING ON; array JAC[1:M + 1,1:N], PR, AID, SOL[1 : N], FU2[1 : M], AUX[2 : 5]; integer array CI[1:N]; boolean procedure LOC FUNCT(M, N, PAR, RV); value M, N; integer M, N; array PAR, RV; begin LOC FUNCT:= TEST THF:= FUNCT(M, N, PAR, RV) and TEST THF; FEVAL:= FEVAL + 1 end LOC FUNCT; ITMAX:= FEVALMAX:= IN[5]; AUX[2]:= N * IN[0]; TIM:= IN[7]; RELTOLPAR:= IN[1] ** 2; ABSTOLPAR:= IN[2] ** 2; ABSTOLRES:= IN[4] ** 2; INRMAX:= IN[6]; DUPVEC(1, N, 0, PR, PAR); if M < N then for I:= 1 step 1 until N do JAC[M + 1, I]:= 0; TEXT:= 4; MIT:= 0; TEST THF:= true; RES2:= STAP:= OUT[5]:= OUT[6]:= OUT[7]:= 0; FUNCT(M, N, PAR, FU2); RN:= VECVEC(1, M, 0, FU2, FU2); OUT[3]:= SQRT(RN); FEVAL:= 1; DAMPING ON:= false; for IT:= 1, IT + 1 while IT <= ITMAX and FEVAL < FEVALMAX do begin OUT[5]:= IT; JACOBIAN(M, N, PAR, FU2, JAC, LOCFUNCT); if not TEST THF then begin TEXT:= 7; "GO TO" FAIL end; LSQORTDEC(JAC, M, N, AUX, AID, CI); if AUX[3] ^= N then begin TEXT:= 5; "GO TO" FAIL end; LSQSOL(JAC, M, N, AID, CI, FU2); DUPVEC(1, N, 0, SOL, FU2); STAP:= VECVEC(1, N, 0, SOL, SOL); RHO:= 2; NORMX:= VECVEC(1, N, 0, PAR, PAR); if STAP > RELTOLPAR * NORMX + ABSTOLPAR or IT = 1 and STAP > 0 then begin for INR:= 0, INR + 1 while if INR = 1 then DAMPING ON or RES2 >= RN else not CONV and (RN <= RES1 or RES2 < RES1) do begin comment DAMPING STOPS WHEN R0 > R1 "AND" R1 <= R2 (BEST RESULT IS X1, R1) WITH X1 = X0 + I * DX, I:= 1, .5, .25, .125, ETC. ; RHO:= RHO / 2; if INR > 0 then begin RES1:= RES2; DUPVEC(1, M, 0, RV, FU2); DAMPING ON:= INR > 1 end; for I:= 1 step 1 until N do PR[I]:= PAR[I] - SOL[I] * RHO; FEVAL:= FEVAL + 1; if not FUNCT(M, N, PR, FU2) then begin TEXT:= 6; "GO TO" FAIL end; RES2:= VECVEC(1, M, 0, FU2, FU2); CONV:= INR >= INRMAX end DAMPING OF STEP VECTOR; if CONV then begin comment RESIDUE CONSTANT; MIT:= MIT + 1; if MIT < TIM then CONV:= false end else MIT:= 0; if INR > 1 then begin RHO:= RHO * 2; ELMVEC(1, N, 0, PAR, SOL, - RHO); RN:= RES1; if INR > 2 then OUT[7]:= IT end else begin DUPVEC(1, N, 0, PAR, PR); RN:= RES2; DUPVEC(1, M, 0, RV, FU2) end; if RN <= ABSTOLRES then begin TEXT:= 1; ITMAX:= IT end else if CONV and INRMAX > 0 then begin TEXT:= 3; ITMAX:= IT end else DUPVEC(1, M, 0, FU2, RV) end ITERATION WITH DAMPING AND TESTS else begin TEXT:= 2; RHO:= 1; ITMAX:= IT end end OF ITERATIONS; LSQINV(JAC, N, AID, CI); for I:= 1 step 1 until N do begin JJINV[I,I]:= JAC[I,I]; for J:= I + 1 step 1 until N do JJINV[I,J]:= JJINV[J,I]:= JAC[I,J] end CALCULATION OF INVERSE MATRIX OF NORMAL EQUATIONS; FAIL : OUT[6]:= SQRT(STAP) * RHO; OUT[2]:= SQRT(RN); OUT[4]:= FEVAL; OUT[1]:= TEXT; OUT[8]:= AUX[3]; OUT[9]:= AUX[5] end GSSNEWTON; eop