2015-02-07 145 views
0

我一直在问,找到下面的功能功能的C++查找根

sin((a*x/(1 + pow(x, 2))) + 1) * atan(b*x - 1/2) + exp(-c*x) * atan(x) 

的根源有两个设定值

  • a=10b=2c=0
  • a=4.5b=2.8的,和c=1

但我没有给开始结束值,我需要找到根。我应该如何继续?

注意:atan()表示tan的反函数。

代码片段:

double f(double x, double a, double b, double c) 
{ 
    return sin((a*x/(1 + pow(x, 2))) + 1) * atan(b*x - 1/2) + exp(-c*x) * atan(x); 
} 

double RootFinder(double f(double, double, double, double), double a, double b, double c, double left, double right, double precision) 
{ 
    double f_left = f(left, a, b, c), now = left + precision, f_right = f(now, a, b, c); 
    while (f_left * f_right > 0 && now < right) 
    { 
     f_left = f_right; 
     now += precision; 
     f_right = now; 
    } 
    return now - precision/2; 
} 
+0

这是一个相当数学问题,我认为它属于另一个网站。但是,并非所有的查找算法都需要一个包含根的区间。所以,要么尝试其他算法,要么挑选一些随机值,如果公式给出不同的符号,则可以平分。 – Axel 2015-02-07 07:50:49

+0

不会。即使如果它似乎是一个数学问题,我被要求用C++实现。所以,我认为它属于这里。你能帮我结构吗? – tmgr 2015-02-07 07:51:56

+0

我不同意。你正在寻找一种算法,你根本没有代码。你必须先选择算法,然后从这里开始寻找:http://math.stackexchange.com/search?q=root+finding。如果您尝试实现算法并陷入困境,那么您可以回到这里,发布您的代码并获得答案。 – Axel 2015-02-07 07:59:30

回答

1

有你实现你的函数的一个bug。

atan(b*x - 1/2) 

术语1/2做整数除法和计算结果为0,这可能不是你想要的。通常,在使用双变量进行算术时使用双字面值。 pow()功能确实需要(double, int)作为它的一个重载,所以你很好。它也有一个(double, double)重载,但如果你的指数实际上是一个整数,那么你不希望这样。

下面是最简单的根查找方法的简单实现 - 平分法(后来我注意到OP使用平分标记,完美)。

#include <iostream> 
#include <cmath> 
#include <random> 

double f(const double x, const double a, const double b, const double c) 
{ 
    return sin((a*x/(1.0 + pow(x, 2))) + 1.0) * atan(b*x - 1.0/2.0) + exp(-c*x) * atan(x); 
} 

double BisectionMethod(
    double f(double, double, double, double), 
    const double a, const double b, const double c, 
    const std::random_device::result_type entropy) 
{ 
    std::mt19937 gen(entropy); 
    static const auto lower_bound = -1.0; 
    static const auto upper_bound = 1.0; 
    std::uniform_real_distribution<> dis(lower_bound, upper_bound); 

    auto pos_pt = dis(gen); 
    auto neg_pt = dis(gen); 

    while (f(pos_pt, a, b, c) < 0.0) 
     pos_pt = dis(gen); 

    while (f(neg_pt, a, b, c) > 0.0) 
     neg_pt = dis(gen); 

    static const auto about_zero_mag = 1E-8; 
    for (;;) 
    { 
     const auto mid_pt = (pos_pt + neg_pt)/2.0; 
     const auto f_mid_pt = f(mid_pt, a, b, c); 
     if (fabs(f_mid_pt) < about_zero_mag) 
      return mid_pt; 

     if (f_mid_pt >= 0.0) 
      pos_pt = mid_pt; 
     else 
      neg_pt = mid_pt; 
    } 
} 

int main() 
{ 
    double a, b, c; 
    std::random_device rd; 
    static const auto entropy = rd(); 

    a =10, b = 2.0, c = 0.0; 
    const auto root1 = BisectionMethod(f, a, b, c, entropy); 
    std::cout << "a = " << a << ", b = " << b << ", c = " << c << std::endl; 
    std::cout << "Found root: (" << root1 << ", " << f(root1, a, b, c) << ")" << std::endl; 

    a =4.5, b = 2.8, c = 1.0; 
    const auto root2 = BisectionMethod(f, a, b, c, entropy); 
    std::cout << "a = " << a << ", b = " << b << ", c = " << c << std::endl; 
    std::cout << "Found root: (" << root2 << ", " << f(root2, a, b, c) << ")" << std::endl; 
} 

Output

g++ -O3 -std=c++11 -Wall -Wextra -pedantic main.cpp -o root && ./root 
a = 10, b = 2, c = 0 
Found root: (0.143042, -2.12425e-09) 
a = 4.5, b = 2.8, c = 1 
Found root: (0.136172, 5.81247e-09) 

输出将每次运行时改变,因为它使用一个RNG。在视觉上,输出看起来是正确的

该代码假定根以-1.0和1.0为界,在您的情况下为true。如果你希望它更一般,那么你需要添加逻辑来处理溢出和检查nans。如果根不在-1.0和1.0之间,这将永远循环。尽管如此,它解决了这个问题中的具体问题,并且是更一般化的开始。

另请注意,您的函数有多个根,并且给定的代码只能找到一个根。

编辑:清理了代码。作为BisectionMethod()的参数添加了entropy,以便重复性好,如果我们正在谈论数值方法,这似乎是可取的。