Метод Дюрана-Кернера для нахождения корней нелинейного уравнения

Меня просят найти корни f(x) = 5x(e^-mod(x))cos(x) + 1 . Ранее я использовал метод Дюрана-Кернера, чтобы найти корни функции x^4 -3x^3 + x^2 + x + 1 с кодом, показанным ниже. Я думал, что могу просто повторно использовать код для поиска корней f(x), но всякий раз, когда я заменяю x^4 -3x^3 + x^2 + x + 1 на f(x), программа выводит nan для всех корней. Что не так с моей реализацией Дюрана-Кернера и как мне изменить ее для работы с f(x)? Буду очень благодарен за любую помощь.

#include <iostream>
#include <complex>
#include <math.h>

using namespace std;

typedef complex<double> dcmplx;

dcmplx f(dcmplx x)
{
    // the function we are interested in
    double a4 = 1;
    double a3 = -3;
    double a2 = 1;
    double a1 = 1;
    double a0 = 1;

    return (a4 * pow(x,4) + a3 * pow(x,3) + a2 * pow(x,2) + a1 * x + a0);
}


int main()
{

dcmplx p(.9,2);
dcmplx q(.1, .5);
dcmplx r(.7,1);
dcmplx s(.3, .5);

dcmplx p0, q0, r0, s0;

int max_iterations = 100;
bool done = false;
int i=0;

while (i<max_iterations && done == false)
{
    p0 = p;
    q0 = q;
    r0 = r;
    s0 = s;


p = p0 - f(p0)/((p0-q)*(p0-r)*(p0-s));
q = q0 - f(q0)/((q0-p)*(q0-r)*(q0-s));
r = r0 - f(r0)/((r0-p)*(r0-q)*(r0-s0));
s = s0 - f(s0)/((s0-p)*(s0-q)*(s0-r));

    // if convergence within small epsilon, declare done
    if (abs(p-p0)<1e-5 && abs(q-q0)<1e-5 && abs(r-r0)<1e-5 && abs(s-s0)<1e-5)
        done = true;

    i++;
}

cout<<"roots are :\n";
cout << p << "\n";
cout << q << "\n";
cout << r << "\n";
cout << s << "\n";
cout << "number steps taken: "<< i << endl;

return 0;
}

Единственное, что я пока изменил, это функция dcmplx f. Я менял его на

dcmplx f(dcmplx x)
{
    // the function we are interested in
    double a4 = 5;

    double a0 = 1;

    return (a4 * x * exp(-x) * cos(x) )+ a0;
}

c++ nonlinear-functions newtons-method polynomial-math
person Ca01an    schedule 14.12.2014    source источник

Ответы (2)


arrow_upward
4
arrow_downward

Используемый вами метод Дюрана-Кернера требует, чтобы функция была непрерывной на интервале, в котором вы работаете.

Здесь мы имеем несоответствие между математическим представлением и ограничениями числовых приложений. Я бы предложил вам построить свою функцию (введя формулу в Google, вы, конечно, получите краткий обзор реальной части). Вы заметите, что:

  • существует бесконечность корней из-за периодичности косинуса.
  • из-за x * exp (-x) абсолютное значение быстро превышает максимальную точность, которую может удерживать число с плавающей запятой.

Чтобы понять последствия для вашего кода, я предлагаю вам проследить разные итерации. Вы заметите, что p, r и s сходятся очень быстро, в то время как q расходится (очевидно, по следу одного из огромных пиков):

  • На 2-й итерации q уже равно 1e74.
  • На 3-й итерации уже больше, чем может хранить двойник.
  • Поскольку q используется при расчете p, r и s, ошибка распространяется на другие члены
  • На 5-й итерации все термины находятся в NAN.
  • Затем он смело продолжается через 100 итераций.

Возможно, вы могли бы заставить его работать, выбрав разные отправные точки. Если нет, вам придется использовать какой-то другой метод и тщательно выбирать межстенные стены, над которыми вы работаете.

person Christophe    schedule 14.12.2014
comment
Любая идея, какой другой метод я мог бы использовать? - person Ca01an; 14.12.2014
comment
Я не математик. Но я бы попробовал простой метод Ньютона-Рафсона (Xi+1 = Xi - f(Xi)/f'(Xi)). Или метод Брента en.wikipedia.org/wiki/Brent%27s_method - person Christophe; 14.12.2014
comment
Я подтверждаю, что метод Брента с алгоритмом, представленным в Википедии, работает хорошо, пока вы можете дать точку выше и точку ниже корня. Я пробовал с двойными вместо сложных, и он нашел корень в 5 итерации - person Christophe; 15.12.2014

arrow_upward
2
arrow_downward

Вы должны были отметить в своей документации по методу Дюрана-Кернера (изобретенному Карлом Вейерштрассом около 1850 г.), что он применяется только к многочленам. Ваша вторая функция далека от многочлена.

Действительно, из-за функции mod ее приходится объявлять как неприятную функцию для численных методов. Большинство из них полагаются на непрерывность заданной функции, т. е. если значение близко к нулю, есть большая вероятность, что поблизости есть корень, а если знак меняется на интервале, то корень на интервале есть. Даже самые простые методы без производных, такие как метод деления пополам или метод Бренца на более сложном конце этого класса, предполагают наличие этих свойств.

person Lutz Lehmann    schedule 27.01.2015