i dont know what happen..when i complie..there's error...
so..what suppose i do??please..
#include<dfp.h>
CMatrix CDfp::operator()(double(*p)(CMatrix&),CMatrix& X0,double e)
{
pfun=p,n=X0.getrow();
X.set(n,1),S.set(n,1);
X=X0, error=e;
CMatrix B(n,n),B0(n,n);
CMatrix TiDu(n,1),TiDu0(n,1);
CMatrix DG(n,1),DGT(1,n),DX(n,1),DXT(1,n);
CMatrix E(n,n);
CDiffer differ;
double a;
double dfm;
int flag(0);
int k(1),c(1);
for(int i=1;i<=n;i++)
B(i,i)=1;
B0=B;
TiDu=differ(pfun,X,error/10);
dfm=TiDu.getmod();
ofstream out(".//result.txt",ios::ate);
for(k;k<=n;k++)
{
S=(-1)*B*TiDu;
a=find_a(TiDu);
X=X+a*S;
TiDu0=TiDu, TiDu=differ(pfun,X,error);
dfm=TiDu.getmod ();
out<<"k="<<c++<<" "<<"[X]'=["<<X<<"] "<<"[TiDu]'=["
<<TiDu<<"] "<<"f(x)="<<pfun(X)<<" a="<<a<<endl;
if(dfm<error)
return X;
else
{
if(k=n)
k=1, B=B0;
else
{
DG=TiDu-TiDu0, DGT=DG.invert();
DX=a*S, DXT=DX.invert();
E=(DX*DXT)/DX.mult(DGT)-(B*DG*DGT*B)/DGT.mult(B*DG);
B=B+E, k+=1;
}
}
}
return X;
}
double CDfp::find_a(CMatrix TiDu)
{
double a=explore(1,TiDu);
return a;
}
double CDfp::f(double a)
{
return pfun(X+a*S);
}
double CDfp::explore0 (double a,double b)
{
double x1,x2;
x1=a+0.382*(b-a);
x2=a+0.618*(b-a);
while((x2-x1)>=error)
{
if(f(x1)<f(x2))
b=x2;
else
a=x1;
x1=a+0.382*(b-a);
x2=a+0.618*(b-a);
}
return f(x1)<f(x2) ? x1:x2;
}
double CDfp::explore(double a0,CMatrix TiDu)
{
double h;
double t=-2*pfun(X)/(TiDu.invert().mult(S));
if(t>0 && t<1/S.getmod())
h=t;
else
h=1/S.getmod();
double a1(a0),a2(a0+h);
double f1,f2;
f1=f(a1), f2=f(a2);
while(f2<f1)
{
h*=2, a1=a2, a2+=h;
f1=f2, f2=f(a2);
if(f2>f1)
return explore0(a1,a2);
}
h/=4;
a2=a1,a1-=h;
f2=f1,f1=f(a1);
while(f2>f1)
{
h*=2, a2=a1, a1-=h;
f2=f1, f1=f(a1);
}
return explore0(a1,a2);
}