code 34440; procedure MARQUARDT(M,N,PAR,G,V,FUNCT,JACOBIAN,IN,OUT); value M,N; integer M,N; array PAR,G,V,IN,OUT; boolean procedure FUNCT; procedure JACOBIAN; begin integer MAXFE,FE,IT,I,J,ERR; real VV,WW,W,MU,RES,FPAR,FPARPRES,LAMBDA,LAMBDAMIN, P,PW,RELTOLRES,ABSTOLRES; array EM[0:7],VAL,B,BB,PARPRES[1:N],JAC[1:M,1:N]; procedure LOCFUNCT(M,N,PAR,G); value M, N; integer M,N; array PAR,G; begin FE:= FE+1; if FE >= MAXFE then ERR:= 1 else if not FUNCT(M,N,PAR,G) then ERR:= 2; if ERR^=0 then goto EXIT end LOCFUNCT; VV:=10; W:=0.5; MU:= 0.01; WW:=(if IN[6]<"-7 then "-8 else "-1*IN[6]); EM[0]:=EM[2]:=EM[6]:=IN[0]; EM[4]:=10*N; RELTOLRES:=IN[3]; ABSTOLRES:=IN[4]**2; MAXFE:=IN[5]; ERR:= 0; FE:= IT:= 1; P:=FPAR:= RES:= 0; PW:=-LN(WW*IN[0])/2.30; if not FUNCT(M,N,PAR,G) then begin ERR:= 3; goto ESCAPE end; FPAR:= VECVEC(1,M,0,G,G); OUT[3]:=SQRT(FPAR); for IT:= 1, IT+1 while FPAR > ABSTOLRES and RES > RELTOLRES*FPAR+ABSTOLRES do begin JACOBIAN(M,N,PAR,G,JAC,LOCFUNCT); I:=QRISNGVALDEC(JAC,M,N,VAL,V,EM); if IT=1 then LAMBDA:= IN[6] * VECVEC(1,N,0,VAL,VAL) else if P =0 then LAMBDA:= LAMBDA*W else P:= 0; for I:=1 step 1 until N do B[I]:=VAL[I]*TAMVEC(1,M,I,JAC,G); L: for I:=1 step 1 until N do BB[I]:=B[I]/(VAL[I]*VAL[I]+LAMBDA); for I:=1 step 1 until N do PARPRES[I]:= PAR[I] - MATVEC(1,N,I,V,BB); LOCFUNCT(M,N,PARPRES,G); FPARPRES:= VECVEC(1,M,0,G,G); RES:=FPAR-FPARPRES; if RES < MU * VECVEC(1,N,0,B,BB) then begin P:= P+1; LAMBDA:= VV * LAMBDA; if P=1 then begin LAMBDAMIN:= WW * VECVEC(1,N,0,VAL,VAL); if LAMBDA<LAMBDAMIN then LAMBDA:= LAMBDAMIN end; if P<PW then goto L else begin ERR:= 4; goto EXIT end; end; DUPVEC(1,N,0,PAR,PARPRES); FPAR:=FPARPRES end ITERATION; EXIT: for I:=1 step 1 until N do MULCOL(1,N,I,I,JAC,V,1/(VAL[I]+IN[0])); for I:=1 step 1 until N do for J:=1 step 1 until I do V[I,J]:= V[J,I]:= MATTAM(1,N,I,J,JAC,JAC); LAMBDA:= LAMBDAMIN:= VAL[1]; for I:= 2 step 1 until N do if VAL[I]>LAMBDA then LAMBDA := VAL[I] else if VAL[I]<LAMBDAMIN then LAMBDAMIN:= VAL[I]; OUT[7]:=(LAMBDA/(LAMBDAMIN+IN[0]))**2; OUT[2]:=SQRT(FPAR); OUT[6]:=SQRT(RES+FPAR)-OUT[2]; ESCAPE: OUT[4]:=FE; OUT[5]:=IT-1; OUT[1]:=ERR end MARQUARDT eop