SCIP in C++11 ― 1.3.4節 その2

Newton法と問題1.40

fixedPointはその1と同じなので略。

問題1.40のテストのひとつは、区間二分法でやった3次方程式と同じ。

-----

const function<double(double)>

differentiate(const function<double(double)>& g)

{

    const double dx(0.0001);

    return([dx,g](const double x){return((g(x+dx)-g(x))/dx);});

}

const function<double(double)>

NewtonTransform(const function<double(double)>& g)

{

    return([g](const double x){return(x-g(x)/differentiate(g)(x));});

}

const double NewtonsMethod

(const function<double(double)>& g,const double guess){

    return(fixedPoint(NewtonTransform(g),guess));

}

const double square(const double x){return(x*x);}

const double cube(const double x){return(x*x*x);}

const function<double(double)>

cubic(const double a,const double b,const double c){

    return([a,b,c](const double x){return(cube(x)+a*square(x)+b*x+c);});

}

const double sqrtNewton(const double x)

{

    return(NewtonsMethod([x](const double y){return(y*y-x);},1.0));

}

   

int main(int argc, char** argv)

{

    cout<<setprecision(16)<<"(d/dx)x^3(5)="<<differentiate(cube)(5)<<endl;

    cout<<setprecision(16)<<"sqrt(2)="<<sqrtNewton(2.0)<<endl;

    double a(0);double b(0);double c(-1);

    cout<<endl<<"Problem 1.40:"<<endl

        <<"zero of x^3+ax^2+bx+c (a="<<a<<", b="<<b<<", c="<<c<<")"<<endl

        <<"x="<<NewtonsMethod(cubic(a,b,c),1.0)<<endl;

    a=0;b=-2;c=-3;

    cout<<endl

        <<"zero of x^3+ax^2+bx+c (a="<<a<<", b="<<b<<", c="<<c<<")"<<endl

        <<"x="<<NewtonsMethod(cubic(a,b,c),1.0)<<endl;

    return(0);

}

----

出力

----

(d/dx)x^3(5)=75.00150000979033

sqrt(2)=1.41421356245306

Problem 1.40:

zero of x^3+ax^2+bx+c (a=0, b=0, c=-1)

x=1

zero of x^3+ax^2+bx+c (a=0, b=-2, c=-3)

x=1.893289196306116