MUG nonlinear curve fit procedure: activate the cells inside

Date: 12-feb-1998 By: Preben Alsholm [ifakpa@pop.dtu.dk] MUG list member

> nonlinfit:=proc(Lx::list(numeric),Ly::list(numeric),eq::name=algebraic,x::name,initparam::set(name=numeric),u::{identical(iter)=posint,name},v::evaln)
local y,n,fx,f,j,N,param,param2,p,S,dp,i,k,pnum,Ls,parseq,df,chg,MAXSTEP; global maxstep;
n:=nops(Lx); if nops(Ly)<>n then ERROR(`The lists must have equal length`) fi;
y:=lhs(eq);
fx:=rhs(eq);
param:= op(select(type,indets(fx) minus {x},name));
param2:=map( t->lhs(t),initparam);
if {param} minus param2 <> {} then ERROR(`Provide starting values for all parameters`) fi;
pnum:=nops([param]);
f:=unapply(fx,param,x);
if nargs>5 and type(args[6],identical(iter)=posint) then N:=rhs(u) else N:=5 fi;
for k to pnum do df[k]:=D[k](f) od;
p[0]:=op(subs(initparam,[param]));
for j from 0 to N do
S:={seq( f(p[j],Lx[i])+add(df[k](p[j],Lx[i])*dp[k],k=1..pnum)=Ly[i],i=1..n)};
Ls:=linalg[leastsqrs](S,{seq(dp[k],k=1..pnum)});
for k in Ls do if type(rhs(k),name) then Ls:=Ls minus {k} union {lhs(k)=0} fi od;
if type(maxstep,procedure) then MAXSTEP:=signum*(maxstep@abs) else MAXSTEP:=s->s fi;
# maxstep determines the absolute value of the damping to be used
chg:=map(MAXSTEP,subs(Ls,[seq(dp[k],k=1..pnum)]));
p[j+1]:=op(chg+[p[j]]);
od;

> parseq:=seq( zip( (s,t)->s=t,[param],[p[j]]),j=0..N);
if nargs>5 and type(args[nargs],name) then if nargs=6 then u:=parseq else v:=parseq fi fi;
subs(parseq[N+1],eq)
end:

Go back