com.actelion.research.chem.forcefield.AbstractForceField Maven / Gradle / Ivy

Go to download
package com.actelion.research.chem.forcefield;

import java.util.ArrayList;

import com.actelion.research.chem.StereoMolecule;

public abstract class AbstractForceField implements ForceField {
	public static final double FUNCTOL = 1e-4;
    public static final double MOVETOL = 1e-7;
    public static final double EPS = 3e-8;
    public static final double TOLX = 4.0*EPS;
    public static final double MAXSTEP = 100.0;
    
	protected ArrayList listeners = new ArrayList();
    protected StereoMolecule mMol;
	protected final int mDim; ;
	protected double[] mPos;
	protected double[] mNewpos;
	protected double[] mGrad;
	protected int[] mFixedAtoms;
	protected double mTotalEnergy;
	protected long mTimeInterval; //time interval for
	protected volatile boolean mIsInterrupted;

	public AbstractForceField(StereoMolecule mol) {
		int implicitHydrogens = 0;
		for(int at=0;at0) {
			throw new IllegalArgumentException("molecule needs explicit hydrogen atoms for force field calculations");
		}
		mMol = mol;
		mDim = 3*mol.getAllAtoms();
        mGrad = new double[mDim];
        mPos = new double[mDim];
        mNewpos = new double[mDim];
        mIsInterrupted = false;
        mTimeInterval = 20;

        // get the atom positions to be placed in the pos array.
        for (int i=0; i 1 ? s / (double)(k - 1) : 0.0);
    }

    @Override
    public int minimise() {
        return minimise(4000, 1e-4, 1e-6);
    }
    
   @Override 
   public void setFixedAtoms(int[] fixedAtoms) {
	   mFixedAtoms = fixedAtoms;
   }
   
 
  public void zeroGradient() {
	   if (mFixedAtoms!=null) {
	   for (int i:mFixedAtoms) {
		   mGrad[3*i] = 0.0;
		   mGrad[3*i+1] = 0.0;
		   mGrad[3*i+2] = 0.0;
	   }
   }
   }

	
    /**
     * Minimise the current molecule.
     *  @param maxIts The maximum number of iterations to run for.
     *  @param gradTol The gradient tolerance.
     *  @param funcTol The energy tolerance.
     *  @return Return code, 0 on success.
     */
    public int minimise(int maxIts, double gradTol, double funcTol) {

        int res = 1;


        for (int i=0; i test)
	                 test = temp;
	             dGrad[i] = mGrad[i];
	         }
	
	         if (test < TOLX) {
	             return 0;
	         }
	
	         // update the gradient:
	         double gradScale = updateGradient();
	         zeroGradient();
	         // is the gradient converged?
	         test = 0.0;
	         double term = Math.max(mTotalEnergy*gradScale, 1.0);
	         for (int i=0; i Math.sqrt(EPS*sumDGrad*sumXi)) {
	             fac = 1.0/fac;
	             double fad = 1.0/fae;
	             for (int i=0; i=mTimeInterval) {
	        	 for(ForceFieldChangeListener listener: listeners) {
	        		 listener.stateChanged();
	        	 }
	        	 t0=t1;
	         }
	     }
	     return 1;
	 }

 /**
  *
  */
	 private int linearSearch(double[] oldPt,
	     double oldVal,
	     double[] dir,
	     double[] newPt,
	     double maxStep) {
	     final int MAX_ITER_LINEAR_SEARCH = 1000;
	     int ret = -1;
	     double [] tmpPt = new double[mDim];
	     double sum = 0.0, slope = 0.0, test = 0.0, lambda = 0.0;
	     double lambda2 = 0.0, lambdaMin = 0.0, tmpLambda = 0.0, val2 = 0.0;
	
	     // get the length of the direction vector:
	     sum = 0.0;
	     for (int i=0; i maxStep)
	         for (int i=0; i= 0.0)
	         return ret;
	
	     test = 0.0;
	     for (int i=0; i test)
	             test=temp;
	     }
	
	     lambdaMin = MOVETOL/test;
	     lambda = 1.0;
	     int it = 0;
	     while (it < MAX_ITER_LINEAR_SEARCH) {
	         if (lambda < lambdaMin) {
	             // the position change is too small.
	             ret = 1;
	             break;
	         }
	
	         for(int i=0; i 0.5*lambda)
	                 tmpLambda = 0.5*lambda;
	         }
	
	         lambda2 = lambda;
	         val2 = mTotalEnergy;
	         lambda = Math.max(tmpLambda, 0.1*lambda);
	         ++it;
	     }
	     // nothing was done
	     for(int i=0; i