How-to apply Newton-Raphson method to find roots of a quintic function

Question

Description

I have developed an algorithm implementing Newton-Raphson method to find a root of a quintic function. The result which I must mirror is 303.6. However, my implementation fails to measure up.

Data

Parameters

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

Equation

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

Where derivative for Y:

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

Solving root for Y

Code

    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);
      }

}

Output

    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

The equation and the derivative has been checked and double-checked, however I am still not getting the desired output

Reference

Newton Method

Answer 1

The code has been modified accordingly. The code provides comments to simplify modification to your specific case.

    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);
    }
}

Answer 2

The derivative of Y^(2/3) is 2/3*Y^(-1/3).

In general it would make the algorithm much more stable to first substitute X=Y^(1/3) so that you have an ordinary polynomial quintic, solve that for X and with the result then set 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

For X positive this is a nice convex function going from negative to positive values. X=0 is a minimum, so not a good initial value. If X=(C/A)^(1/5), then the first and last term in F(X) cancel and the second is positive, which indicates a good starting value.

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;
  }

}

with results

 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