1

我已经有一段时间没有出现在这个非常好的论坛上了。我正在学习数值分析课程,我被要求编写二分法,这是我的代码

/*
 * Bisection.cpp
 *
 *  Created on: 08/10/2012
 *  Author: BRabbit27
 *  École Polytechnique Fédérale de Lausanne - M.Sc. CSE
 */

#include <cmath>
#include <iostream>

using namespace std;

double functionA(double x) {
    return sin(2.0 * x) - 1.0 + x;
}

double functionB(double x) {
    return 2.0 * x / (1.0 + x / 1.5);
}

double bisectionMethod(double (*function)(double), double a, double b, double tolerance) {
    double x;
    double f;
    double error = tolerance + 1;
    int step = 0;
    double fa = (*function)(a);
    double fb = (*function)(b);
    //Check the conditions of a root in the given interval
    if (a < b) {
        if (fa * fb < 0) {
            while (error > tolerance) {
                step++;

                x = (a + b) / 2.0;
                f = (*function)(x);

                if (f == 0) {
                    cout << "Root found in x = " << x;
                    return x;
                } else if (f * fa > 0) {
                    a = x;
                } else if (f * fa < 0) {
                    b = x;
                }
                error = (b - a) / pow(2.0, (double) step + 1);
            }
            cout << "Root found in x = " << x;
            return x;
        } else {
            cout << "There not exist a root in that interval." << endl;
            return -1;
        }
    } else {
        cout << "Mandatory \"a < b\", verify." << endl;
        return -1;
    }
}

int main(int argc, char *argv[]){
    bisectionMethod(functionA, -3.0, 3.0, 10.0e-7);
}

我唯一的问题是,当x = 0.354492时找到了根,而真正的根在x=1/3中,所以实际上我的双精度或容差不好。我不知道如何改进此代码以获得更好的结果。任何想法?

4

1 回答 1

2

真正的根不是x = 1/3!它是 3.52288

罪(2.0 * x) - 1.0 + x = 0

罪(2.0 * x)= 1 - x

1 - 1/3 = 2/3

罪(2/3)!= 2/3


你对宽容的定义让我觉得很奇怪。容忍度应该是x你愿意接受的范围吧?嗯,这很容易:b - a > tolerance

double bisectionMethod(double (*function)(double), double a, double b, double tolerance) {
    if (a > b) {
        cout << "Mandatory \"a < b\", verify." << endl;
        return -1;
    }
    double fa = function(a);
    double fb = function(b);
    if (fa * fb > 0) {
        cout << "No change of sign - bijection not possible" << endl;
        return -1;
    }

    do {    
        double x = (a + b) / 2.0;
        double f = (*function)(x);
 
        if (f * fa > 0) {
            a = x;
        } else {
            b = x;
        }
    } while (b - a > tolerance);

    cout << "Root found in x = " << x;
    return x;
}
于 2012-10-08T18:33:06.100 回答