Hi, Isn't it possible to use BigDecimal or similar (e.g. http://www.apfloat.org)? I know that some code might need to be rewritten.
Cheers, Mikkel. 2010/8/15 Phil Steitz <[email protected]>: > Jbb wrote: >> Thanks Phil (and Luc). >> >> I understand the point, but is there any known workaround? >> Is there anyway to get more precision on the root finder solver taken into >> account that calculation time (speed of the algorithm) is not important in >> this case? >> Anyway to exchange speed for accuracy? > > The only thing that I can think of is to recast the function so that > its values near the root are distinguishable from zero. If the > function is returning values indistinguishable from zero near the > root, there is nothing that the solver can do to distinguish the one > you are looking for. I did try changing the function definition to > eliminate the pow and sqrt, but that did not help much. If you have > access to R, you can see the problem by executing the following: > > x <- seq( -2*pi - .00005, -2*pi + .00005, length = 100) > y <- x * x - ((1 - cos(x)) * (1 - cos(x))) - ((2*pi - sin(x)) * > (2*pi - sin(x))) > plot(x, y, axes=T, type="b", pch="*", xlab=" ", ylab=" ") > > The y vector is your function with sqrt and pow removed. You can > see from the plot and also by just examining the y values that the > function is not too well-behaved near the root and appears constant > at zero over (-2pi - .00001, -2pi + .00001). > > You can see the same thing by tabulating values using your Java code: > > AlfaFunctionClass alfaFunction = new AlfaFunctionClass(x, y, z); > final double answer = -2d * Math.PI; > final double diff = .00005; > final int count = 100; > final double delta = diff / count; > double root = answer - diff / 2; > for (int i = 0; i < count; i++) { > System.out.println( > "point: " + root + > " value: " + alfaFunction.value(root) + > " arg error: " + (root + 2*Math.PI) + > " close to zero: " + (Math.abs(alfaFunction.value(root)) <1E-300) > ); > root += delta; > } > > You will see that the above displays "true" for "close to zero" > throughout the interval described above. Replacing "<1E-300" with > "==0d" yields the same result. So bottom line is that from the > perspective of the JVM, the function vanishes identically over the > an interval of length approximately 1E-5 around -2pi. So the only > way you are going to be able to distinguish the "true" root > numerically is to transform the function so that values near the > root are distinguishable from zero. > > Phil > > >> >> Thanks again and best regards, >> >> Jbb >> >> El 15/08/2010, a las 05:47, Phil Steitz <[email protected]> escribió: >> >>> Juan Barandiaran wrote: >>>> Hi, I'm working in a theoretical physics problem in which I have to find >>>> the >>>> roots of a function for every point in x,y,z. >>>> I can bracket quite precisely where every root will be found and I know the >>>> exact solution for many points: >>>> For every point (x= 1, y= n*2*PI, z= 0) the root is in -n*2*PI. >>>> >>>> I have been using the newBrentSolver function, but to my surprise, even >>>> when >>>> I set the Accuracy settings to very >>>> small values and the MaxIterationCount to Integer.MAX_VALUE I get an error >>>> of 1.427E-5 (too much). >>>> >>>> The function to be solved is alfa in >>>> : >>>> alfa+Math.sqrt(Math.pow(x-Math.cos(alfa),2.0)+Math.pow(y-Math.sin(alfa),2.0)+z*z); >>>> >>>> Any hint of how could I get more precision on a root finder solver? >>>> I attach the code of a short Test Case program that shows the error I'm >>>> getting. >>> The problem is in the function evaluation. The function vanishes >>> numerically (i.e., returns a value with absolute value less than >>> 1E-90) throughout a neighborhood of width 1E-5 around the root you >>> are trying to find. The solver will (correctly, given its contract) >>> return the first value that it finds within this interval. >>> >>> Phil >>>> Thanks and best regards, >>>> >>>> JBB >>>> >>>> /***************************************************************************************************************************************** >>>> >>>> ELECTRON DENSITY OF ENERGY: ALFA Retarded Position Calculation >>>> SpinningParticles Model >>>> >>>> @author Martin Rivas (SpinningParticles theory) & Juan Barandiaran >>>> (Numerical Analysis & Calculation Program) >>>> @web http://spinningparticles.com/ >>>> >>>> Test case to try different solvers from the Apache Commons-Math library, >>>> showing 1.427E-5 (too much) error in the resulting roots. >>>> >>>> Compiled in JAVA with Eclipse SDK 3.6: >>>> http://www.eclipse.org/downloads/packages/eclipse-classic-360/heliosr >>>> and with the Apache Commons-Math library 2.1 < >>>> http://commons.apache.org/math/> >>>> >>>> ******************************************************************************************************************************************/ >>>> >>>> >>>> package testCaseRootFindingSolversPackage; >>>> >>>> import org.apache.commons.math.analysis.UnivariateRealFunction; >>>> import org.apache.commons.math.analysis.solvers.*; >>>> import org.apache.commons.math.ConvergenceException; >>>> import org.apache.commons.math.FunctionEvaluationException; >>>> >>>> public class TestCaseRootFindingSolvers { >>>> >>>> public static class AlfaFunctionClass implements UnivariateRealFunction { >>>> //Alfa is the angle of the point dependent retarded position where the >>>> electron charge was when it produced the field that influences the >>>> evaluated >>>> x,y,z point >>>> //To calculate it we must find the root of its equation between >>>> -sqr((x^2+y^2+z^2))-1 and -sqr((x^2+y^2+z^2))+1 as the correct root must be >>>> negative (retarded) >>>> //and bracketed around the module of the position vector >>>> >>>> private double x; >>>> private double y; >>>> private double z; >>>> >>>> public AlfaFunctionClass(double x, double y, double z) { >>>> this.x = x; >>>> this.y = y; >>>> this.z = z; >>>> } >>>> // this is the method specified by interface UnivariateRealFunction >>>> public double value(double alfa) { >>>> // here, we have to evaluate the function that calculates the retarded >>>> position angle alfa. >>>> return >>>> alfa+Math.sqrt(Math.pow(x-Math.cos(alfa),2.0)+Math.pow(y-Math.sin(alfa),2.0)+z*z); >>>> } >>>> } >>>> >>>> public static class EnergyDensityFunctionAlfaClass { >>>> public double evaluate(double x, double y, double z) { >>>> double alfa; >>>> UnivariateRealFunction alfaFunction = new AlfaFunctionClass(x,y,z); >>>> UnivariateRealSolverFactory factory = >>>> UnivariateRealSolverFactory.newInstance(); >>>> UnivariateRealSolver solver = factory.newBrentSolver(); //Fails when the >>>> result in the two bracketing points has the same sign >>>> // UnivariateRealSolver solver = factory.newBisectionSolver(); >>>> >>>> //Trying to get the maximal precision with the available settings >>>> solver.setAbsoluteAccuracy(1E-90); >>>> solver.setMaximalIterationCount(Integer.MAX_VALUE); >>>> solver.setRelativeAccuracy(1E-90); >>>> solver.setFunctionValueAccuracy(1E-90); >>>> try { >>>> alfa = solver.solve(alfaFunction, -1.0*Math.sqrt(x*x+y*y+z*z)-1.0, >>>> -1.0*Math.sqrt(x*x+y*y+z*z)+1.0); >>>> } catch (ConvergenceException e) { >>>> // TODO Auto-generated catch block >>>> e.printStackTrace(); >>>> alfa= 0.0; >>>> } catch (FunctionEvaluationException e) { >>>> // TODO Auto-generated catch block >>>> e.printStackTrace(); >>>> alfa= 0.0; >>>> } >>>> System.out.println("getFunctionValue="+String.valueOf(solver.getFunctionValue())+"??? >>>> Can't be cero with 1,42E-5 error?????"); >>>> System.out.println("getResult="+String.valueOf(solver.getResult())); >>>> System.out.println("getAbsoluteAccuracy="+String.valueOf(solver.getAbsoluteAccuracy())); >>>> >>>> System.out.println("getRelativeAccuracy="+String.valueOf(solver.getRelativeAccuracy())); >>>> System.out.println("getFunctionValueAccuracy="+String.valueOf(solver.getFunctionValueAccuracy())); >>>> return alfa; >>>> } >>>> } >>>> public static void main(String args[]) { >>>> // Variables >>>> double alfa=0.0; >>>> double error; >>>> EnergyDensityFunctionAlfaClass energyDensityFunction= new >>>> EnergyDensityFunctionAlfaClass(); >>>> //Example point, although we have to evaluate the function for every point >>>> in space. >>>> double x=1.0; >>>> double n=1.0; >>>> double y=n*2.0*Math.PI; >>>> double z=0.0; >>>> alfa= energyDensityFunction.evaluate(x, y, z); >>>> System.out.println(); >>>> System.out.println("Alfa= Root calculated with Apache Solver= >>>> "+String.valueOf(alfa)); >>>> error= alfa+n*2*Math.PI; >>>> System.out.println(); >>>> System.out.println("Error= Difference between the Alfa Root calculated with >>>> Apache Solver and n*2*PI (exact root for any integer n)= >>>> "+String.valueOf(error)); >>>> } // End of main >>>> } //End >>>> >>> >>> --------------------------------------------------------------------- >>> To unsubscribe, e-mail: [email protected] >>> For additional commands, e-mail: [email protected] >>> >> >> --------------------------------------------------------------------- >> To unsubscribe, e-mail: [email protected] >> For additional commands, e-mail: [email protected] >> > > > --------------------------------------------------------------------- > To unsubscribe, e-mail: [email protected] > For additional commands, e-mail: [email protected] > > --------------------------------------------------------------------- To unsubscribe, e-mail: [email protected] For additional commands, e-mail: [email protected]
