Below is the syntax highlighted version of BohrRadius.java
from §9.6 Numerical Linear Algebra.
/****************************************************************************** * Compilation: javac BohrRadius.java * Execution: java BohrRadius * * Comptues the function f(r) = (1 - 3r/4 + r^2 / 8 - r^3 / 192)^2 e^(-r/2) * and its first and second derivatives. * The sample client runs Newton's method to finds the roots and local * minima and maxima. The Bohr radius function has three zeros at 1.872, * 6.611, and 15.518. These are also global minima. It has two local * maxima at 3.420 and 22.676. * * It illustrates the use of a "callback" using Java interfaces. * * % java BohrRadius * Possible roots: * f(1.8716444550481592) = 1.1817103296889357E-29 * f(6.610814578664537) = 6.529313277927562E-31 * f(6.6108145786645425) = 4.630205537810115E-31 * f(15.517540966287259) = 5.38926090372337E-33 * f(1.8716444550481608) = 9.872670880655912E-30 * * Possible optimum: * f(1.8716444550481757) = 3.191866873333672E-33 * f(3.419589235198309) = 0.017526400937842416 * f(6.610814578664558) = 1.6278066343863564E-32 * f(15.51754096628727) = 0.0 * f(22.676316353129153) = 0.0018496340643479514 * ******************************************************************************/ public class BohrRadius implements Function { // f(r) public double eval(double x) { double temp = 1 - 3*x/4 + x*x/8 - x*x*x/192; double exp = Math.exp(-x/2); return temp * temp * exp; } // f'(r) = first derivative public double deriv(double x) { double temp = 1 - 3*x/4 + x*x/8 - x*x*x/192; double temp1 = - 3.0/4 + x/4 - x*x/64; double exp = Math.exp(-x/2); return 2*temp*temp1*exp - temp*temp*exp/2; } // f''(r) = second derivative public double deriv2(double x) { double temp = 1 - 3*x/4 + x*x/8 - x*x*x/192; double temp1 = - 3.0/4 + x/4 - x*x/64; double temp2 = 1.0/4 - x/32; double exp = Math.exp(-x/2); return 2*temp1*temp1*exp - 2*temp*exp*temp1 + 2*temp*temp2*exp + temp*temp*exp/4; } // sample client public static void main(String[] args) { Function f = new BohrRadius(); double root, optimum; StdOut.println("Possible roots: "); root = Newton.root(f, 0.0); StdOut.println("f(" + root + ") = " + f.eval(root)); root = Newton.root(f, 4.0); StdOut.println("f(" + root + ") = " + f.eval(root)); root = Newton.root(f, 5.0); StdOut.println("f(" + root + ") = " + f.eval(root)); root = Newton.root(f, 13.0); StdOut.println("f(" + root + ") = " + f.eval(root)); root = Newton.root(f, 22.0); StdOut.println("f(" + root + ") = " + f.eval(root)); StdOut.println(); StdOut.println("Possible optimum: "); optimum = Newton.optimum(f, 0.0); StdOut.println("f(" + optimum + ") = " + f.eval(optimum)); optimum = Newton.optimum(f, 4.0); StdOut.println("f(" + optimum + ") = " + f.eval(optimum)); optimum = Newton.optimum(f, 5.0); StdOut.println("f(" + optimum + ") = " + f.eval(optimum)); optimum = Newton.optimum(f, 13.0); StdOut.println("f(" + optimum + ") = " + f.eval(optimum)); optimum = Newton.optimum(f, 22.0); StdOut.println("f(" + optimum + ") = " + f.eval(optimum)); StdOut.println(); } }