用于 SIN(X) + N 和 COS(X) + N 的模板化 Newton-Raphson 算法，其中 N 在 [0, 2] 范围内

• 下载源代码 - 14.3 KB

引言

您好，感谢您的光临。如果您在这里，您就知道 Newton-Raphson 算法的作用。简而言之，它找到 f(x) 形式方程的解

f(x) = ax + c

通过迭代使用导数，直到结果在 y 轴上达到 0。有一个 x 的初始值来启动迭代（种子）。在这种情况下，这些函数是 SIN 和 COS，只是为了让它更有趣（而且我认为更简单）。如果您想了解更多关于 Newton-Raphson 的信息，网上有很多信息。请使用您最喜欢的 Google 搜索引擎搜索它！ ;）。

背景

这段代码是作为 C++ 模板测试提供给我的，也用来展示对基本微积分的一些理解，我想也许有人会觉得它有用。

经过一些研究，我发现用于启动迭代的“种子”的选择必须足够接近所需的根。对于非三角函数，建议选择一个接近导致 f(x) = 0 的值的“种子”。

然而，对于三角函数，它们的周期性将使它们始终在 0 到 2PI 之间变化。换句话说，迭代的结果（下一次迭代的种子）必须减小到圆周上的一个值（以弧度为单位），以确定我们是否在根上或足够接近根。

找到函数的根

为了知道迭代的结果是否近似于 0，我们将相应的三角函数（SIN 或 COS）应用于迭代的结果，以确定该值是否为 0 或接近 0。此值存储在变量 *delta* 中。

如果 *delta* 在接近 0 的“可接受”值范围内，则可以停止迭代。为了实现这一点，需要设置一个可接受的最小误差阈值。选择 1E-5 作为停止迭代的可接受误差阈值。

因此，如果 *delta* 小于误差，则停止迭代并显示结果，从而避免对于不收敛于 0 的函数的无限迭代。设置了最大迭代次数（由静态变量 size_t MAXITER 给出）。否则，对于非收敛图，程序可能会进入无限循环而永远找不到结果。

使用代码

这是 main() 中代码的示例

switch(n) 
{ 
    default:
    case 0:
    {
        NR <Sine> ns(x);
        ns.Process();
        break;
    }
    case 1:
    {
        NR <Sine,1> ns(x);
        ns.Process();
        break;
    }
    case 2:
    {
        NR <Sine,2> ns(x);
        ns.Process(); 
        break;
    }
}

正如您所看到的，要实例化的模板需要在运行时才知道。换句话说，如果您尝试这样做

NR <Sine, n> ns(x);

...编译器会抱怨。如果有人知道如何实现这一点，我乐于接受建议并学习。

执行计算的代码是

// interface class to homogenize the 
// behavior of our  function classes
struct Funct
{
    virtual double Function( double x, int n = 0) = 0;
    virtual bool CalculateDiff(double fX) = 0;
    template < double f(double) >
    bool CalculateDiff(double fX)
    {
        double delta(f(fX));
        // f can be SIN or COS functions from math.h header
#ifdef _DEBUG
        std::cout << "Delta : " << delta << std::endl;
#endif
        if (0.0 > delta)
            delta *= -1.0;
        return delta < Error;
    }

};
// class Sine: performs calculations for SIN(X) function
struct Sine : public Funct
{
    // Function sin
    double Function( double x, int n = 0);
    bool CalculateDiff(double fX);
};
// class Cosine: performs calculations for COS(X) function
struct Cosine : public Funct
{
    // Function sin
    double Function( double x, int n = 0);
    bool CalculateDiff(double fX);
};

正如您所看到的，我们为每个函数使用一个类，使用相同的接口。为了强制执行此操作，我们从一个抽象类（或协议类，或接口）派生函数类。我们的类函数的实现如下

// Newton-Raphson for sin
double Sine::Function( double x, int n)
{
    return (sin(x) + n)/cos(x);
}
////////////////////////////////////////////////////////
// Calculates differences for sine
bool Sine::CalculateDiff(double fX)
{
    return Funct::CalculateDiff<sin>(fX);
}
////////////////////////////////////////////////////////
// Newton-Raphson for cosine
double Cosine::Function( double x, int n)
{
    return (cos(x)+n)/-sin(x);
}
////////////////////////////////////////////////////////
// Calculate differences for Cosine
bool Cosine::CalculateDiff(double fX)
{
    return Funct::CalculateDiff<cos>(fX);
}

最后，允许我们利用上述所有魔术的包装器是

// template class to work as a wrapper
// for the NewtonRapshon algorithm
template < class T , int N = 0>
class NR
{
    double _x; // seed to start iterations from. 
    double _fX; // function to approximate to 0. 
    // to mark if the iterations where successful
    // after the max number of iterations
    bool _success;
    size_t _i; // iterator

    // Copy constructor and assignment made
    // private to disallow copy and assignment
    NR(const NR& nr);
    NR& operator=(const NR& nr);
    // methods
    inline void PrintCalculations();
    inline void PrintResults();
    inline void Function();
public:
    static size_t MAXITER; // Maximum number of iterations
    explicit NR(double x) : _x(x), _fX(0.0) , _success(false), _i(0)
    {
    }
    void Process();
};

template  < class T, int N > size_t NR < T, N > ::MAXITER(50);

// execute function
template  < class T , int N>
void NR<T, N>::Function()
{
    T t;
    _fX = _x - t.Function(_x,N);
    _success = t.CalculateDiff(_fX);
}
// Print results
template  < class T, int N >
void NR<T, N>::PrintResults()
{
    using namespace std;
    if (_success)
        cout << "The root is: [" << _fX << "]. Result reached in " 
             << _i+1 << " Iterations." << endl;
    else
        cout << "** Result NOT reached in " << _i 
             << " Iterations." << endl;
}

template  < class T, int N >
void NR<T,N>::PrintCalculations()
{
    std::cout << _i+1 << "] x: " << _x << " - (" << sin(_x) 
              << "/" << cos(_x) << ") = f(x) : " << _fX << std::endl;
}
// Executes algorithm
// Iterates until: 
// a) result is found 
// or 
//b) MAXITER has been reached without a result found
template  < class T, int N >
void NR<T,N>::Process()
{
    for (; _i < NR::MAXITER ; _i++)
    {
        Function();
#ifdef _DEBUG
        PrintCalculations();
#endif
        if (_success)
        {
            break;
        }
        _x = _fX;
    }
    PrintResults();
}

包含的代码已使用 MS Visual Studio 2005 和 GNU C++ 编译器 3.4.4 编写和测试。要在 Linux 上编译，请使用

c++ -D_DEBUG  -Wall -pedantic -o NewRaph NewRaph.cpp StrHelper.cpp

...如果你想看到迭代和 *delta*，或者

c++ -Wall -pedantic -o NewRaph NewRaph.cpp StrHelper.cpp 

...如果你不想。

示例执行

下面是一个不带参数执行程序的示例

NewRaph
Usage: NewRaph SEED {S|C|B} [SHIFT] where:
SEED: Integer. Seed to start the iterations from.
{S|C|B}: S to calculate SIN(x), C for COS(X), B for Both
[SHIFT]: n value to establish the function f(x) + n. n varies from 0 to 2

关注点

对于 SIN(X) 和 COS(X)，该算法非常快地收敛到 f(x) = 0。对于 SIN(X) + N 和 COS(X) + N，该算法收敛到 -1。

可能有一种更好的方法来使用模板来实现这一点（使用模板元编程或类似的东西），但我不知道。欢迎反馈。

历史

• 版本 1。