用于计算逆的牛顿方法

Newton method for computing an inverse

本文关键字：方法计算用于更新时间：2023-10-16

这是遵循我在此线程中提出的问题：缺少vtable的链接错误

我定义了一个类"函数"和另外两个从"函数"继承的类"多项式"和"仿射"。

class function {
        public:
            function(){};
            virtual function* clone()const=0;
            virtual float operator()(float x)const=0; //gives the image of a number by the function
            virtual function* derivative()const=0;
            virtual float inverse(float y)const=0; 
            virtual ~function(){}
        };
        class polynomial : public function {
        protected:
            int degree;
        private:
            float *coefficient;
        public:
            polynomial(int d);
            virtual~polynomial();
            virtual function* clone()const;
            int get_degree()const;
            float operator[](int i)const; //reads coefficient number i
            float& operator[](int i); //updates coefficient number i
            virtual float operator()(float x)const; 
            virtual function* derivative()const;
            virtual float inverse(float y)const; 
        };
        class affine : public polynomial {
            int a; 
            int b; 
            //ax+b
    public:
            affine(int d,float a_, float b_);
            function* clone()const;
            float operator()(float x)const;
            function* derivative()const;
            float inverse(float y)const;
            ~affine(){}
        };

多项式中的方法逆似乎不能正常工作。它基于牛顿方法，该方法应用于固定 y（我们正在计算逆元素）和当前多项式 f 的函数 x->f（x）-y。

float polynomial::inverse(float y)const
{
    int i=0;
    float x0=1;
    function* deriv=derivative();
    float x1=x0+(y-operator()(x0))/(deriv->operator()(x0));
    while(i<=100 && abs(x1-x0)>1e-5)
    {
        x0=x1;
        x1=x0+(y-operator()(x0))/(deriv->operator()(x0));
        i++;
    }
    if(abs(x1-x0)<=1e-5) 
    {
        //delete deriv; //I get memory problems when I uncomment this line
        return x1;   
    }
    else    
    {
        cout<<"Maximum iteration reached in polynomial method 'inverse'"<<endl;
        //delete deriv; //same here
        return -1;
    }
}
double polynomial::operator()(double x)const
{
    double value=0;
    for(int i=0;i<=degree;i++) value+=coefficient[i]*pow(x,i);
    return value;
}
polynomial* polynomial::derivative()const
{
    if(degree==0)
    {
        return new affine(0,0,0);
    }
    polynomial* deriv=new polynomial(degree-1);
    for(int i=0;i<degree;i++)
        deriv[i]=(i+1)*coefficient[i+1];   
    return deriv;
}

我用 p：x->x^3 测试此方法：

#include "function.h"
int main(int argc, const char * argv[])
{
        polynomial p(3);
        for(int i=0;i<=2;i++) p[i]=0;
        p[3]=1;
        cout<<"27^(1/3)="<<p.inverse(27);

        return 0;
}

即使我输入 10，000 而不是 100，此脚本也会输出27^(1/3)=Maximum iteration reached in polynomial method 'inverse' -1。我在互联网上读过一些文章，似乎这是一种计算逆数的常用方法。

ABS函数原型是：int abs(int)因此，像abs(x1-x0)<=1e-5这样的测试不会像您预期的那样运行;您将整数与浮点数进行比较。在这种情况下，浮点数将转换为 int，因此它与 abs(x1-x0)<=0 相同

这可能就是您没有得到预期结果的原因 - 我建议再添加一些打印输出以深入了解事情。

嗯，问题出在方法"导数"上。我没有使用我重新定义的"operator[]"，而是使用了"->系数[]"，主脚本在p.inverse(27)中运行良好（只有 14 次迭代）。我只是用deriv->coefficient[i]=(i+1)*coefficient[i+1];替换了deriv[i]=(i+1)*coefficient[i+1];

检查此代码：

#include<iostream>
#include<cmath>
#include<math.h>
using namespace std;
void c_equation(int choose, double x);
void Processes(double x, double fx1, double fdx1, int choose);
void main()
{
    int choose,choose2;
    double x;
    system("color A");
    cout << "                                                             " << endl;
    cout << "=============================================================" << endl;
    cout << "Choose Equation : " << endl;
    cout << "_____________________________________" << endl;
    cout << "1- x-2sin(x)" << endl;
    cout << "2- x^2 + 10 cos(x)" << endl;
    cout << "3- e^x - 3x^2" << endl;
    cout << "                                                    " << endl;
    cin >> choose;
    cout << "If you have values  press 1/ random press 2 :" << endl;
    cin >> choose2;
    if (choose2 == 1)
    {
        cout << "                                         " << endl;
        cout << "Enter Xo : " << endl;
        cin >> x;
        c_equation(choose, x);
    }
    else if (choose2 == 2)
    {
        x = rand() % 20;
        cout << "Xo = " << x << endl;
        c_equation(choose, x);
        choose2 = NULL;
    }
    else
    {
        cout << "Worng Choice !! " << endl;
        choose = NULL;
        choose2 = NULL;
        main();
    }
}
void c_equation(int choose, double x)
{
    double fx;
    double fdx;
    double fddx;
    double result;
    if (choose == 1)
    {
        fx = x - 2 * sin(x);
        fdx = 1 - 2 * cos(x);
        fddx = 2 * sin(x);
        result = abs((fx * fddx) / pow(fdx, 2));
    }
    else if (choose == 2)
    {
        fx = pow(x, 2) + 10 * cos(x);
        fdx = 2 * x - 10 * sin(x);
        fddx = 2 - 10 * cos(x);
        result = abs((fx * fddx) / pow(fdx, 2));
    }
    else if (choose == 3)
    {
        fx = exp(x) - 3 * pow(x, 2);
        fdx = exp(x) - 6 * x;
        fddx = exp(x) - 6;
        result = abs((fx * fddx) / pow(fdx, 2));
    }
    else
    {
        cout << " " << endl;
    }
    //------------------------------------------------------------
    if (result < 1)
    {
        cout << "True Equation :) " << endl;
        Processes(x, fx, fdx , choose);
    }
    else
    {
        system("cls");
        cout << "False Equation !!" << endl;
        choose = NULL;
        x = NULL;
        main();
    }
}
void Processes(double x, double fx, double fdx , int choose)
{
    double xic;
    for (int i = 0; i < 3; i++)
    {
        xic = x - (fx / fdx);
        cout << "                                          " << endl;
        cout << "Xi = " << x << "   " << "F(Xi) = " << fx   << "   " << "  F'(Xi) = " << fdx   << "  " << "  Xi+1 = " << xic << endl;
        x = xic;
        if (choose == 1)
        {
            fx = xic - 2 * sin(xic);
            fdx = 1 - 2 * cos(xic);
        }
        else if (choose == 2)
        {
            fx = pow(xic, 2) + 10 * cos(xic);
            fdx = 2 * xic - 10 * sin(xic);
        }
        else if (choose == 3)
        {
            fx = exp(xic) - 3 * pow(xic, 2);
            fdx = exp(xic) - 6 * xic;
        }
    }
}