Index: src/Dynamics/LinearInterpolationBetweenSteps.hpp
===================================================================
--- src/Dynamics/LinearInterpolationBetweenSteps.hpp	(revision fe6f202b08255305d13ec735eddd6dc5bc8cdcf4)
+++ src/Dynamics/LinearInterpolationBetweenSteps.hpp	(revision e355762f9447acb9c234eadd3df6adc314e5dee7)
@@ -59,8 +59,9 @@
     atom *Sprinter = NULL;
     if (!MapByIdentity) {
+      std::map<atom *, atom*> PermutationMap;
       molecule::atomSet atoms_list;
       copy(atoms.begin(), atoms.end(), atoms_list.begin());
-      MinimiseConstrainedPotential Minimiser(atoms_list);
-      Minimiser(PermutationMap, startstep, endstep, IsAngstroem);
+      MinimiseConstrainedPotential Minimiser(atoms_list, PermutationMap);
+      Minimiser(startstep, endstep, IsAngstroem);
     } else {
       // TODO: construction of enumeration goes here
Index: src/Dynamics/MinimiseConstrainedPotential.cpp
===================================================================
--- src/Dynamics/MinimiseConstrainedPotential.cpp	(revision fe6f202b08255305d13ec735eddd6dc5bc8cdcf4)
+++ src/Dynamics/MinimiseConstrainedPotential.cpp	(revision e355762f9447acb9c234eadd3df6adc314e5dee7)
@@ -39,6 +39,9 @@
 
 
-MinimiseConstrainedPotential::MinimiseConstrainedPotential(molecule::atomSet &_atoms) :
-  atoms(_atoms)
+MinimiseConstrainedPotential::MinimiseConstrainedPotential(
+    molecule::atomSet &_atoms,
+    std::map<atom*, atom *> &_PermutationMap) :
+  atoms(_atoms),
+  PermutationMap(_PermutationMap)
 {}
 
@@ -46,12 +49,7 @@
 {}
 
-double MinimiseConstrainedPotential::operator()(atom **&PermutationMap, int _startstep, int _endstep, bool IsAngstroem)
+double MinimiseConstrainedPotential::operator()(int _startstep, int _endstep, bool IsAngstroem)
 {
   double Potential, OldPotential, OlderPotential;
-  PermutationMap = new atom *[atoms.size()];
-  DistanceList = new DistanceMap *[atoms.size()];
-  DistanceIterators = new DistanceMap::iterator[atoms.size()];
-  DoubleList = new int[atoms.size()];
-  StepList = new DistanceMap::iterator[atoms.size()];
   int round;
   atom *Sprinter = NULL;
@@ -59,8 +57,12 @@
 
   // set to zero
-  for (size_t i=0;i<atoms.size();i++) {
-    PermutationMap[i] = NULL;
-    DoubleList[i] = 0;
-  }
+  PermutationMap.clear();
+  DoubleList.clear();
+  for (molecule::atomSet::const_iterator iter = atoms.begin(); iter != atoms.end(); ++iter) {
+    DistanceList[*iter].clear();
+  }
+  DistanceList.clear();
+  DistanceIterators.clear();
+  DistanceIterators.clear();
 
   /// Minimise the potential
@@ -79,5 +81,5 @@
   DoLog(1) && (Log() << Verbose(1) << "Making the PermutationMap injective ... " << endl);
   MakeInjectivePermutation();
-  delete[](DoubleList);
+  DoubleList.clear();
 
   // argument minimise the constrained potential in this injective PermutationMap
@@ -92,18 +94,19 @@
       iter = atoms.begin();
       for (; iter != atoms.end(); ++iter) {
+        CalculateDoubleList();
         PrintPermutationMap();
-        Sprinter = DistanceIterators[(*iter)->getNr()]->second;   // store initial partner
-        Strider = DistanceIterators[(*iter)->getNr()];  //remember old iterator
-        DistanceIterators[(*iter)->getNr()] = StepList[(*iter)->getNr()];
-        if (DistanceIterators[(*iter)->getNr()] == DistanceList[(*iter)->getNr()]->end()) {// stop, before we run through the list and still on
-          DistanceIterators[(*iter)->getNr()] == DistanceList[(*iter)->getNr()]->begin();
+        Sprinter = DistanceIterators[(*iter)]->second;   // store initial partner
+        Strider = DistanceIterators[(*iter)];  //remember old iterator
+        DistanceIterators[(*iter)] = StepList[(*iter)];
+        if (DistanceIterators[(*iter)] == DistanceList[(*iter)].end()) {// stop, before we run through the list and still on
+          DistanceIterators[(*iter)] == DistanceList[(*iter)].begin();
           break;
         }
-        //Log() << Verbose(2) << "Current Walker: " << *(*iter) << " with old/next candidate " << *Sprinter << "/" << *DistanceIterators[(*iter)->getNr()]->second << "." << endl;
+        //Log() << Verbose(2) << "Current Walker: " << *(*iter) << " with old/next candidate " << *Sprinter << "/" << *DistanceIterators[(*iter)]->second << "." << endl;
         // find source of the new target
         molecule::atomSet::const_iterator runner = atoms.begin();
         for (; runner != atoms.end(); ++runner) { // find the source whose toes we might be stepping on (Walker's new target should be in use by another already)
-          if (PermutationMap[(*runner)->getNr()] == DistanceIterators[(*iter)->getNr()]->second) {
-            //Log() << Verbose(2) << "Found the corresponding owner " << *(*runner) << " to " << *PermutationMap[(*runner)->getNr()] << "." << endl;
+          if (PermutationMap[(*runner)] == DistanceIterators[(*iter)]->second) {
+            //Log() << Verbose(2) << "Found the corresponding owner " << *(*runner) << " to " << *PermutationMap[(*runner)] << "." << endl;
             break;
           }
@@ -111,13 +114,14 @@
         if (runner != atoms.end()) { // we found the other source
           // then look in its distance list for Sprinter
-          Rider = DistanceList[(*runner)->getNr()]->begin();
-          for (; Rider != DistanceList[(*runner)->getNr()]->end(); Rider++)
+          Rider = DistanceList[(*runner)].begin();
+          for (; Rider != DistanceList[(*runner)].end(); Rider++)
             if (Rider->second == Sprinter)
               break;
-          if (Rider != DistanceList[(*runner)->getNr()]->end()) { // if we have found one
-            //Log() << Verbose(2) << "Current Other: " << *(*runner) << " with old/next candidate " << *PermutationMap[(*runner)->getNr()] << "/" << *Rider->second << "." << endl;
+          if (Rider != DistanceList[(*runner)].end()) { // if we have found one
+            //Log() << Verbose(2) << "Current Other: " << *(*runner) << " with old/next candidate " << *PermutationMap[(*runner)] << "/" << *Rider->second << "." << endl;
             // exchange both
-            PermutationMap[(*iter)->getNr()] = DistanceIterators[(*iter)->getNr()]->second; // put next farther distance into PermutationMap
-            PermutationMap[(*runner)->getNr()] = Sprinter;  // and hand the old target to its respective owner
+            PermutationMap[(*iter)] = DistanceIterators[(*iter)]->second; // put next farther distance into PermutationMap
+            PermutationMap[(*runner)] = Sprinter;  // and hand the old target to its respective owner
+            CalculateDoubleList();
             PrintPermutationMap();
             // calculate the new potential
@@ -126,13 +130,13 @@
             if (Potential > OldPotential) { // we made everything worse! Undo ...
               //Log() << Verbose(3) << "Nay, made the potential worse: " << Potential << " vs. " << OldPotential << "!" << endl;
-              //Log() << Verbose(3) << "Setting " << *(*runner) << "'s source to " << *DistanceIterators[(*runner)->getNr()]->second << "." << endl;
+              //Log() << Verbose(3) << "Setting " << *(*runner) << "'s source to " << *DistanceIterators[(*runner)]->second << "." << endl;
               // Undo for Runner (note, we haven't moved the iteration yet, we may use this)
-              PermutationMap[(*runner)->getNr()] = DistanceIterators[(*runner)->getNr()]->second;
+              PermutationMap[(*runner)] = DistanceIterators[(*runner)]->second;
               // Undo for Walker
-              DistanceIterators[(*iter)->getNr()] = Strider;  // take next farther distance target
-              //Log() << Verbose(3) << "Setting " << *(*iter) << "'s source to " << *DistanceIterators[(*iter)->getNr()]->second << "." << endl;
-              PermutationMap[(*iter)->getNr()] = DistanceIterators[(*iter)->getNr()]->second;
+              DistanceIterators[(*iter)] = Strider;  // take next farther distance target
+              //Log() << Verbose(3) << "Setting " << *(*iter) << "'s source to " << *DistanceIterators[(*iter)]->second << "." << endl;
+              PermutationMap[(*iter)] = DistanceIterators[(*iter)]->second;
             } else {
-              DistanceIterators[(*runner)->getNr()] = Rider;  // if successful also move the pointer in the iterator list
+              DistanceIterators[(*runner)] = Rider;  // if successful also move the pointer in the iterator list
               DoLog(3) && (Log() << Verbose(3) << "Found a better permutation, new potential is " << Potential << " vs." << OldPotential << "." << endl);
               OldPotential = Potential;
@@ -148,7 +152,7 @@
           }
         } else {
-          PermutationMap[(*iter)->getNr()] = DistanceIterators[(*iter)->getNr()]->second; // new target has no source!
+          PermutationMap[(*iter)] = DistanceIterators[(*iter)]->second; // new target has no source!
         }
-        StepList[(*iter)->getNr()]++; // take next farther distance target
+        StepList[(*iter)]++; // take next farther distance target
       }
     } while (++iter != atoms.end());
@@ -157,9 +161,4 @@
 
 
-  /// free memory and return with evaluated potential
-  for (int i=atoms.size(); i--;)
-    DistanceList[i]->clear();
-  delete[](DistanceList);
-  delete[](DistanceIterators);
   return ConstrainedPotential();
 };
@@ -167,12 +166,7 @@
 void MinimiseConstrainedPotential::FillDistanceList()
 {
-  for (int i=atoms.size(); i--;) {
-    DistanceList[i] = new DistanceMap;    // is the distance sorted target list per atom
-    DistanceList[i]->clear();
-  }
-
   for (molecule::atomSet::const_iterator iter = atoms.begin(); iter != atoms.end(); ++iter) {
     for (molecule::atomSet::const_iterator runner = atoms.begin(); runner != atoms.end(); ++runner) {
-      DistanceList[(*iter)->getNr()]->insert( DistancePair((*iter)->getPositionAtStep(startstep).distance((*runner)->getPositionAtStep(endstep)), (*runner)) );
+      DistanceList[(*iter)].insert( DistancePair((*iter)->getPositionAtStep(startstep).distance((*runner)->getPositionAtStep(endstep)), (*runner)) );
     }
   }
@@ -182,9 +176,9 @@
 {
   for (molecule::atomSet::const_iterator iter = atoms.begin(); iter != atoms.end(); ++iter) {
-    StepList[(*iter)->getNr()] = DistanceList[(*iter)->getNr()]->begin();    // stores the step to the next iterator that could be a possible next target
-    PermutationMap[(*iter)->getNr()] = DistanceList[(*iter)->getNr()]->begin()->second;   // always pick target with the smallest distance
-    DoubleList[DistanceList[(*iter)->getNr()]->begin()->second->getNr()]++;            // increase this target's source count (>1? not injective)
-    DistanceIterators[(*iter)->getNr()] = DistanceList[(*iter)->getNr()]->begin();    // and remember which one we picked
-    DoLog(2) && (Log() << Verbose(2) << **iter << " starts with distance " << DistanceList[(*iter)->getNr()]->begin()->first << "." << endl);
+    StepList[(*iter)] = DistanceList[(*iter)].begin();    // stores the step to the next iterator that could be a possible next target
+    PermutationMap[(*iter)] = DistanceList[(*iter)].begin()->second;   // always pick target with the smallest distance
+    DoubleList[DistanceList[(*iter)].begin()->second]++;            // increase this target's source count (>1? not injective)
+    DistanceIterators[(*iter)] = DistanceList[(*iter)].begin();    // and remember which one we picked
+    DoLog(2) && (Log() << Verbose(2) << **iter << " starts with distance " << DistanceList[(*iter)].begin()->first << "." << endl);
   }
 };
@@ -201,42 +195,56 @@
   }
   while ((Potential) > PenaltyConstants[2]) {
+    CalculateDoubleList();
     PrintPermutationMap();
     iter++;
     if (iter == atoms.end()) // round-robin at the end
       iter = atoms.begin();
-    if (DoubleList[DistanceIterators[(*iter)->getNr()]->second->getNr()] <= 1)  // no need to make those injective that aren't
+    if (DoubleList[DistanceIterators[(*iter)]->second] <= 1)  // no need to make those injective that aren't
       continue;
     // now, try finding a new one
     Potential = TryNextNearestNeighbourForInjectivePermutation((*iter), Potential);
   }
-  for (int i=atoms.size(); i--;) // now each single entry in the DoubleList should be <=1
-    if (DoubleList[i] > 1) {
+  for (molecule::atomSet::const_iterator iter = atoms.begin(); iter != atoms.end(); ++iter) {
+    // now each single entry in the DoubleList should be <=1
+    if (DoubleList[*iter] > 1) {
       DoeLog(0) && (eLog()<< Verbose(0) << "Failed to create an injective PermutationMap!" << endl);
       performCriticalExit();
     }
+  }
   DoLog(1) && (Log() << Verbose(1) << "done." << endl);
 };
 
+unsigned int MinimiseConstrainedPotential::CalculateDoubleList()
+{
+  for (molecule::atomSet::const_iterator iter = atoms.begin(); iter != atoms.end(); ++iter)
+    DoubleList[*iter] = 0;
+  unsigned int doubles = 0;
+  for (molecule::atomSet::const_iterator iter = atoms.begin(); iter != atoms.end(); ++iter)
+      DoubleList[ PermutationMap[*iter] ]++;
+  for (molecule::atomSet::const_iterator iter = atoms.begin(); iter != atoms.end(); ++iter)
+    if (DoubleList[*iter] > 1)
+      doubles++;
+  if (doubles >0)
+    DoLog(2) && (Log() << Verbose(2) << "Found " << doubles << " Doubles." << endl);
+  return doubles;
+};
+
 void MinimiseConstrainedPotential::PrintPermutationMap() const
 {
   stringstream zeile1, zeile2;
-  const int AtomCount = atoms.size();
-  int *DoubleList = new int[AtomCount];
-  for(int i=0;i<AtomCount;i++)
-    DoubleList[i] = 0;
   int doubles = 0;
   zeile1 << "PermutationMap: ";
   zeile2 << "                ";
-  for (int i=0;i<AtomCount;i++) {
-    DoubleList[PermutationMap[i]->getNr()]++;
-    zeile1 << i << " ";
-    zeile2 << PermutationMap[i]->getNr() << " ";
-  }
-  for (int i=0;i<AtomCount;i++)
-    if (DoubleList[i] > 1)
-    doubles++;
+  for (molecule::atomSet::const_iterator iter = atoms.begin(); iter != atoms.end(); ++iter) {
+    zeile1 << (*iter)->getName() << " ";
+    zeile2 << (PermutationMap[*iter])->getName() << " ";
+  }
+  for (molecule::atomSet::const_iterator iter = atoms.begin(); iter != atoms.end(); ++iter) {
+    std::map<atom *, unsigned int>::const_iterator value_iter = DoubleList.find(*iter);
+    if (value_iter->second > (unsigned int)1)
+      doubles++;
+  }
   if (doubles >0)
     DoLog(2) && (Log() << Verbose(2) << "Found " << doubles << " Doubles." << endl);
-  delete[](DoubleList);
 //  Log() << Verbose(2) << zeile1.str() << endl << zeile2.str() << endl;
 };
@@ -250,5 +258,5 @@
   for (molecule::atomSet::const_iterator iter = atoms.begin(); iter != atoms.end(); ++iter) {
     // first term: distance to target
-    Runner = PermutationMap[(*iter)->getNr()];   // find target point
+    Runner = PermutationMap[(*iter)];   // find target point
     tmp = ((*iter)->getPositionAtStep(startstep).distance(Runner->getPositionAtStep(endstep)));
     tmp *= IsAngstroem ? 1. : 1./AtomicLengthToAngstroem;
@@ -269,17 +277,17 @@
 {
   double Potential = 0;
-  DistanceMap::iterator NewBase = DistanceIterators[Walker->getNr()];  // store old base
+  DistanceMap::iterator NewBase = DistanceIterators[Walker];  // store old base
   do {
     NewBase++;  // take next further distance in distance to targets list that's a target of no one
-  } while ((DoubleList[NewBase->second->getNr()] != 0) && (NewBase != DistanceList[Walker->getNr()]->end()));
-  if (NewBase != DistanceList[Walker->getNr()]->end()) {
-    PermutationMap[Walker->getNr()] = NewBase->second;
+  } while ((DoubleList[NewBase->second] != 0) && (NewBase != DistanceList[Walker].end()));
+  if (NewBase != DistanceList[Walker].end()) {
+    PermutationMap[Walker] = NewBase->second;
     Potential = fabs(ConstrainedPotential());
     if (Potential > OldPotential) { // undo
-      PermutationMap[Walker->getNr()] = DistanceIterators[Walker->getNr()]->second;
+      PermutationMap[Walker] = DistanceIterators[Walker]->second;
     } else {  // do
-      DoubleList[DistanceIterators[Walker->getNr()]->second->getNr()]--;  // decrease the old entry in the doubles list
-      DoubleList[NewBase->second->getNr()]++;    // increase the old entry in the doubles list
-      DistanceIterators[Walker->getNr()] = NewBase;
+      DoubleList[DistanceIterators[Walker]->second]--;  // decrease the old entry in the doubles list
+      DoubleList[NewBase->second]++;    // increase the old entry in the doubles list
+      DistanceIterators[Walker] = NewBase;
       OldPotential = Potential;
       DoLog(3) && (Log() << Verbose(3) << "Found a new permutation, new potential is " << OldPotential << "." << endl);
@@ -293,5 +301,5 @@
   double result = 0.;
   for (molecule::atomSet::const_iterator iter = atoms.begin(); iter != atoms.end(); ++iter) {
-    if ((PermutationMap[Walker->getNr()] == PermutationMap[(*iter)->getNr()]) && (Walker->getNr() < (*iter)->getNr())) {
+    if ((PermutationMap[Walker] == PermutationMap[(*iter)]) && (Walker < (*iter))) {
   //    atom *Sprinter = PermutationMap[Walker->nr];
   //        Log() << Verbose(0) << *Walker << " and " << *(*iter) << " are heading to the same target at ";
@@ -317,9 +325,9 @@
       break;
     // determine normalized trajectories direction vector (n1, n2)
-    Sprinter = PermutationMap[Walker->getNr()];   // find first target point
+    Sprinter = PermutationMap[Walker];   // find first target point
     trajectory1 = Sprinter->getPositionAtStep(endstep) - Walker->getPositionAtStep(startstep);
     trajectory1.Normalize();
     Norm1 = trajectory1.Norm();
-    Sprinter = PermutationMap[(*iter)->getNr()];   // find second target point
+    Sprinter = PermutationMap[(*iter)];   // find second target point
     trajectory2 = Sprinter->getPositionAtStep(endstep) - (*iter)->getPositionAtStep(startstep);
     trajectory2.Normalize();
@@ -329,5 +337,5 @@
       tmp = Walker->getPositionAtStep(startstep).distance((*iter)->getPositionAtStep(startstep));
     } else if (Norm1 < MYEPSILON) {
-      Sprinter = PermutationMap[Walker->getNr()];   // find first target point
+      Sprinter = PermutationMap[Walker];   // find first target point
       trajectory1 = Sprinter->getPositionAtStep(endstep) - (*iter)->getPositionAtStep(startstep);
       trajectory2 *= trajectory1.ScalarProduct(trajectory2); // trajectory2 is scaled to unity, hence we don't need to divide by anything
@@ -335,5 +343,5 @@
       tmp = trajectory1.Norm();  // remaining norm is distance
     } else if (Norm2 < MYEPSILON) {
-      Sprinter = PermutationMap[(*iter)->getNr()];   // find second target point
+      Sprinter = PermutationMap[(*iter)];   // find second target point
       trajectory2 = Sprinter->getPositionAtStep(endstep) - Walker->getPositionAtStep(startstep);  // copy second offset
       trajectory1 *= trajectory2.ScalarProduct(trajectory1); // trajectory1 is scaled to unity, hence we don't need to divide by anything
@@ -402,5 +410,5 @@
   DoLog(1) && (Log() << Verbose(1) << "Calculating forces and adding onto ForceMatrix ... " << endl);
   for(molecule::atomSet::const_iterator iter = atoms.begin(); iter != atoms.end(); ++iter) {
-    atom *Sprinter = PermutationMap[(*iter)->getNr()];
+    atom *Sprinter = PermutationMap[(*iter)];
     // set forces
     for (int i=NDIM;i++;)
Index: src/Dynamics/MinimiseConstrainedPotential.hpp
===================================================================
--- src/Dynamics/MinimiseConstrainedPotential.hpp	(revision fe6f202b08255305d13ec735eddd6dc5bc8cdcf4)
+++ src/Dynamics/MinimiseConstrainedPotential.hpp	(revision e355762f9447acb9c234eadd3df6adc314e5dee7)
@@ -27,5 +27,16 @@
 {
 public:
-  MinimiseConstrainedPotential(molecule::atomSet &_atoms);
+  /** Constructor.
+   *
+   * @param _atoms set of atoms to operate on
+   * \param _PermutationMap on return: mapping between the atom label of the initial and the final configuration
+   * @return
+   */
+  MinimiseConstrainedPotential(molecule::atomSet &_atoms, std::map<atom*, atom *> &_PermutationMap);
+
+  /** Destructor.
+   *
+   * @return
+   */
   ~MinimiseConstrainedPotential();
 
@@ -46,6 +57,4 @@
    *     right direction).
    *  -# Finally, we calculate the potential energy and return.
-   * \param *out output stream for debugging
-   * \param **PermutationMap on return: mapping between the atom label of the initial and the final configuration
    * \param startstep current MD step giving initial position between which and \a endstep we perform the constrained MD (as further steps are always concatenated)
    * \param endstep step giving final position in constrained MD
@@ -57,11 +66,7 @@
    * \bug this all is not O(N log N) but O(N^2)
    */
-  double operator()(atom **&PermutationMap, int startstep, int endstep, bool IsAngstroem);
+  double operator()(int startstep, int endstep, bool IsAngstroem);
 
   /** Evaluates the (distance-related part) of the constrained potential for the constrained forces.
-   * \param *out output stream for debugging
-   * \param startstep current MD step giving initial position between which and \a endstep we perform the constrained MD (as further steps are always concatenated)
-   * \param endstep step giving final position in constrained MD
-   * \param **PermutationMap mapping between the atom label of the initial and the final configuration
    * \param *Force ForceMatrix containing force vectors from the external energy functional minimisation.
    * \todo the constant for the constrained potential distance part is hard-coded independently of the hard-coded value in MinimiseConstrainedPotential()
@@ -77,36 +82,29 @@
   int startstep; //!< start configuration (MDStep in atom::trajectory)
   int endstep; //!< end configuration (MDStep in atom::trajectory)
-  atom **PermutationMap; //!< gives target ptr for each atom, array of size molecule::AtomCount (this is "x" in \f$ V^{con}(x) \f$ )
-  DistanceMap **DistanceList; //!< distance list of each atom to each atom
-  DistanceMap::iterator *StepList; //!< iterator to ascend through NearestNeighbours \a **DistanceList
-  int *DoubleList; //!< count of which sources want to move to this target, basically the injective measure (>1 -> not injective)
-  DistanceMap::iterator *DistanceIterators; //!< marks which was the last picked target as injective candidate with smallest distance
+  std::map<atom*, atom *> &PermutationMap; //!< gives target ptr for each atom, array of size molecule::AtomCount (this is "x" in \f$ V^{con}(x) \f$ )
+  std::map<atom *, DistanceMap> DistanceList; //!< distance list of each atom to each atom
+  std::map<atom *, DistanceMap::iterator> StepList; //!< iterator to ascend through NearestNeighbours \a **DistanceList
+  std::map<atom *, unsigned int> DoubleList; //!< count of which sources want to move to this target, basically the injective measure (>1 -> not injective)
+  std::map<atom *, DistanceMap::iterator> DistanceIterators; //!< marks which was the last picked target as injective candidate with smallest distance
   bool IsAngstroem; //!< whether coordinates are in angstroem (true) or bohrradius (false)
   double *PenaltyConstants; //!<  penalty constant in front of each term
 
   /** \f$O(N^2)\f$ operation of calculation distance between each atom pair and putting into DistanceList.
-   * \param *mol molecule to scan distances in
-   * \param &Params constrained potential parameters
    */
   void FillDistanceList();
 
   /** Initialize lists.
-   * \param *out output stream for debugging
-   * \param *mol molecule to scan distances in
-   * \param &Params constrained potential parameters
    */
   void CreateInitialLists();
 
   /** Permutes \a **&PermutationMap until the penalty is below constants[2].
-   * \param *out output stream for debugging
-   * \param *mol molecule to scan distances in
-   * \param &Params constrained potential parameters
    */
   void MakeInjectivePermutation();
 
+  /** Calculates the number of doubles in PermutationMap.
+   */
+  unsigned int CalculateDoubleList();
+
   /** Print the current permutation map.
-   * \param *out output stream for debugging
-   * \param &Params constrained potential parameters
-   * \param AtomCount number of atoms
    */
   void PrintPermutationMap() const;
@@ -121,6 +119,4 @@
    * \f$n_3\f$ is the normal vector to both directions. \f$C\f$ would be the minimum distance between the two lines.
    * \sa molecule::MinimiseConstrainedPotential(), molecule::VerletForceIntegration()
-   * \param *out output stream for debugging
-   * \param &Params constrained potential parameters
    * \return potential energy
    * \note This routine is scaling quadratically which is not optimal.
@@ -130,9 +126,6 @@
 
   /** Try the next nearest neighbour in order to make the permutation map injective.
-   * \param *out output stream for debugging
-   * \param *mol molecule
    * \param *Walker atom to change its target
    * \param &OldPotential old value of constraint potential to see if we do better with new target
-   * \param &Params constrained potential parameters
    */
   double TryNextNearestNeighbourForInjectivePermutation(atom *Walker, double &OldPotential);
@@ -140,6 +133,4 @@
   /** Penalizes atoms heading to same target.
    * \param *Walker atom to check against others
-   * \param *mol molecule with other atoms
-   * \param &Params constrained potential parameters
    * \return \a penalty times the number of equal targets
    */
@@ -148,6 +139,4 @@
   /** Penalizes long trajectories.
    * \param *Walker atom to check against others
-   * \param *mol molecule with other atoms
-   * \param &Params constraint potential parameters
    * \return penalty times each distance
    */
Index: src/molecule_dynamics.cpp
===================================================================
--- src/molecule_dynamics.cpp	(revision fe6f202b08255305d13ec735eddd6dc5bc8cdcf4)
+++ src/molecule_dynamics.cpp	(revision e355762f9447acb9c234eadd3df6adc314e5dee7)
@@ -85,10 +85,9 @@
   if (configuration.DoConstrainedMD) {
     // calculate forces and potential
-    atom **PermutationMap = NULL;
-    MinimiseConstrainedPotential Minimiser(atoms);
+    std::map<atom*, atom *> PermutationMap;
+    MinimiseConstrainedPotential Minimiser(atoms, PermutationMap);
     //double ConstrainedPotentialEnergy =
-    Minimiser(PermutationMap, configuration.DoConstrainedMD, 0, configuration.GetIsAngstroem());
+    Minimiser(configuration.DoConstrainedMD, 0, configuration.GetIsAngstroem());
     Minimiser.EvaluateConstrainedForces(&Force);
-    delete[](PermutationMap);
   }
 
