package etomica.virial;

import etomica.api.IBox;
import etomica.api.IMoleculeList;
import etomica.api.IPotentialMaster;
import etomica.api.IRandom;
import etomica.api.ISimulation;
import etomica.integrator.mcmove.MCMoveMolecule;
import etomica.space.ISpace;
import etomica.space.IVectorRandom;


/**
 * @author kofke
 *
 * Extension of MCMoveAtom that does trial in which several atom positions are
 * perturbed.  However, position of first atom is never altered.  
 */
public class MCMoveClusterMoleculeMulti extends MCMoveMolecule {

    private static final long serialVersionUID = 1L;
    private final MeterClusterWeight weightMeter;

    public MCMoveClusterMoleculeMulti(ISimulation sim, IPotentialMaster potentialMaster,
    		                          ISpace _space) {
    	this(potentialMaster,sim.getRandom(), _space, 1.0);
    }
    
    /**
     * Constructor for MCMoveAtomMulti.
     * @param parentIntegrator
     * @param nAtoms number of atoms to move in a trial.  Number of atoms in
     * box should be at least one greater than this value (greater
     * because first atom is never moved)
     */
    public MCMoveClusterMoleculeMulti(IPotentialMaster potentialMaster,
            IRandom random, ISpace _space, double stepSize) {
        super(potentialMaster,random, _space, stepSize,Double.POSITIVE_INFINITY);
        weightMeter = new MeterClusterWeight(potential);
    }

    public void setBox(IBox p) {
        super.setBox(p);
        weightMeter.setBox(p);
        translationVectors = new IVectorRandom[box.getMoleculeList().getMoleculeCount()-1];
        for (int i=0; i<box.getMoleculeList().getMoleculeCount()-1; i++) {
            translationVectors[i] = (IVectorRandom)space.makeVector();
        }
    }
    
    //note that total energy is calculated
    public boolean doTrial() {
        uOld = weightMeter.getDataAsScalar();
        IMoleculeList moleculeList = box.getMoleculeList();
        for(int i=1; i<moleculeList.getMoleculeCount(); i++) {
            translationVectors[i-1].setRandomCube(random);
            translationVectors[i-1].TE(stepSize);
            groupTranslationVector.E(translationVectors[i-1]);
            moveMoleculeAction.actionPerformed(moleculeList.getMolecule(i));
        }
        ((BoxCluster)box).trialNotify();
        return true;
    }
	
    public void rejectNotify() {
        IMoleculeList moleculeList = box.getMoleculeList();
        for(int i=1; i<moleculeList.getMoleculeCount(); i++) {
            groupTranslationVector.E(translationVectors[i-1]);
            groupTranslationVector.TE(-1);
            moveMoleculeAction.actionPerformed(moleculeList.getMolecule(i));
        }
        ((BoxCluster)box).rejectNotify();
        if (weightMeter.getDataAsScalar() == 0) {
            throw new RuntimeException("oops oops, reverted to illegal configuration");
        }
    }

    public void acceptNotify() {
        if (uNew == 0) {
            throw new RuntimeException("oops, accepted illegal configuration");
        }
        ((BoxCluster)box).acceptNotify();
        if (weightMeter.getDataAsScalar() == 0) {
            throw new RuntimeException("oops oops, accepted illegal configuration");
        }
    }
    
    public double getB() {
        return 0.0;
    }
    
    public double getA() {
        uNew = weightMeter.getDataAsScalar();
        return uNew/uOld;
    }
	
    private IVectorRandom[] translationVectors;
}

@@ -632,16 +632,99 @@
 htMeter;
+%0A    private IVectorRandom%5B%5D translationVectors;%0A    protected int%5B%5D constraintMap;
 %0A%0A    pu
@@ -1652,18 +1652,16 @@
 eCount()
--1
 %5D;%0A     
@@ -1719,18 +1719,16 @@
 eCount()
--1
 ; i++) %7B
@@ -1801,32 +1801,379 @@
 or();%0A        %7D%0A
+        if (constraintMap == null) %7B%0A            constraintMap = new int%5Bbox.getMoleculeList().getMoleculeCount()%5D;%0A            for (int i=0; i%3CconstraintMap.length; i++) %7B%0A                constraintMap%5Bi%5D = i;%0A            %7D%0A        %7D%0A    %7D%0A    %0A    public void setConstraintMap(int%5B%5D newConstraintMap) %7B%0A        constraintMap = newConstraintMap;%0A
     %7D%0A    %0A    /
@@ -2279,32 +2279,152 @@
 DataAsScalar();%0A
+//        if (uOld == 0) %7B%0A//            throw new RuntimeException(%22oops, initial configuration unhappy%22);%0A//        %7D%0A
         IMolecul
@@ -2522,32 +2522,102 @@
 Count(); i++) %7B%0A
+            int tv = constraintMap%5Bi%5D;%0A            if (tv == i) %7B%0A    
             tran
@@ -2631,19 +2631,18 @@
 Vectors%5B
-i-1
+tv
 %5D.setRan
@@ -2662,32 +2662,36 @@
 m);%0A            
+   
+ 
 translationVecto
@@ -2693,19 +2693,18 @@
 Vectors%5B
-i-1
+tv
 %5D.TE(ste
@@ -2711,16 +2711,30 @@
 pSize);%0A
+            %7D%0A
         
@@ -2773,35 +2773,34 @@
 nslationVectors%5B
-i-1
+tv
 %5D);%0A            
@@ -3128,33 +3128,41 @@
 nslationVector.E
-(
+a1Tv1(-1,
 translationVecto
@@ -3168,55 +3168,25 @@
 ors%5B
-i-1%5D);%0A            groupTranslationVector.TE(-1
+constraintMap%5Bi%5D%5D
 );%0A 
@@ -3493,32 +3493,34 @@
 cceptNotify() %7B%0A
+//
         if (uNew
@@ -3520,32 +3520,34 @@
 f (uNew == 0) %7B%0A
+//
             thro
@@ -3602,32 +3602,34 @@
 onfiguration%22);%0A
+//
         %7D%0A      
@@ -3668,153 +3668,8 @@
 ();%0A
-        if (weightMeter.getDataAsScalar() == 0) %7B%0A            throw new RuntimeException(%22oops oops, accepted illegal configuration%22);%0A        %7D%0A
     
@@ -3844,54 +3844,6 @@
 %7D%0A%09%0A
-    private IVectorRandom%5B%5D translationVectors;%0A
 %7D%0A

allow molecules to be constrained to share a position with another molecule (they get translated together)