Как применить метод Ньютона-Рафсона для поиска корней функции пятого числа

Описание

Я разработал алгоритм, реализующий метод Ньютона-Рафсона для нахождения корня функции пятого числа. Результат, который я должен отразить, это 303.6. Однако моя реализация не соответствует действительности.

Данные

Параметры

g = 9.81; 
Ds = 0.198; 
uj = 805.9; 
W = 0.0557;

Уравнение

0.024*((gDs/uj^2)^(1/3))(Y^(5/3)) + 0.2*(Y^(2/3)) - ((2.85/W)^(2/3)) = 0

Где производная для Y:

(0.04*d^(1/3)⋅g(1/3)⋅y^(2/3)) / u(2/3) + 2/15*y^(1/3)

Решение корня для Y

Код

    import java.lang.*;
    public class InvokeNewton {
    public static void main(String argv[]) {
    double del = 1e-5,
    double xx = 0 ;
    double dx =0, 
    double x= Math.PI/2;
    int k = 0;
    while (Math.abs(xx-x) > del && k<10 && f(x)!=0) {
      dx = f(x)/d(x);
      xx=x;
      x =x - dx;
      k++;

    System.out.println("Iteration number: " + k);
    System.out.println("Root obtained: " + x);
    System.out.println("Estimated error: " + Math.abs(xx-x));
    }
  }

    // Method to provide function f(x)

      public static double f(double x) {
        return 0.024*(Math.pow(g * Ds / Math.pow(uj, 2.0),(1.0/3.0)) *                       (Math.pow(Y,5.0/3.0))+ 0.2*(Math.pow(Y,2.0/3.0)) - (Math.pow((2.85/W)(2.0/3.0))));
      }

    // Method to provide the derivative f'(x).

      public static double d(double x) {
        return (0.04*Math.pow(Ds,1.0/3.0)*Math.pow(Y,2.0/3.0)) / Math.pow*uj,2.0/3.0) + 2 /    15*Math.pow(Y,1.0/3.0);
      }

}

Выход

    Iteration number: 1
Root obtained: 3.65373153496716
Estimated error: 2.0829352081722634
Iteration number: 2
Root obtained: 5.2246000232674215
Estimated error: 1.5708684883002615
Iteration number: 3
Root obtained: 6.618389759316356
Estimated error: 1.3937897360489346
Iteration number: 4
Root obtained: 7.906164279270034
Estimated error: 1.287774519953678
Iteration number: 5
Root obtained: 9.119558352547333
Estimated error: 1.213394073277299
Iteration number: 6
Root obtained: 10.27633029334909
Estimated error: 1.1567719408017574
Iteration number: 7
Root obtained: 11.387769167896339
Estimated error: 1.1114388745472485
Iteration number: 8
Root obtained: 12.461641418739712
Estimated error: 1.0738722508433725
Iteration number: 9
Root obtained: 13.503592201954325
Estimated error: 1.041950783214613
Iteration number: 10
Root obtained: 14.517895007865569
Estimated error: 1.0143028059112442

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

Справка

Метод Ньютона


person e.doroskevic    schedule 06.12.2013    source источник
comment
Какой результат вы получаете? Это может помочь помочь.   -  person prmottajr    schedule 06.12.2013
comment
Строка import java.lang.*; не нужна. Все в java.lang импортируется по умолчанию.   -  person Viktor Seifert    schedule 06.12.2013
comment
В функции d помогает ли замена 2 / 15 на 2.0/15.0?   -  person doctorlove    schedule 06.12.2013
comment
@doctorlove Я пробовал это раньше, это не имеет никакого значения, спасибо, что заметили это.   -  person e.doroskevic    schedule 06.12.2013
comment
Вы также можете изменить double x= Math.PI/2; удвоить x= Math.PI/2.0; Но опять же, публикация результата, который вы получаете, поможет.   -  person prmottajr    schedule 06.12.2013
comment
В вашем опубликованном коде много переменных, которые я не вижу, например uj, Ds,...   -  person doctorlove    schedule 06.12.2013
comment
@prmottajr я добавил вывод, который я получаю до 10-й итерации.   -  person e.doroskevic    schedule 06.12.2013
comment
@doctorlove переменные доступны в схеме проблемы. В целях тестирования я просто использовал значения этих переменных.   -  person e.doroskevic    schedule 06.12.2013
comment
Если вы знаете, что корень 303.6, почему бы не начать ближе - ваше решение только растет, и 10 шагов недостаточно.   -  person doctorlove    schedule 06.12.2013
comment
del это дельта? Возможно, размер вашего шага довольно мал.   -  person Viktor Seifert    schedule 06.12.2013
comment
@doctorlove Я хотел бы оставить его равным 1 или, возможно, использовать корень альфа / бета в качестве моего первоначального предположения, где альфа - это значение первого члена, а бета-значение третьего члена.   -  person e.doroskevic    schedule 06.12.2013
comment
@Viktor Seifert, это правильно, я использовал del для delta.   -  person e.doroskevic    schedule 06.12.2013
comment
Кажется, я приближаюсь к ответу. Я опубликую результат и исходный код метода Ньютона, когда закончу.   -  person e.doroskevic    schedule 06.12.2013
comment
Тогда я думаю, что вам нужно намного больше итераций, чем 10. Попробуйте пару сотен итераций и посмотрите, приближается ли результат к ожидаемому результату.   -  person Viktor Seifert    schedule 06.12.2013
comment
Ваша производная должна иметь y^-1/3 (обратите внимание на отрицательный результат).   -  person Teepeemm    schedule 06.12.2013


Ответы (2)


Производная Y^(2/3) равна 2/3*Y^(-1/3).

В общем, было бы намного более стабильным алгоритм, если бы сначала заменить X = Y ^ (1/3), чтобы у вас была обычная полиномиальная пятерка, решить ее для X и с результатом затем установить Y = X ^ 3.

A=0.024*(g*Ds/uj^2)^(1/3)), B=0.2, C= ((2.85/W)^(2/3))

0=F(X)=A*X^5+B*X^2-C
F'(X)=4A*X^4+2B*X

Для положительного X это хорошая выпуклая функция, идущая от отрицательных значений к положительным. X=0 — это минимум, поэтому не самое хорошее начальное значение. Если X=(C/A)^(1/5), то первый и последний члены в F(X) отменяются, а второй положителен, что указывает на хорошее начальное значение.

public class InvokeNewton {
  // Define the given constants
  static double g = 9.81, Ds = 0.198, uj = 805.9, W = 0.0557;
  // Combine the constants into to coefficients of the polynomial
  // Executed during class creation
  static double A=0.024*Math.pow(g*Ds/(uj*uj), 1.0/3), B=0.2, C= Math.pow(2.85/W,2.0/3);
  /* The original problem asks for the solution of 
   * A*y^(5/3)+B*y^(2/3)-C, introduce x=y^(1/3) <=> y=x^3
   * The equation to solve now is f(x)=Ax^5+Bx^2-C=0, 
   * a polynomial with derivative f'(x)=5Ax^4+2Bx
   */


  public static void main(String argv[]) {
    // set the precision with some buffer to the full 1e-16
    double del = 1e-12;
    double xx = 0 ;
    double dx =0;
    /* -------------------
     * initial point obtained by solving Ax^5-C=0
     * Now f(0)<0, f(x)>0 and f is convex monotone increasing
     * for x>0, so Newtons method produces a decreasing 
     * sequence of points quadratically converging to
     * the root of the equation.
     */
    double x= Math.pow(C/A,1.0/5);
    int k = 0;
    while (Math.abs(xx-x) > del && k<20) {
      // Standard Newton method
      dx = f(x)/d(x);
      xx=x;
      x =x - dx;
      k++;

      System.out.println("Iteration number: " + k);
      System.out.print("Root obtained: " + x);
      // the number of interest is y=x^3
      System.out.println(" solution Y="+Math.pow(x,3));
      System.out.println("Estimated error: " + Math.abs(xx-x));
    }
  }

  // Method to provide function f(x)

  public static double f(double x) {
    return A*Math.pow(x,5)+B*Math.pow(x,2)-C;
  }

  // Method to provide the derivative f'(x).

  public static double d(double x) { 
    return 5*A*Math.pow(x,4)+2*B*x;
  }

}

с результатами

 Iteration number: 1
 Root obtained: 7.127382504549578 solution Y=362.06804746760736
 Estimated error: 1.192264780096913
 Iteration number: 2
 Root obtained: 6.7530415629870015 solution Y=307.962806003808
 Estimated error: 0.37434094156257647
 Iteration number: 3
 Root obtained: 6.722154364534209 solution Y=303.75640454448063
 Estimated error: 0.030887198452792752
 Iteration number: 4
 Root obtained: 6.72196107618877 solution Y=303.73020272815836
 Estimated error: 1.932883454385248E-4
 Iteration number: 5
 Root obtained: 6.721961068677178 solution Y=303.7302017099299
 Estimated error: 7.511592237108289E-9
 Iteration number: 6
 Root obtained: 6.721961068677178 solution Y=303.7302017099299
 Estimated error: 0.0
person Lutz Lehmann    schedule 13.12.2013
comment
не могли бы вы прокомментировать свой код для большей ясности? Я рассмотрю его и оценю его для ответа. Ваши усилия очень ценятся. - person e.doroskevic; 13.12.2013

Код был соответствующим образом изменен. Код содержит комментарии для упрощения внесения изменений в ваш конкретный случай.

    public class RootFinder {

    // Method to be called to calculate the root 

    public void InvokeNewton(){
        double del = 1e-10;                                  // Delta or precision; 
        double xx = 0 ;                                      // Storage for previous root

        double dx = 0;                                       // Storage to hold a derivative of a pre-defined function;
        double x= Math.sqrt(beta/alpha);                     // Initial guess;

        while (Math.abs(xx-x) >= del &&  f(x)!= 0.0) {       // Math.abs(xx-x) - Estimated error;
            dx = f(x)/d(x);                                  // Derivative
            xx=x;                                            // New xx value set to previous root for compersion;
            x = x - dx;                                      // Root obtained;
        }
    }// Method to provide function f(x)
                          // !!! Provide your function bellow !!!
    public double f(double x) {
        return (alpha * (Math.pow(x,5.0/3.0)) + 0.2*(Math.pow(x,2.0/3.0)) - beta);
    }

    // Method to provide the derivative f'(x).
                          // !!!Provide the prime derivative of your function bellow!!!
    public double d(double x) {
        return (0.04*Math.pow(jf.getSourceDiameter(),1.0/3.0)*Math.pow(x,2.0/3.0)) / Math.pow(jf.getJetVelocity(),2.0/3.0) + 2.0 /  15.0*Math.pow(x,1.0/3.0);
    }
}
person e.doroskevic    schedule 06.12.2013