Java Code Examples for cern.colt.matrix.DoubleMatrix1D#set()

void uniformVoxelSample(int s, DoubleMatrix1D[] sampRel, DoubleMatrix1D[] sampAbs, double[] trueVal,
        ObjectiveFunction of, double[] relMin, double[] relMax){
    //Draw sample # s for generateSamples, filling in trueVal, sampRel, and sampAbs
    //Generate vector relative to minimum
    DoubleMatrix1D dx = DoubleFactory1D.dense.make(numDOFs);
    //and absolute
    DoubleMatrix1D x = DoubleFactory1D.dense.make(numDOFs);

    for(int dof=0; dof<numDOFs; dof++){
        double top = relMax[dof];
        double bottom = relMin[dof];

        dx.set(dof, bottom + Math.random()*(top-bottom));
        x.set(dof, center.get(dof)+dx.get(dof));
    }

    trueVal[s] = of.getValue(x) - minE;
    
    sampRel[s] = dx;
    sampAbs[s] = x;
}
 

public ArrayList<Double> allTermValues(){
    //values of all epic terms at current curDOFVals
    ArrayList<Double> ans = new ArrayList<>();
    
    if(curDOFVals==null){
        throw new RuntimeException("ERROR: Trying to evaluate an EPICEnergyFunction "
                + "before assigning it to a vector of DOF values");
    }
    
    for(int termNum=0; termNum<terms.size(); termNum++){
        EPoly term = terms.get(termNum);
        
        DoubleMatrix1D DOFValsForTerm = DoubleFactory1D.dense.make(term.numDOFs);
        for(int DOFCount=0; DOFCount<term.numDOFs; DOFCount++)
            DOFValsForTerm.set( DOFCount, curDOFVals.get(termDOFs.get(termNum).get(DOFCount)) );
        
        double termVal = term.evaluate(DOFValsForTerm, includeMinE, useSharedMolec);
        ans.add(termVal);
    }
    
    return ans;
}
 

@Override
public double getEnergy() {
    if(curDOFVals==null){
        throw new RuntimeException("ERROR: Trying to evaluate an EPICEnergyFunction "
                + "before assigning it to a vector of DOF values");
    }
    
    double E = 0;
    for(int termNum=0; termNum<terms.size(); termNum++){
        EPoly term = terms.get(termNum);
        
        DoubleMatrix1D DOFValsForTerm = DoubleFactory1D.dense.make(term.numDOFs);
        for(int DOFCount=0; DOFCount<term.numDOFs; DOFCount++)
            DOFValsForTerm.set( DOFCount, curDOFVals.get(termDOFs.get(termNum).get(DOFCount)) );
        
        double termVal = term.evaluate(DOFValsForTerm, includeMinE, useSharedMolec);
        E += termVal;
    }
    
    return E;
}
 

DoubleMatrix1D getScaleAdaptive(DoubleMatrix1D sp){
    //Get an appropriate scale for sampling around some starting point sp
    /*
    //initially, we'll go about a tenth of the way across the voxel
    //in each dimension
    samplingScale = DOFmax.copy();
    samplingScale.assign(DOFmin,Functions.minus);
    samplingScale.assign( Functions.mult(0.1) );*/
    
    //Let's aim for a scaling based on what dimensions allow more flexibility
    //we'll let the scale in each dimension be ~ the distance we can go from our starting point (current x)
    //in that dimension (sum of farthest we can go in positive and negative directions)
    DoubleMatrix1D ans = DoubleFactory1D.dense.make(numDOFs);
    
    for(int dim=0; dim<numDOFs; dim++){
        
        DoubleMatrix1D y = sp.copy();
        //we'll just try and estimate distances within a factor of two
        double upDist = 1e-6;//how far up we can go in this dimension
        double downDist = 1e-6;//how far down
        
        do {
            upDist*=2;
            y.set(dim, sp.get(dim)+upDist);
        } while(checkValidPt(y));
        
        do {
            downDist*=2;
            y.set(dim, sp.get(dim)-downDist);
        } while(checkValidPt(y));
        
        
        ans.set(dim,(upDist+downDist)/6);//divide by 6 so we can make moves within the allowed region
        //though tuneScale should handle uniform scaling of samplingScale better
        //if not adapativeScale..
    }
    
    return ans;
}
 

private List<Minimizer.Result> sampleDensely(ParametricMolecule pmol, MoleculeObjectiveFunction f, int numSamplesPerDof) {

		int numDims = pmol.dofBounds.size();
		int[] dims = new int[numDims];
		Arrays.fill(dims, numSamplesPerDof);

		int numSamples = numSamplesPerDof;
		for (int i=1; i<numDims; i++) {
			numSamples *= numSamplesPerDof;
		}
		numSamples += 1;
		List<Minimizer.Result> samples = new ArrayList<>(numSamples);

		// start with the minimized center point
		samples.add(new SimpleCCDMinimizer(f).minimizeFromCenter());

		// sample points from a dense regular grid
		for (int[] p : new MathTools.GridIterable(dims)) {

			DoubleMatrix1D x = DoubleFactory1D.dense.make(numDims);
			for (int d=0; d<numDims; d++) {
				double min = pmol.dofBounds.getMin(d);
				double max = pmol.dofBounds.getMax(d);
				double xd = min + (max - min)*p[d]/(dims[d] - 1);
				x.set(d, xd);
			}

			samples.add(new Minimizer.Result(x, f.getValue(x)));
		}

		return samples;
	}
 

DoubleMatrix1D nextCandidate(){
    DoubleMatrix1D y = x.copy();
    for(int dof=0; dof<numDOFs; dof++){
        y.set(dof, x.get(dof)+random.nextGaussian()*samplingScale.get(dof) );
    }
    
    return y;
}
 

@Override
public Minimizer.Result minimizeFromCenter() {

	DoubleMatrix1D x = DoubleFactory1D.dense.make(dihedrals.size());
	for (Dihedral dihedral : dihedrals) {
		x.set(dihedral.d, (dihedral.xdmax + dihedral.xdmin)/2);
	}

	return minimizeFrom(x);
}
 

default DoubleMatrix1D getDOFsCenter() {
	int n = getNumDOFs();
	DoubleMatrix1D x = DoubleFactory1D.dense.make(n);
	for (int d=0; d<n; d++) {
		double xdmin = getConstraints()[0].get(d);
		double xdmax = getConstraints()[1].get(d);
		x.set(d, (xdmin + xdmax)/2);
	}
	return x;
}
 

DoubleMatrix1D getAnglesInRange(DoubleMatrix1D u){
    //Given a vector of this minimizer's DOFs, identify those that are angles
    //add multiples of 360 if possible to get them in the [DOFmin,DOFmax] range
    //angles are assumed to be bounded
    //either on both max and min side or not at all
    //this is for preparing initial values
    
    //we will move each out-of-range angle to the first equivalent angle after DOFmin
    DoubleMatrix1D ans = u.copy();
    
    for(int dof=0; dof<numDOFs; dof++){
        
        if( objFcn.isDOFAngle(dof) ){
        
            if(ans.get(dof)<DOFmin.get(dof)-0.001){
                ans.set(dof,
                        DOFmin.get(dof) + (ans.get(dof)-DOFmin.get(dof))%(360) + 360 );
            }
            else if(ans.get(dof)>DOFmax.get(dof)+0.001){
                ans.set(dof,
                        DOFmin.get(dof) + (ans.get(dof)-DOFmin.get(dof))%(360) );
            }
        }
    }

    return ans;
}
 

protected void scaleEigenvalues(DoubleMatrix1D eigenvalues) {
    double sum = 0.0;
    for (int i = 0; i < eigenvalues.cardinality(); ++i) {
        sum += eigenvalues.get(i);
    }

    double mean = -sum / eigenvalues.cardinality();

    for (int i = 0; i < eigenvalues.cardinality(); ++i) {
        eigenvalues.set(i, eigenvalues.get(i) / mean);
    }
}
 

private static ParametricMolecule addDofs(ParametricMolecule pmol, List<BoundedDof> newDofs) {

		// expand the dofs
		List<DegreeOfFreedom> combinedDofs = new ArrayList<>(pmol.dofs);
		for (BoundedDof bdof : newDofs) {
			combinedDofs.add(bdof.dof);
		}

		// expand the dof bounds
		DoubleMatrix1D mins = DoubleFactory1D.dense.make(pmol.dofs.size() + newDofs.size());
		DoubleMatrix1D maxs = mins.copy();
		for (int i=0; i<pmol.dofs.size(); i++) {
			mins.set(i, pmol.dofBounds.getMin(i));
			maxs.set(i, pmol.dofBounds.getMax(i));
		}
		for (int i=0; i<newDofs.size(); i++) {
			mins.set(i + pmol.dofs.size(), newDofs.get(i).min);
			maxs.set(i + pmol.dofs.size(), newDofs.get(i).max);
		}

		// build a new pmol
		return new ParametricMolecule(
			pmol.mol,
			combinedDofs,
			new ObjectiveFunction.DofBounds(new DoubleMatrix1D[] { mins, maxs })
		);
	}
 

public void getCenter(DoubleMatrix1D out) {
	for (int d=0; d<size(); d++) {
		out.set(d, getCenter(d));
	}
}
 

public void clamp(DoubleMatrix1D x) {
	for (int d=0; d<size(); d++) {
		x.set(d, clamp(d, x.get(d)));
	}
}
 

public void clampDelta(DoubleMatrix1D x, DoubleMatrix1D delta) {
	for (int d=0; d<size(); d++) {
		delta.set(d, clampDelta(d, x.get(d), delta.get(d)));
	}
}
 

@Override
public Minimizer.Result minimizeFrom(DoubleMatrix1D startx) {

	int n = f.getNumDOFs();
	DoubleMatrix1D herex = startx.copy();
	DoubleMatrix1D nextx = startx.copy();

	// ccd is pretty simple actually
	// just do a line search along each dimension until we stop improving
	// we deal with cycles by just capping the number of iterations
	
	// get the current objective function value
	double herefx = f.getValue(herex);
	
	for (int iter=0; iter<MaxIterations; iter++) {
		
		// update all the dofs using line search
		for (int d=0; d<n; d++) {
			
			LineSearcher lineSearcher = lineSearchers.get(d);
			if (lineSearcher != null) {
				
				// get the next x value for this dof
				double xd = nextx.get(d);
				xd = lineSearcher.search(xd);
				nextx.set(d, xd);
			}
		}
		
		// how much did we improve?
		double nextfx = f.getValue(nextx);
		double improvement = herefx - nextfx;
		
		if (improvement > 0) {
			
			// take the step
			herex.assign(nextx);
			herefx = nextfx;
			
			if (improvement < ConvergenceThreshold) {
				break;
			}
			
		} else {
			break;
		}
	}

	// update the protein conf, one last time
	f.setDOFs(herex);

	return new Minimizer.Result(herex, herefx);
}
 

double[] sampleResid(double[][] freeDOFVoxel){
    //draw a sample from the voxel and measure its constr resid and free-DOF resid
    
    //draw free DOFs for sample...
    DoubleMatrix1D sampFreeDOFs = DoubleFactory1D.dense.make(numFreeDOFs);
    
    for(int freeDOF=0; freeDOF<numFreeDOFs; freeDOF++){
        double voxWidth = freeDOFVoxel[1][freeDOF] - freeDOFVoxel[0][freeDOF];
        sampFreeDOFs.set( freeDOF, freeDOFVoxel[0][freeDOF] + Math.random()*voxWidth );
    }
    
    //OK now do full DOFs, placing them in NCoord and CACoord
    int fullDOFCount = 0;
    DoubleMatrix1D fullDOFVals = DoubleFactory1D.dense.make(numFullDOFs);
    
    for(int resNum=1; resNum<numRes; resNum++){
        
        for(int dim=0; dim<3; dim++){
            NCoord[resNum][dim] = evalFullDOF(fullDOFCount,sampFreeDOFs);
            fullDOFVals.set(fullDOFCount, NCoord[resNum][dim]);
            fullDOFCount++;
        }
        
        if(resNum==numRes-1)//no CA variables
            break;
        
        for(int dim=0; dim<3; dim++){
            CACoord[resNum][dim] = evalFullDOF(fullDOFCount,sampFreeDOFs);
            fullDOFVals.set(fullDOFCount, CACoord[resNum][dim]);
            fullDOFCount++;
        }
        
        
    }
    
    //Once N and CA in place, can calc C'.  Use plane projection, to match constr in jac
    for(int resNum=0; resNum<numRes-1; resNum++){
        CCoord[resNum]= pepPlanes[resNum].calcCCoords(CACoord[resNum], NCoord[resNum+1],
                CACoord[resNum+1], true);
    }
    
    //OK now handle add up constraint resids!
    ArrayList<Double> sampConstraintVals = calcConstraintVals();
    
    double constrResid = 0;
    for(int c=0; c<numFullDOFs-numFreeDOFs; c++){
        double dev = sampConstraintVals.get(c) - targetConstraintVals.get(c);
        constrResid += dev*dev;
    }
    
    constrResid /= (numFullDOFs-numFreeDOFs);//normalize resid by # of constraints
    
    
    DoubleMatrix1D freeDOFsCheck = Algebra.DEFAULT.mult(freeDOFMatrix, fullDOFVals);
    freeDOFsCheck.assign(freeDOFCenter, Functions.minus);
    freeDOFsCheck.assign(sampFreeDOFs, Functions.minus);//calc deviation in free DOFs
    double freeDOFResid = freeDOFsCheck.zDotProduct(freeDOFsCheck) / numFreeDOFs;
    
    return new double[] {constrResid, freeDOFResid};
}
 

DoubleMatrix1D getCurDOFVals(){
    DoubleMatrix1D ans = DoubleFactory1D.dense.make(numDOFs);
    for(int dof=0; dof<numDOFs; dof++)
        ans.set(dof, mms.getCurValueOfDOF(dof));
    return ans;
}
 

private void printFitTests(EPICFitter fitter, RCTuple RCList, double minEnergy,
        MoleculeModifierAndScorer mof, DoubleMatrix1D bestDOFVals, ArrayList<EPoly> series){
    //Do some tests on fit performance, and print the results
    int numDOFs = fitter.numDOFs;

    //TESTING FITS
    System.out.println("RCs: "+RCList.stringListing());
    System.out.println("Minimum energy: "+minEnergy);

    double testScales[] = new double[] { 0.01, 0.5, 5, 100 };//100
    int samplesPerScale = 3;



    double relMax[] = new double[numDOFs];//maximum shifts of degrees of freedom relative to minimum point (startVec)
    double relMin[] = new double[numDOFs];
    DoubleMatrix1D constr[] = mof.getConstraints();
    for(int dof=0; dof<numDOFs; dof++){
        relMax[dof] = constr[1].get(dof) - bestDOFVals.get(dof);
        relMin[dof] = constr[0].get(dof) - bestDOFVals.get(dof);
    }


    for(double scale : testScales){
        for(int s=0; s<samplesPerScale; s++){

            //Generate vector relative to minimum
            double dx[] = new double[numDOFs];
            //and absolute
            DoubleMatrix1D sampAbs = DoubleFactory1D.dense.make(numDOFs);
            for(int dof=0; dof<numDOFs; dof++){
                double top = Math.min(relMax[dof], scale);
                double bottom = Math.max(relMin[dof], -scale);
                dx[dof] = bottom + Math.random()*(top-bottom);
                sampAbs.set(dof, bestDOFVals.get(dof)+dx[dof]);
            }

            double trueVal = mof.getValue(sampAbs) - minEnergy;

            System.out.print("TEST: scale="+scale+" dx=");
            for(int dof=0; dof<numDOFs; dof++)
                System.out.print(dx[dof]+" ");

            System.out.print("TRUE="+trueVal+" FIT=");

            for(EPoly b : series){
                if(b!=null)
                    System.out.print(b.evaluate(sampAbs,false,false)+" ");
            }

            System.out.println();
        }
    }
}
 

protected void negateEigenvalues(DoubleMatrix1D eigenvalues) {
    for (int i = 0; i < eigenvalues.cardinality(); i++) {
        eigenvalues.set(i, -eigenvalues.get(i));
    }
}
 

protected void inverseNegateEigenvalues(DoubleMatrix1D eigenvalues) {
    for (int i = 0; i < eigenvalues.cardinality(); i++) {
        eigenvalues.set(i, -1.0 / eigenvalues.get(i));
    }
}