Цикл программы метода Ньютона (на C) работает бесконечно

Этот код (прикрепленный к сообщению) на C использует метод Ньютона-Рафсона для поиска корней многочлена на определенном интервале.

Этот код отлично работает для некоторых полиномов, таких как x^3 + x^2 + x + 1, но время выполнения становится бесконечным для некоторых полиномов, таких как x^3 - 6*x^2 + 11*x - 6. То есть этот код отлично работает для многочленов, имеющих один или нулевой корень в введенном интервале, но если присутствует более одного корня, он работает бесконечно.

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

#include<stdio.h>
#include<math.h>
#include<stdlib.h>
#include<ctype.h>

int check(float num)                          //just a function to check for  the correct input
{
    char c;
    scanf("%c",&c);

    if(isalpha((int)c))
        printf("you entered an alphabet\n");

    else
        printf("you entered a character, please retry\n");

    return 0;
}

float func(float *p,int order,double x)                     //calculates the value of the function required in the formmula in main
{
    double fc=0.0;
    int i;

    for(i=0;i<=order;i++)
    {
        fc=fc+(double)(*p)*pow(x,(double)i);
        p++;
    }

    return fc;
}

float derv(float *q,int order,double x)               //calculates the derivative of the function required in the formmula in main
{     
    double dv=0.0,i;

    for(i=1;i<=order;i++)
    {
        dv=dv+(double)(*q)*(pow(x,(double)(i-1)))*(double)i;
        q++;
    }

    return dv;
}


int main()
{
    float coeff[1000];
    int order,count=0,i,j=0;
    char ch;
    float a,b;
    double val[5];

    printf("roots of polynomial using newton and bisection method\n");
    printf("enter the order of the equation\n");

    while(scanf("%d",&order)!=1)
    {
        printf("invalid input.please retry\n");
        while(getchar()!='\n'){}          
    }     

    printf("enter the cofficients\n");

    for(i=0;i<=order;i++)
    {
        printf("for x^%d  :",order-i);
        printf("\n");

        while(scanf("%f",&coeff[i])!=1)
        {
            check(coeff[i]);
        }   
    }

    while(getchar()!='\n'){}                                 //this clears the buffer accumulated upto pressing enter

    printf("the polynomial you entered is :\n");

    for(i=0;i<=order;i++)
    {
        printf(" %fx^%d ",coeff[i],order-i);
    }

    printf("\n");

    //while(getchar()!='\n'){};

    /* fflush(stdout);
    fflush(stdin);*/

    printf("plese enter the interval domain [a,b]\n");
    printf("enter a and b:\n");
    scanf("%f %f",&a,&b);

    while(getchar()!='\n'){}

    printf("the entered interval is [%f,%f]",a,b);

    fflush(stdout);
    fflush(stdin);

    //this array is used to choose a different value of x to apply newton's formula recurcively in an interval to scan it roperly for 3 roots

    val[0]=(double)b;       
    val[1]=(double)a;
    val[2]=(double)((a+b)/2);

    double t,x=val[0],x1=0.0,roots[10];

    while(1)
    {

        t=x1;
        x1=(x-(func(&coeff[0],order,x)/derv(&coeff[0],order,x)));           //this is the newton's formula

        x=x1;

        if((fabs(t - x1))<=0.0001 && count!=0)
        {
            roots[j]=x;
            j++;
            x=val[j];   // every time a root is encountered this stores the root in roots array and runs the loop again with different value of x to find other roots
            t=0.0;
            x1=0.0;
            count=(-1);

            if(j==3)
                break;
        }

        count++;
    }

    printf("the roots are = \n");

    int p=0;

    for(j=0;j<3;j++)
    {
        if(j==0 && roots[j]>=a && roots[j]<=b)
        {
            printf("  %f  ",roots[j]);
            p++;
        }

        if(fabs(roots[j]-roots[j-1])>0.5 && j!=0 && roots[j]>=a && roots[j]<=b)
        {
            printf(" %f  ",roots[j]);
            p++;
        }
    }

    if(p==0)
        printf("Sorry,no roots are there in this interval \n");

    return 0;
}

person Priyanshu Khandelwal    schedule 26.11.2016    source источник
comment
В статье Википедии о методе Ньютона описаны причины, по которым он может не сходится. Предположительно, вы столкнулись с некоторыми из этих причин.   -  person Jonathan Leffler    schedule 27.11.2016


Ответы (1)


Вы неправильно вычисляете функцию или производную, потому что сохраняете коэффициенты в обратном порядке, но не учитываете это.

Когда вы распечатываете уравнение, вы учитываете его, печатая order-i:

printf(" %fx^%d ",coeff[i],order-i);

Итак, вам нужно сделать то же самое в func:

fc=fc+(double)(*p)*pow(x,(double)(order-i));

и derv:

dv=dv+(double)(*q)*(pow(x,(double)((order-i)-1)))*(double)(order-i);

Причина, по которой он работает для многочленов, таких как x^3 + x^2 + x + 1, заключается в том, что в этом примере все коэффициенты одинаковы, поэтому не будет иметь значения, читаете ли вы массив вперед или назад.

Кроме того, вам может потребоваться учитывать другие причины, по которым метод может не сойтись, как упоминал Джонатон Леффлер в комментарии. Вы можете установить максимальное количество итераций для своего цикла и просто прерваться, если оно превысит максимум.

Хороший способ отладки чего-то подобного (помимо использования отладчика, конечно) — добавить несколько дополнительных операторов printf для отображения вычисляемых значений. Вы можете проверить работоспособность вывода, введя уравнение в поиск Google, это даст интерактивный график функции.

person samgak    schedule 26.11.2016