Menu
▾
▴
[Cdk-commits] CVS: cdk/src/org/openscience/cdk/modeling/forcefield ConjugateGradientMethod.java,1.6,1.7 ForceField.java,1.6,1.7 GeometricMinimizer.java,1.8,1.9 LineSearch.java,1.3,1.4 PotentialFunction.java,1.3,1.4
Update of /cvsroot/cdk/cdk/src/org/openscience/cdk/modeling/forcefield In directory sc8-pr-cvs1.sourceforge.net:/tmp/cvs-serv5592/src/org/openscience/cdk/modeling/forcefield Modified Files: ConjugateGradientMethod.java ForceField.java GeometricMinimizer.java LineSearch.java PotentialFunction.java Log Message: Corrected ASCII Errors, rearranged methods in Geometric Minimizer Index: ConjugateGradientMethod.java =================================================================== RCS file: /cvsroot/cdk/cdk/src/org/openscience/cdk/modeling/forcefield/ConjugateGradientMethod.java,v retrieving revision 1.6 retrieving revision 1.7 diff -u -r1.6 -r1.7 --- ConjugateGradientMethod.java 18 Jan 2005 11:02:56 -0000 1.6 +++ ConjugateGradientMethod.java 21 Jan 2005 12:58:58 -0000 1.7 @@ -13,7 +13,7 @@ * */ public class ConjugateGradientMethod { - double µk = 0; + double uk = 0; GVector vk = new GVector(3); GVector vkminus1 = new GVector(3); @@ -33,17 +33,17 @@ * uk = gk gk / gk-1 gk-1 * */ - public void setµk(GVector xkminus1, GVector xk, PotentialFunction forceFieldFunction) { + public void setuk(GVector xkminus1, GVector xk, PotentialFunction forceFieldFunction) { GVector temporalVector = new GVector(forceFieldFunction.gradientInPoint(xk)); - µk = temporalVector.dot(temporalVector); + uk = temporalVector.dot(temporalVector); temporalVector.set(forceFieldFunction.gradientInPoint(xkminus1)); - µk = µk / temporalVector.dot(temporalVector); - //System.out.println("µk = " + µk); + uk = uk / temporalVector.dot(temporalVector); + //System.out.println("uk = " + uk); return; } /** - * vk=-gk + µk vk-1 + * vk=-gk + uk vk-1 * * @param gk Gradient at coordinates Xk * @param iterNumber Iteration number @@ -54,7 +54,7 @@ vkminus1.set(vk); vk.set(gk); vk.scale(-1); - vkminus1.scale(µk); + vkminus1.scale(uk); vk.add(vkminus1); //System.out.println("vector vk : " + vk); } Index: ForceField.java =================================================================== RCS file: /cvsroot/cdk/cdk/src/org/openscience/cdk/modeling/forcefield/ForceField.java,v retrieving revision 1.6 retrieving revision 1.7 diff -u -r1.6 -r1.7 --- ForceField.java 18 Jan 2005 10:57:31 -0000 1.6 +++ ForceField.java 21 Jan 2005 12:58:58 -0000 1.7 @@ -18,28 +18,28 @@ * Constructor for the ForceField object */ public ForceField() {} - - + + /** - * Read the 3d coordinates of atoms in a molecule from its atom container + * Get the 3xN coordinates vector of a molecule of N atoms from its atom container * *@param molecule molecule store in an AtomContainer - *@return GVector with 3N coordinates (N: atom numbers) + *@return GVector with 3xN coordinates (N: atom numbers) */ public GVector getCoordinatesVector(AtomContainer molecule) { - System.out.println("molecule: " + molecule.toString()); - System.out.println("Atoms number = " + molecule.getAtomCount()); + //System.out.println("molecule: " + molecule.toString()); + //System.out.println("Atoms number = " + molecule.getAtomCount()); GVector point0 = new GVector(3 * (molecule.getAtomCount())); - System.out.println("point0 = " + point0); - + //System.out.println("point0 = " + point0); + Atom thisAtom = new Atom(); int j = 0; for (int i = 0; i < molecule.getAtomCount(); i++) { thisAtom = molecule.getAtomAt(i); - System.out.println("thisAtom = " + thisAtom); - System.out.println("thisAtom.getPoint3d() = " + thisAtom.getPoint3d()); + //System.out.println("thisAtom = " + thisAtom); + //System.out.println("thisAtom.getPoint3d() = " + thisAtom.getPoint3d()); j = 3 * i; point0.setElement(j, thisAtom.getX3d()); @@ -47,40 +47,50 @@ point0.setElement(j + 2, thisAtom.getZ3d()); } - System.out.println("Atoms coordinates vector: " + point0); + //System.out.println("Atoms coordinates vector: " + point0); return point0; } - + + /** - * Read the coordinates vector of atoms in a molecule from its atom container + * Get the set of N coordinates 3d of a molecule of N atoms from its atom container * *@param molecule molecule store in an AtomContainer - *@return GVector with 3N coordinates (N: atom numbers) + *@return Vector with the N coordinates 3d of a molecule of N atoms */ public Vector getPoint3dCoordinates(AtomContainer molecule) { - System.out.println("molecule: " + molecule.toString()); - System.out.println("Atoms number = " + molecule.getAtomCount()); + //System.out.println("molecule: " + molecule.toString()); + //System.out.println("Atoms number = " + molecule.getAtomCount()); Vector point3dCoordinates = new Vector(); Atom thisAtom = new Atom(); for (int i = 0; i < molecule.getAtomCount(); i++) { thisAtom = molecule.getAtomAt(i); - System.out.println("thisAtom = " + thisAtom); - System.out.println("thisAtom.getPoint3d() = " + thisAtom.getPoint3d()); + //System.out.println("thisAtom = " + thisAtom); + //System.out.println("thisAtom.getPoint3d() = " + thisAtom.getPoint3d()); point3dCoordinates.add( new Point3d(thisAtom.getPoint3d()) ); //Point3d ia = (Point3d)point3dCoordinates.get(i); //System.out.println(i + "a = " + ia); } - System.out.println("Atoms 3d coordinates : " + point3dCoordinates); + //System.out.println("Atoms 3d coordinates : " + point3dCoordinates); return point3dCoordinates; } - + + + /** + * Calculate 3d distance between two atoms in one molecule from its N coordinates 3d + * + *@param atoms3dCoordinates Vector with the N coordinates 3d + *@param atomNum1 Atom position in the 3xN coordinates vector (from 0 to N-1) for the first atom. + *@param atomNum2 Atom position in the 3xN coordinates vector (from 0 to N-1) for the second atom. + *@return Distance between the two atoms in the molecule. + */ public double get3dDistanceBetweenTwoAtoms(Vector atoms3dCoordinates, int atomNum1, int atomNum2) { Point3d atom13dCoord = (Point3d)atoms3dCoordinates.get(atomNum1); @@ -88,21 +98,32 @@ double atomsDistance = 0; atomsDistance = atom13dCoord.distance(atom23dCoord); - System.out.println("atomsDistance = " + atomsDistance); + //System.out.println("atomsDistance = " + atomsDistance); return atomsDistance; } - public double get3dDistanceBetweenTwoAtomOfTwoDiffMolecules(GVector atomsCoordinatesVector1, GVector atomsCoordinatesVector2, int atomNumM1, int atomNumM2) { + + + /** + * Calculate 3d distance between two atoms from two different 3xN coordinate vectors. + * + *@param atomsCoordinatesVector1 3xN coordinates from the first molecule. + *@param atomsCoordinatesVector2 3xN coordinates from the second molecule. + *@param atomNumM1 Atom position in the first molecule. + *@param atomNumM2 Atom position in the second molecule. + *@return Distance between the two atoms. + */ + public double calculate3dDistanceBetweenTwoAtomOfTwoDiff3xNCoordinates(GVector atomsCoordinatesVector1, GVector atomsCoordinatesVector2, int atomNumM1, int atomNumM2) { double atomsDistance = 0; double difference = 0; for (int j = 0; j < 3; j++) { - difference = atomsCoordinatesVector2.getElement(atomNumM2*3+j) - atomsCoordinatesVector2.getElement(atomNumM1*3+j); + difference = atomsCoordinatesVector2.getElement(atomNumM2*3+j) - atomsCoordinatesVector1.getElement(atomNumM1*3+j); difference = Math.pow(difference, 2); atomsDistance = atomsDistance + difference; } atomsDistance = Math.sqrt(atomsDistance); - System.out.println("atomsDistance = " + atomsDistance); + //System.out.println("atomsDistance = " + atomsDistance); return atomsDistance; } } Index: GeometricMinimizer.java =================================================================== RCS file: /cvsroot/cdk/cdk/src/org/openscience/cdk/modeling/forcefield/GeometricMinimizer.java,v retrieving revision 1.8 retrieving revision 1.9 diff -u -r1.8 -r1.9 --- GeometricMinimizer.java 18 Jan 2005 10:57:31 -0000 1.8 +++ GeometricMinimizer.java 21 Jan 2005 12:58:58 -0000 1.9 @@ -14,80 +14,99 @@ */ public class GeometricMinimizer { + int SDMaximumIteration = 0; + int CGMaximumIteration = 0; + int NRMaximumIteration = 0; + + double CGconvergenceCriterion = 0.001; + double SDconvergenceCriterion = 0.001; + double NRconvergenceCriterion = 0.001; + GVector kCoordinates = new GVector(3); GVector kplus1Coordinates = new GVector(3); + int AtomsNumber =1; + int dimension = 1; + GVector gradient = new GVector(3); - GVector minimumCoordinates = new GVector(3); + + GVector steepestDescentsMinimum = new GVector(3); + GVector conjugateGradientMinimum = new GVector(3); + GVector newtonRaphsonMinimum = new GVector(3); + int iterationNumber = 0; + int iterationNumberBefore = -1; boolean convergence = false; + double RMSD = 0; + double RMS = 0; + double d = 0; + ForceField forceFieldObject = new ForceField(); + +/* double B11 = 0; double B12 = 0; double B21 = 0; double B22 = 0; double B = 0; - double RMSD = 0; - double RMS = 0; - double d = 0; double difference = 0; - int dimension = 1; double temp = 0; - +*/ /** * Constructor for the GeometricMinimizer object */ - public GeometricMinimizer() { } - /** - * Gets the minimumCoordinates attribute of the GeometricMinimizer object - * - *@return The minimumCoordinates value - */ - public GVector getMinimumCoordinates() { - return minimumCoordinates; + public void initializeMinimizationParameters(GVector initialCoord) { + + dimension = initialCoord.getSize(); + AtomsNumber = dimension/3; + kCoordinates.setSize(dimension); + kCoordinates.set(initialCoord); + //System.out.println("Coordinates at iteration 1: X1 = " + kCoordinates); + gradient.setSize(dimension); + kplus1Coordinates.setSize(kCoordinates.getSize()); + + convergence = false; + iterationNumber = 0; + iterationNumberBefore = -1; + + return; } /** * To check convergence */ - public void checkConvergence() { - + public void checkConvergence(double convergenceCriterion) { + RMS = 0; RMS = gradient.dot(gradient); RMS = RMS / dimension; RMS = Math.sqrt(RMS); - ////ln(""); //System.out.print("RMS = " + RMS); - if (RMS < 0.001) { + + if (RMS < convergenceCriterion) { convergence = true; //System.out.println("RMS convergence"); } RMSD = 0; - for (int i = 0; i < dimension; i++) { - d = 0; - for (int j = 0; j < 3; j++) { - difference = kplus1Coordinates.getElement(i + j) - kCoordinates.getElement(i + j); - difference = Math.pow(difference, 2); - d = d + difference; - } - d = Math.sqrt(d); - i = i + 2; + for (int i = 0; i < AtomsNumber; i++) { + d = forceFieldObject.calculate3dDistanceBetweenTwoAtomOfTwoDiff3xNCoordinates(kplus1Coordinates, kCoordinates, i, i); RMSD = RMSD + Math.pow(d, 2); } RMSD = RMSD / dimension; RMSD = Math.sqrt(RMSD); //System.out.print("RMSD = " + RMSD); - if (RMSD < 0.001) { + if (RMSD < convergenceCriterion) { convergence = true; //System.out.println("RMSD convergence"); } +/* if (iterationNumber == 1) { B11 = kplus1Coordinates.norm() / kCoordinates.norm(); } @@ -113,7 +132,7 @@ //System.out.println("Superlinear convergence"); } else { - if (B < 0.000001) { + if (B < convergenceCriterion) { convergence = true; //System.out.println("Linear convergence"); } @@ -145,44 +164,55 @@ B = B22 - B21; B = Math.abs(B); //System.out.print("B = " + B); - if (B < 0.0000001) { + if (B < convergenceCriterion) { convergence = true; //System.out.println("Quadratic convergence"); } } +*/ + return; + } + + + public void setkplus1Coordinates(GVector direction, double stepSize) { + kplus1Coordinates.set(direction); + kplus1Coordinates.scale(stepSize); + kplus1Coordinates.add(kCoordinates); return; } - public void setConvergenceParameters(){ + /** + * Set convergence parameters for Steepest Decents Method + * + * @param changeSDMaximumIteration Maximum number of iteration for steepest descents method. + * @param changeSDConvergenceCriterion Convergence criterion for steepest descents method. + */ + public void setConvergenceParametersForSDM(int changeSDMaximumIteration, double changeSDConvergenceCriterion){ + SDMaximumIteration = changeSDMaximumIteration; + SDconvergenceCriterion = changeSDConvergenceCriterion; return; } + /** - * Optimize the potential energy function + * Minimize the potential energy function using steepest descents method * - * @param SDMaximumIteration Maximum number of iteration for steepest descents method - * @param CGMaximumIteration Maximum number of iteration for conjugate gradient method - * @param NRMaximumIteration Maximum number of iteration for Newton-Raphson method * @param forceField The potential function to be used */ - public void energyMinimization(GVector initialCoordinates, int SDMaximumIteration, int CGMaximumIteration, int NRMaximumIteration, PotentialFunction forceField) { + public void steepestDescentsMinimization(GVector initialCoordinates, PotentialFunction forceField) { - dimension = initialCoordinates.getSize(); - kCoordinates.setSize(initialCoordinates.getSize()); - kCoordinates.set(initialCoordinates); - //System.out.println("Coordinates at iteration 1: X1 = " + kCoordinates); - gradient.setSize(initialCoordinates.getSize()); + initializeMinimizationParameters(initialCoordinates); gradient.set(forceField.gradientInPoint(kCoordinates)); //System.out.println("gradient at iteration 1 : g1 = " + gradient); - kplus1Coordinates.setSize(kCoordinates.getSize()); - minimumCoordinates.setSize(kCoordinates.getSize()); - minimumCoordinates.set(kCoordinates); - int iterationNumberBefore = -1; + + steepestDescentsMinimum.setSize(kCoordinates.getSize()); + steepestDescentsMinimum.set(kCoordinates); + LineSearch ls = new LineSearch(); - - System.out.println(""); - System.out.println("FORCEFIELDTESTS steepestDescentTest"); + + //System.out.println(""); + //System.out.println("FORCEFIELDTESTS steepestDescentTest"); SteepestDescentsMethod sdm = new SteepestDescentsMethod(kCoordinates); @@ -190,114 +220,221 @@ iterationNumber += 1; iterationNumberBefore += 1; - System.out.println(""); - System.out.println(""); - System.out.println("SD Iteration number: " + iterationNumber); + //System.out.println(""); + //System.out.println("SD Iteration number: " + iterationNumber); if (iterationNumber != 1) { kCoordinates.set(kplus1Coordinates); } //System.out.println("Search direction: "); sdm.setSk(gradient, iterationNumber); + //System.out.println(""); //System.out.println("Start the line search: "); ls.bracketingTheMinimum(kCoordinates, sdm.sk, forceField, iterationNumberBefore); - kplus1Coordinates.set(sdm.sk); - kplus1Coordinates.scale(ls.arbitraryStepSize); - kplus1Coordinates.add(kCoordinates); - System.out.print("x" + iterationNumber + " = " + kplus1Coordinates + " "); - System.out.println("f(x" + iterationNumber + ") = " + forceField.functionInPoint(kplus1Coordinates)); + setkplus1Coordinates(sdm.sk, ls.arbitraryStepSize); + //System.out.print("x" + iterationNumber + " = " + kplus1Coordinates + " "); + //System.out.println("f(x" + iterationNumber + ") = " + forceField.functionInPoint(kplus1Coordinates)); + //System.out.println("Parabolic interpolation: "); ls.parabolicInterpolation(); - kplus1Coordinates.set(sdm.sk); - kplus1Coordinates.scale(ls.parabolMinimumLambda); - kplus1Coordinates.add(kCoordinates); - System.out.print("x" + iterationNumber + " = " + kplus1Coordinates + " "); - System.out.println("f(x" + iterationNumber + ") = " + forceField.functionInPoint(kplus1Coordinates)); - //kplus1Coordinates.set(sdm.setNewCoordinates(kCoordinates, sdm.getArbitraryStepSize())); + setkplus1Coordinates(sdm.sk, ls.parabolMinimumLambda); + //System.out.print("x" + iterationNumber + " = " + kplus1Coordinates + " "); + //System.out.println("f(x" + iterationNumber + ") = " + forceField.functionInPoint(kplus1Coordinates)); gradient.set(forceField.gradientInPoint(kplus1Coordinates)); - checkConvergence(); - System.out.println("convergence: " + convergence); - System.out.println(""); - System.out.println("f(x" + iterationNumberBefore + ") = " + forceField.functionInPoint(kCoordinates)); - System.out.println("f(x" + iterationNumber + ") = " + forceField.functionInPoint(kplus1Coordinates)); + + checkConvergence(SDconvergenceCriterion); + //System.out.println("convergence: " + convergence); + + //System.out.println(""); + //System.out.println("f(x" + iterationNumberBefore + ") = " + forceField.functionInPoint(kCoordinates)); + //System.out.println("f(x" + iterationNumber + ") = " + forceField.functionInPoint(kplus1Coordinates)); } - minimumCoordinates.set(kplus1Coordinates); - System.out.println("The SD minimum energy is at: " + minimumCoordinates); + steepestDescentsMinimum.set(kplus1Coordinates); + // System.out.println("The SD minimum energy is at: " + steepestDescentsMinimum); - convergence = false; - iterationNumber = 0; - iterationNumberBefore = -1; - - System.out.println(""); - System.out.println("FORCEFIELDTESTS ConjugatedGradientTest"); + return; + } + + + /** + * Gets the steepestDescentsMinimum attribute of the GeometricMinimizer object + * + *@return The minimumCoordinates value + */ + public GVector getSteepestDescentsMinimum() { + return steepestDescentsMinimum; + } + + + /** + * Set convergence parameters for Conjugate Gradient Method + * + * @param changeCGMaximumIteration Maximum number of iteration for conjugated gradient method + * @param changeCGConvergenceCriterion Convergence criterion for conjugated gradient method + */ + public void setConvergenceParametersForCGM(int changeCGMaximumIteration, double changeCGConvergenceCriterion){ + CGMaximumIteration = changeCGMaximumIteration; + CGconvergenceCriterion = changeCGConvergenceCriterion; + return; + } + + + /** + * Minimize the potential energy function using conjugate gradient method + * + * @param forceField The potential function to be used + */ + public void conjugateGradientMinimization(GVector initialCoordinates, PotentialFunction forceField) { + + initializeMinimizationParameters(initialCoordinates); + gradient.set(forceField.gradientInPoint(kCoordinates)); + //System.out.println("gradient at iteration 1 : g1 = " + gradient); + + conjugateGradientMinimum.setSize(kCoordinates.getSize()); + conjugateGradientMinimum.set(kCoordinates); + LineSearch ls = new LineSearch(); - ConjugateGradientMethod cgm = new ConjugateGradientMethod(minimumCoordinates); - kCoordinates.set(minimumCoordinates); + //System.out.println(""); + //System.out.println("FORCEFIELDTESTS ConjugatedGradientTest"); + + ConjugateGradientMethod cgm = new ConjugateGradientMethod(kCoordinates); while ((iterationNumber < CGMaximumIteration) & (convergence == false)) { iterationNumber += 1; iterationNumberBefore += 1; - System.out.println(""); - System.out.println(""); - System.out.println("CG Iteration number: " + iterationNumber); + //System.out.println(""); + //System.out.println(""); + //System.out.println("CG Iteration number: " + iterationNumber); if (iterationNumber != 1) { - cgm.setµk(kCoordinates, kplus1Coordinates, forceField); + cgm.setuk(kCoordinates, kplus1Coordinates, forceField); kCoordinates.set(kplus1Coordinates); } //System.out.println("Search direction: "); cgm.setvk(gradient, iterationNumber); + //System.out.println(""); //System.out.println("Start the line search: "); ls.bracketingTheMinimum(kCoordinates, cgm.vk, forceField, iterationNumberBefore); - kplus1Coordinates.set(cgm.vk); - kplus1Coordinates.scale(ls.arbitraryStepSize); - kplus1Coordinates.add(kCoordinates); - System.out.print("x" + iterationNumber + " = " + kplus1Coordinates + " "); - System.out.println("f(x" + iterationNumber + ") = " + forceField.functionInPoint(kplus1Coordinates)); - System.out.println("Parabolic interpolation: "); + setkplus1Coordinates(cgm.vk, ls.arbitraryStepSize); + //System.out.print("x" + iterationNumber + " = " + kplus1Coordinates + " "); + //System.out.println("f(x" + iterationNumber + ") = " + forceField.functionInPoint(kplus1Coordinates)); + + //System.out.println("Parabolic interpolation: "); ls.parabolicInterpolation(); - kplus1Coordinates.set(cgm.vk); - kplus1Coordinates.scale(ls.parabolMinimumLambda); - kplus1Coordinates.add(kCoordinates); - System.out.print("x" + iterationNumber + " = " + kplus1Coordinates + " "); - System.out.println("f(x" + iterationNumber + ") = " + forceField.functionInPoint(kplus1Coordinates)); + setkplus1Coordinates(cgm.vk, ls.parabolMinimumLambda); + //System.out.print("x" + iterationNumber + " = " + kplus1Coordinates + " "); + //System.out.println("f(x" + iterationNumber + ") = " + forceField.functionInPoint(kplus1Coordinates)); gradient.set(forceField.gradientInPoint(kplus1Coordinates)); - checkConvergence(); - System.out.println("convergence: " + convergence); - System.out.println(""); - System.out.println("f(x" + iterationNumberBefore + ") = " + forceField.functionInPoint(kCoordinates)); - System.out.println("f(x" + iterationNumber + ") = " + forceField.functionInPoint(kplus1Coordinates)); + + checkConvergence(CGconvergenceCriterion); + //System.out.println("convergence: " + convergence); + + //System.out.println(""); + //System.out.println("f(x" + iterationNumberBefore + ") = " + forceField.functionInPoint(kCoordinates)); + //System.out.println("f(x" + iterationNumber + ") = " + forceField.functionInPoint(kplus1Coordinates)); } - minimumCoordinates.set(kplus1Coordinates); - System.out.println("The CG minimum energy is at: " + minimumCoordinates); + conjugateGradientMinimum.set(kplus1Coordinates); + //System.out.println("The CG minimum energy is at: " + conjugateGradientMinimum); + + return; + } + + + /** + * Gets the conjugatedGradientMinimum attribute of the GeometricMinimizer object + * + *@return The minimumCoordinates value + */ + public GVector getConjugateGradientMinimum() { + return conjugateGradientMinimum; + } + + + /** + * Set convergence parameters for Newton-Raphson Method + * + * @param changeNRMaximumIteration Maximum number of iteration for Newton-Raphson method + * @param changeNRConvergenceCriterion Convergence criterion for Newton-Raphson method + */ + public void setConvergenceParametersForNRM(int changeNRMaximumIteration, double changeNRConvergenceCriterion){ + NRMaximumIteration = changeNRMaximumIteration; + NRconvergenceCriterion = changeNRConvergenceCriterion; + return; + } + + + /** + * Minimize the potential energy function using the Newton-Raphson method + * + * @param forceField The potential function to be used + */ + public void newtonRaphsonMinimization(GVector initialCoordinates, PotentialFunction forceField) { + + initializeMinimizationParameters(initialCoordinates); + gradient.set(forceField.gradientInPoint(kCoordinates)); + //System.out.println("gradient at iteration 1 : g1 = " + gradient); + + newtonRaphsonMinimum.setSize(kCoordinates.getSize()); + newtonRaphsonMinimum.set(kCoordinates); + + //System.out.println(""); + //System.out.println("FORCEFIELDTESTS NewtonRaphsonTest"); + + NewtonRaphsonMethod nrm = new NewtonRaphsonMethod(kCoordinates); + GMatrix hessian = new GMatrix(dimension,dimension); + + while ((iterationNumber < NRMaximumIteration) & (convergence == false)) { + + iterationNumber += 1; + iterationNumberBefore += 1; + //System.out.println(""); + //System.out.println("NR Iteration number: " + iterationNumber); + + if (iterationNumber != 1) { + kCoordinates.set(kplus1Coordinates); + } + hessian.set(forceField.hessianInPoint(kCoordinates)); + //System.out.println("hessian = " + hessian); + nrm.gradientPerInverseHessian(gradient,hessian); + //System.out.println("GradientPerInverseHessian = " + nrm.getGradientPerInverseHessian()); + + setkplus1Coordinates(nrm.getGradientPerInverseHessian(), -1); + //System.out.print("x" + iterationNumber + " = " + kplus1Coordinates + " "); + //System.out.println("f(x" + iterationNumber + ") = " + forceField.functionInPoint(kplus1Coordinates)); + + + gradient.set(forceField.gradientInPoint(kplus1Coordinates)); + + checkConvergence(NRconvergenceCriterion); + //System.out.println("convergence: " + convergence); + + //System.out.println(""); + //System.out.println("f(x" + iterationNumberBefore + ") = " + forceField.functionInPoint(kCoordinates)); + //System.out.println("f(x" + iterationNumber + ") = " + forceField.functionInPoint(kplus1Coordinates)); + } + newtonRaphsonMinimum.set(kplus1Coordinates); + //System.out.println("The NR minimum energy is at: " + newtonRaphsonMinimum); - /* - * convergence = false; - * iterationNumber = 0; - * System.out.println(""); - * System.out.println("FORCEFIELDTESTS NewtonRaphsonTest"); - * NewtonRaphsonMethod nrm = new NewtonRaphsonMethod(minimumCoordinates); - * while ((iterationNumber < NRMaximumIteration) & (convergence == false)) { - * iterationNumber += 1; - * System.out.println(""); - * System.out.println("NR Iteration number: " + iterationNumber); - * if iterationNumber =! 1 { - * sdm.setArbitraryStepSize(kCoordinates, kplus1Coordinates,forceField); - * } - * kCoordinates.set(kplus1Coordinates); - * checkConvergence(); - * } - * minimumCoordinates.set(kplus1Coordinates); - * System.out.println("The NR minimum energy is at: " + minimumCoordinates); - */ return; } + + + /** + * Gets the newtonRaphsonMinimum attribute of the GeometricMinimizer object + * + *@return The newtonRaphsonMinimum value + */ + public GVector getNewtonRaphsonMinimum() { + return newtonRaphsonMinimum; + } + + } Index: LineSearch.java =================================================================== RCS file: /cvsroot/cdk/cdk/src/org/openscience/cdk/modeling/forcefield/LineSearch.java,v retrieving revision 1.3 retrieving revision 1.4 diff -u -r1.3 -r1.4 --- LineSearch.java 18 Jan 2005 11:02:57 -0000 1.3 +++ LineSearch.java 21 Jan 2005 12:58:58 -0000 1.4 @@ -146,9 +146,9 @@ } } arbitraryStepSize = lambdab; - System.out.println(""); - System.out.println("ArbitraryStep: "); - System.out.println("arbitraryStepSize = " + arbitraryStepSize); + //System.out.println(""); + //System.out.println("ArbitraryStep: "); + //System.out.println("arbitraryStepSize = " + arbitraryStepSize); return; } @@ -157,7 +157,7 @@ parabolMinimumLambda = fxa * (Math.pow(lambdac,2) - Math.pow(lambdab,2)) + fxb * (Math.pow(lambdaa,2) - Math.pow(lambdac,2)) + fxc * (Math.pow(lambdab,2) - Math.pow(lambdaa,2)); parabolMinimumLambda = parabolMinimumLambda / (fxa * (lambdac-lambdab) + fxb * (lambdaa-lambdac) + fxc * (lambdab-lambdaa)); parabolMinimumLambda = 0.5 * parabolMinimumLambda; - System.out.println("parabolMinimumLambda = " + parabolMinimumLambda); + //System.out.println("parabolMinimumLambda = " + parabolMinimumLambda); return parabolMinimumLambda; } Index: PotentialFunction.java =================================================================== RCS file: /cvsroot/cdk/cdk/src/org/openscience/cdk/modeling/forcefield/PotentialFunction.java,v retrieving revision 1.3 retrieving revision 1.4 diff -u -r1.3 -r1.4 --- PotentialFunction.java 7 Jan 2005 10:08:11 -0000 1.3 +++ PotentialFunction.java 21 Jan 2005 12:58:58 -0000 1.4 @@ -4,6 +4,7 @@ import java.lang.*; import java.util.*; import javax.vecmath.*; +import Jama.*; import org.openscience.cdk.*; /** @@ -14,8 +15,8 @@ */ public interface PotentialFunction { - double functionInWishedPoint = 0; // Function value in a point - GVector gradientInWishedPoint = new GVector(3); // Gradient value in a point + double functionInWishedPoint = 0; //Function value in a point + GVector gradientInWishedPoint = new GVector(3); //Gradient value in a point // GVector getGradientInWishedPoint(); @@ -39,6 +40,11 @@ GVector gradientInPoint(GVector point); -} - + /** + * Evaluate the hessian of the potential energy function in a given point + * + * @return Hessian value in the wished point + */ + GMatrix hessianInPoint(GVector point); +} |