linesearchforthewolfeconditions.java

化学图形处理软件

				}
				if (linearFunctionDerivativeInAlphaj * (alphaHigh-alphaLow) >= 0) {		
					alphaHigh = alphaLow;
					linearFunctionInAlphaHigh = linearFunctionInAlphaLow;
					//linearFunctionDerivativeInAlphaHigh = linearFunctionDerivativeInAlphaLow;
				}
				alphaLow = alphaj;
				linearFunctionInAlphaLow = linearFunctionInAlphaj;
				linearFunctionDerivativeInAlphaLow = linearFunctionDerivativeInAlphaj;
				dfInAlphaLow = dfInAlphaj;
			}
			
			//logger.debug("AlphaLow = " + alphaLow + ", AlphaHigh = " + alphaHigh);
			//logger.debug("linearFunctionInAlphaLow = " + linearFunctionInAlphaLow + ", linearFunctionInAlphaHigh = " + linearFunctionInAlphaHigh);
			functionEvaluationNumber = functionEvaluationNumber + 1;
			//logger.debug("functionEvaluationNumber = " + functionEvaluationNumber);

			if ((functionEvaluationNumber == 10) | (Math.abs(linearFunctionInAlphaHigh - linearFunctionInAlphaLow) <= 0.000001) | (Math.abs(alphaLow - alphaHigh) <= 0.000000000001)) {
				//logger.debug("ZOOM WAS TERMINATE EARLIER");
				/*System.out.println("functionEvaluationNumber = " + functionEvaluationNumber + 
						", Math.abs(linearFunctionInAlphaHigh - linearFunctionInAlphaLow) = " + Math.abs(linearFunctionInAlphaHigh - linearFunctionInAlphaLow) + 
						", Math.abs(alphaLow - alphaHigh) = " + Math.abs(alphaLow - alphaHigh));*/
				this.alphaOptimum = alphaLow;
				this.linearFunctionInAlphaOptimum = linearFunctionInAlphaLow;
				this.dfOptimum = dfInAlphaLow;
				
				//logger.debug("(functionEvaluationNumber == 10) | (Math.abs(linearFunctionInAlphaHigh - linearFunctionInAlphaLow) <= 0.000001) | (Math.abs(alphaLow - alphaHigh) <= 0.0000001)");
				//logger.debug("zoom end -> this.alphaOptimum = " + this.alphaOptimum); 
				//logger.debug("zoom end -> this.linearFunctionInAlphaOptimum = " + this.linearFunctionInAlphaOptimum); 
				break;
			}

		
		} while ((Math.abs(linearFunctionInAlphaHigh - linearFunctionInAlphaLow) > 0.000001) 
					& (functionEvaluationNumber < 10) 
					& (Math.abs(alphaLow - alphaHigh) > 0.000000000001));
		
		//logger.debug("zoom end");
		return;
	}

	
	 /*
	 * Cubic interpolation in the interval [a,b] known to contain desirable step length 
	 * and given two previous step length estimates in this interval.
	 *
	 *@param 	alphai										Previous step length.									
	 *@param 	linearFunctionInAlphai						Function value at the previous step length alphai.
	 *@param 	linearFunctionDerivativeInAlphai			Derivative at the previous step length alphai.
	 *@param 	alphaiMinus1								Previous step length.
	 *@param 	linearFunctionInAlphaiMinus1	 			Function value at the previous step length alphaiMinus1.
	 *@param 	linearFunctionDerivativeInAlphaiMinus1	 	Derivative value at the previous step length alphaiMinus1.
	 *@param 	a											Inferior value of the interval [a,b].
	 *@param 	b											Superior value of the interval [a,b].	
 	 *	  
	 * @return												Cubic interpolation in the interval [a,b]
	 */
	public double cubicInterpolation(double alphai, double linearFunctionInAlphai, double linearFunctionDerivativeInAlphai, 
									double alphaiMinus1, double linearFunctionInAlphaiMinus1, double linearFunctionDerivativeInAlphaiMinus1, 
									double a, double b) {
		
		//logger.debug("The interval [" + a + ", " + b + "] contains acceptable step lengths.");
		
		if (alphai < alphaiMinus1) {
			this.alphaTemporal = alphai;
			this.linearFunctionInAlphaTemporal = linearFunctionInAlphai;
			this.linearFunctionDerivativeInAlphaTemporal = linearFunctionDerivativeInAlphai;
			alphai = alphaiMinus1;
			linearFunctionInAlphai = linearFunctionInAlphaiMinus1;
			linearFunctionDerivativeInAlphai = linearFunctionDerivativeInAlphaiMinus1;
			alphaiMinus1 = this.alphaTemporal;
			linearFunctionInAlphaiMinus1 = this.linearFunctionInAlphaTemporal;
			linearFunctionDerivativeInAlphaiMinus1 = this.linearFunctionDerivativeInAlphaTemporal;
		}
		
		this.d1 = linearFunctionDerivativeInAlphaiMinus1 + linearFunctionDerivativeInAlphai - 3 * ((linearFunctionInAlphaiMinus1 - linearFunctionInAlphai)/(alphaiMinus1 - alphai));
		//logger.debug("d1 = " + d1);
		
		//logger.debug("linearFunctionDerivativeInAlphaiMinus1 = " + linearFunctionDerivativeInAlphaiMinus1);
		//logger.debug("linearFunctionDerivativeInAlphai = " + linearFunctionDerivativeInAlphai);
		
		this.d2 = Math.sqrt(Math.abs(Math.pow(d1,2) - linearFunctionDerivativeInAlphaiMinus1 * linearFunctionDerivativeInAlphai));
		//logger.debug("d2 = " + d2);
		
		this.alphaiplus1 = alphai-(alphai-alphaiMinus1) * ((linearFunctionDerivativeInAlphai + d2 - d1) / (linearFunctionDerivativeInAlphai - linearFunctionDerivativeInAlphaiMinus1 + 2 * d2));
		
		//logger.debug("alphaiplus1 = " + alphaiplus1);
		
		if (alphaiplus1 < a) {alphaiplus1 = a;}
		if (alphaiplus1 > b) {alphaiplus1 = b;}
		
		//logger.debug("alphaiplus1 = " + alphaiplus1);
		
		if (Math.abs(alphaiplus1 - alphai) < 0.000000001) {
			/*System.out.println("We reset alphaiplus1 = (alphaiMinus1 + alphai) / 2, because alphaiplus1 = " + alphaiplus1 + " is too close to its predecessor " +
					"alphaiMinus1 = " + alphaiMinus1); */
			alphaiplus1 = (alphaiMinus1 + alphai) / 2;
		} else {if (alphaiplus1 < (alphai - 9 * (alphai-alphaiMinus1) / 10)) {
			//logger.debug("We reset alphaiplus1 = (alphaiMinus1 + alphai) / 2, because alphaiplus1 = " + alphaiplus1 + " is 	too much smaller than alphai = " + alphai); 
			alphaiplus1 = (alphaiMinus1 + alphai) / 2;;
			}
		}
	
		return alphaiplus1;
	}
	

	 /*
	 *	The aim is to find a value of alpha that satisfies the sufficient decrease condition, without being too small.
	 *	The procedures generate a value alphai such that is not too much smaller than its predecesor alphai-1. 
	 *  The interpolation in the first is quadratic but if the sufficient decrease condition is not satisfied 
	 *  then the interpolation is cubic.
 	 *
	 * @param alphaLow								Among all step lengths generated so far and satisfying the sufficient decrease condition, the one giving the smallest function value.
	 * @param linearFunctionInAlphaLow				Energy function value at alphaLow.
	 * @param linearFunctionDerivativeInAlphaLow	Derivative value at alphaLow.
	 * @param alphaHigh								AlphaHigh is chosen so that linearFunctionDerivativeInAlphaj * (alphaHigh-alphaLow) < 0
	 * @param linearFunctionInAlphaHigh				Energy function value at alphaHigh.
	 * @return										Value of alpha that satisfies the sufficient decrease condition, without being too small.
	 */
	private double interpolation(double alphaLow, double linearFunctionInAlphaLow, double linearFunctionDerivativeInAlphaLow, 
								double alphaHigh, double linearFunctionInAlphaHigh) {
		
		double minAlpha = Math.min(alphaLow, alphaHigh);
		double alphaDiff = Math.abs(alphaHigh - alphaLow);
		double alphaInterpolation;

		//logger.debug("We form a quadratic approximation to the linear function");
		double alpha1 = -1 * ((linearFunctionDerivativeInAlphaLow * Math.pow(alphaDiff,2)) / (2 * (linearFunctionInAlphaHigh - linearFunctionInAlphaLow - linearFunctionDerivativeInAlphaLow * alphaDiff)));
		
		//logger.debug("The value alpha1 = " + alpha1 + ", is the minimizer of this quadratic function");
		
		if ((alpha1 > alphaDiff) | (Math.abs(alpha1 - alphaDiff) < 0.000000001)) {
			if (alpha1 < 1E-7) {}
			else {
				/*System.out.println("We reset alpha1 = alphaDiff / 2, because alphaInterpolation = " + alpha1 + " is too close to its predecessor " +
						"alphaiMinus1 = " + alphaDiff); */
				alpha1 = alphaDiff / 2;
			}
		} else {
			if ((alpha1 < 0) & (alpha1 < (alphaDiff - 9 * alphaDiff / 10))) {
				if (alpha1 < 1E-7) {}
				else {
					//logger.debug("We reset alphai = alphaiMinus1 / 2, because alphaInterpolation = " + alpha1 + " is 	too much smaller than alphaiMinus1 = " + alphaDiff); 
					alpha1 = alphaDiff / 2;
				}
			}
		}

		//logger.debug("alpha1 = " + alpha1);

		alphaInterpolation = minAlpha + alpha1;
		this.linearFunctionAlphaInterpolation = this.evaluateEnergyFunction(alphaInterpolation);
		//logger.debug("alphaInterpolation = " + alphaInterpolation);
		//logger.debug("linearFunctionAlphaInterpolation = " + this.linearFunctionAlphaInterpolation);
		if (this.linearFunctionAlphaInterpolation <= this.linearFunctionInAlpha0 + this.c1 * (alphaInterpolation) * this.linearFunctionDerivativeInAlpha0) {
			//logger.debug("The sufficient decrease condition is satisfied at alpha1 and we termine the interpolation");
		}
		else {
			//double alphaiMinus2;
			//double alphaiMinus1 = alphaDiff;
			//double linearFunctionInAlphaiMinus2;
			//double linearFunctionInAlphaiMinus1 = linearFunctionInAlphaHigh;
			double alphai; // = alpha1;
			//double linearFunctionInAlphai = this.linearFunctionAlphaInterpolation;
				
			//do {
				//alphaiMinus2 = alphaiMinus1;
				//alphaiMinus1 = alphai;
				//linearFunctionInAlphaiMinus2 = linearFunctionInAlphaiMinus1;
				//linearFunctionInAlphaiMinus1 = linearFunctionInAlphai;
					
				//logger.debug("We construct a cubic function that interpolates the fours pieces of information");	
				a = 1/(Math.pow(alphaDiff,2) * Math.pow(alpha1, 2) * (alpha1-alphaDiff));
				b = a;
				a = a * (Math.pow(alphaDiff,2) * (this.linearFunctionAlphaInterpolation - linearFunctionInAlphaLow - linearFunctionDerivativeInAlphaLow * alpha1) 
						+ (-Math.pow(alpha1,2)) * (linearFunctionInAlphaHigh - linearFunctionInAlphaLow - linearFunctionDerivativeInAlphaLow * alphaDiff));
				b = b * (- Math.pow(alphaDiff,3) * (this.linearFunctionAlphaInterpolation - linearFunctionInAlphaLow - linearFunctionDerivativeInAlphaLow * alpha1) 
						+ Math.pow(alpha1,3) * (linearFunctionInAlphaHigh - linearFunctionInAlphaLow - linearFunctionDerivativeInAlphaLow * alphaDiff));
				
				//logger.debug("a = " + a);
				//logger.debug("b = " + b);
				
				alphai = (-b + Math.sqrt(Math.pow(b,2) - 3 * a * linearFunctionDerivativeInAlphaLow)) / (3 * a);
				//logger.debug("alphai = " + alphai);
				
				if (Math.abs(alphai - alpha1) < 0.000000001) {
					/*System.out.println("We reset alphai = alpha1 / 2, because alphaInterpolation = " + alphai + " is too close to its predecessor " +
							"alpha1 = " + alpha1); */
					alphai = alpha1 / 2;
				} else {if (alphai < (alpha1 - 9 * alpha1 / 10)) {
					//logger.debug("We reset alphai = alpha1 / 2, because alphaInterpolation = " + alphai + " is 	too much smaller than alpha1 = " + alpha1); 
					alphai = alpha1 / 2;
					}
				}
			
				alphaInterpolation = minAlpha + alphai;
				this.linearFunctionAlphaInterpolation = this.evaluateEnergyFunction(alphaInterpolation);
				//logger.debug("alphaInterpolation = " + alphaInterpolation);
				//logger.debug("linearFunctionAlphaInterpolation = " + this.linearFunctionAlphaInterpolation);
				//functionEvaluationNumber = functionEvaluationNumber + 1;
				
			/*} while (((linearFunctionInAlphai > linearFunctionInAlphaLow + this.c1 * (alphaLow + alphai) * linearFunctionDerivativeInAlphaLow) & (functionEvaluationNumber < 5)) 
					| ((linearFunctionInAlphai - this.linearFunctionAlphaInterpolation) < 0.00000001) | ((alphai - alpha1) < 0.00000001));*/
				
				
		}
			
		return alphaInterpolation;
	}
	

	/**Evaluate the energy function from an alpha value, using the current coordinates and the current direction.
	 * 
	 * @param alpha	
	 * @return			Energy function value.
	 */
	private double evaluateEnergyFunction(double alpha) {
		//logger.debug("alpha= " + alpha);
		this.xAlpha.set(this.x);
		//logger.debug("xAlpha = " + xAlpha);
		GVector directionStep = direction;
		//logger.debug("directionStep = " + directionStep);
		xAlpha.scaleAdd(alpha, directionStep, xAlpha);
		//logger.debug("xAlpha = " + xAlpha);
		double fxAlpha = pf.energyFunction(xAlpha);
		//logger.debug("fxAlpha = " + fxAlpha);
		return fxAlpha;
	}
	
	
	/**Evaluate the gradient of the energy function from an alpha value, 
	 * using the current coordinates and the current direction.
	 * 
	 * @param alpha		Alpha value for the one-dimensional problem generate from the current coordinates and the current direction.
	 * @return				Gradient of the energy function at alpha. 
	 */
	private GVector evaluateEnergyFunctionDerivative(double alpha) {
		//logger.debug("alpha= " + alpha);
		this.xAlpha.set(this.x);
		//logger.debug("xAlpha = " + xAlpha);
		GVector directionStep = direction;
		//logger.debug("directionStep = " + directionStep);
		xAlpha.scaleAdd(alpha, directionStep, xAlpha);
		//logger.debug("xAlpha = " + xAlpha);
		pf.setEnergyGradient(xAlpha);
		GVector dfxAlpha = pf.getEnergyGradient();
		//logger.debug("dfxAlpha = " + dfxAlpha);
		return dfxAlpha;
	}

	
	/**
	 *	From the interval [a, b] that bracket the minimum, evaluates the energy function at an intermediate point x 
	 *  and obtain a new, smaller bracketing interval, either (a,x) or (x,b).	 
	 *
	 * @param lambdaMin		a
	 * @param flambdaMin		Energy function at a.
	 * @param lambdaMax		b
	 * @param flambdaMax		Energy function at b.
	 * @return					An intermediate point x
	 */
	/*private double goldenSectionMethod(double lambdaMin, double flambdaMin, double lambdaMax, double flambdaMax) {

		//logger.debug("Golden Section Search");
		double goldenLambda;
		
		double lambda1 = lambdaMin;
		double lambda4 = lambdaMax;
		double lambda2 = lambda1 + 0.3819660112 * (lambda4 - lambda1);
		double lambda3 = lambda1 + 0.6180339887498948482 * (lambda4 - lambda1);
		
		//double flambda1 = flambdaMin;
		double flambda2 = evaluateEnergyFunction(lambda2);
		double flambda3 = evaluateEnergyFunction(lambda3);
		//double flambda4 = flambdaMax;
		
		//logger.debug("lambda1 = " + lambda1 + ", flambda1 = " + flambda1);
		//logger.debug("lambda2 = " + lambda2 + ", flambda2 = " + flambda2);
		//logger.debug("lambda3 = " + lambda3 + ", flambda3 = " + flambda3);
		//logger.debug("lambda4 = " + lambda4 + ", flambda4 = " + flambda4);

		if (flambda2 < flambda3) {		//we can bracket the minimum by the interval [lambda1, lambda3]
			goldenLambda = lambda2;
			linearFunctionGoldenAlpha = flambda2;
		}	
		else {			//we can bracket the minimum by the interval [lambda2, lambda4]
			goldenLambda = lambda3;
			linearFunctionGoldenAlpha = flambda3;
		}

		return goldenLambda;
	}*/


}