Reading individual sites from moleculeProperties dictionary entries, finished.

applications/solvers/molecularDynamics/mdFoam/mdFoam.C

@@ -43,6 +43,10 @@ int main(int argc, char *argv[])
 
     potential pot(mesh);
 
+    const pairPotentialList& pairPot(pot.pairPotentials());
+
+    Info<< pairPot.energy(0, 0, 0.45e-9) << endl;
+
     moleculeCloud molecules(mesh, pot);
 
     Info << "\nStarting time loop\n" << endl;

applications/solvers/molecularDynamics/mdWaterTest/mdWaterTest.C

@@ -85,6 +85,11 @@ int main(int argc, char *argv[])
           + 2.0 * mSites[0] * (pSites[0].y() * pSites[0].y() + pSites[0].x() * pSites[0].x())
     );
 
+    Info<< m
+        << nl << qSites
+        << nl << pSites
+        << nl << I << endl;
+
     vector p1(0, 0, 0);
 
     vector v1(100, 0, 0);

src/lagrangian/molecularDynamics/molecule/molecule/molecule.C

@@ -40,8 +40,8 @@ Foam::tensor Foam::molecule::rotationTensor(scalar deltaT) const
         tensor
         (
             1, 0, 0,
-            0, Foam::cos(phi1), Foam::sin(phi1),
-            0, -Foam::sin(phi1), Foam::cos(phi1)
+            0, Foam::cos(phi1), -Foam::sin(phi1),
+            0, Foam::sin(phi1), Foam::cos(phi1)
         )
     );
 
@@ -51,9 +51,9 @@ Foam::tensor Foam::molecule::rotationTensor(scalar deltaT) const
     (
         tensor
         (
-            Foam::cos(phi2), 0, -Foam::sin(phi2),
+            Foam::cos(phi2), 0, Foam::sin(phi2),
             0, 1, 0,
-            Foam::sin(phi2), 0, Foam::cos(phi2)
+            -Foam::sin(phi2), 0, Foam::cos(phi2)
         )
     );
 
@@ -63,13 +63,13 @@ Foam::tensor Foam::molecule::rotationTensor(scalar deltaT) const
     (
         tensor
         (
-            Foam::cos(phi3), Foam::sin(phi3), 0,
-            -Foam::sin(phi3), Foam::cos(phi3), 0,
+            Foam::cos(phi3), -Foam::sin(phi3), 0,
+            Foam::sin(phi3), Foam::cos(phi3), 0,
             0, 0, 1
         )
     );
 
-    return (U1.T() & U2.T() & U3.T() & U2.T() & U1.T());
+    return (U1 & U2 & U3 & U2 & U1);
 }
 
 

src/lagrangian/molecularDynamics/molecule/molecule/moleculeI.H

@@ -67,8 +67,6 @@ inline Foam::molecule::constantProperties::constantProperties
     mass_(readScalar(dict.lookup("mass")))
 {
     checkSiteListSizes();
-
-    Info<< nl << "siteIds not set yet." << endl;
 }
 
 

src/lagrangian/molecularDynamics/molecule/moleculeCloud/moleculeCloud.C

@@ -54,6 +54,19 @@ void Foam::moleculeCloud::buildConstProps()
 
     const List<word>& allSiteIdNames(pot_.allSiteIdNames());
 
+    IOdictionary moleculePropertiesDict
+    (
+        IOobject
+        (
+            "moleculeProperties",
+            mesh_.time().constant(),
+            mesh_,
+            IOobject::MUST_READ,
+            IOobject::NO_WRITE,
+            false
+        )
+    );
+
     forAll(idList, i)
     {
         const word& id(idList[i]);
@@ -83,14 +96,6 @@ void Foam::moleculeCloud::buildConstProps()
         constProp = molecule::constantProperties(molDict);
 
         constProp.siteIds() = siteIds;
-
-        Info<< constPropList_[i].siteReferencePositions()
-            << nl << constPropList_[i].siteCharges()
-            << nl << constPropList_[i].siteIds()
-            << nl << constPropList_[i].momentOfInertia()
-            << nl << constPropList_[i].mass()
-            << nl << constPropList_[i].nSites()
-            << endl;
     }
 }
 

src/lagrangian/molecularDynamics/potential/Make/files

@@ -15,6 +15,7 @@ $(pairPotential)/derived/maitlandSmith/maitlandSmith.C
 $(pairPotential)/derived/azizChen/azizChen.C
 $(pairPotential)/derived/exponentialRepulsion/exponentialRepulsion.C
 $(pairPotential)/derived/coulomb/coulomb.C
+$(pairPotential)/derived/noInteraction/noInteraction.C
 
 energyScalingFunction = energyScalingFunction
 

src/lagrangian/molecularDynamics/potential/potential/potential.C

@@ -26,27 +26,9 @@ License
 
 #include "potential.H"
 
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
-Foam::potential::potential(const polyMesh& mesh)
-:
-    mesh_(mesh)
-{}
-
-
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-Foam::potential::~potential()
-{}
-
-
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
-
-
-void Foam::potential::potential::readPotentialDict
-(
-    const List<word>& allSiteIdNames
-)
+void Foam::potential::potential::readPotentialDict()
 {
     Info<< nl <<  "Reading potential dictionary:" << endl;
 
@@ -79,9 +61,9 @@ void Foam::potential::potential::readPotentialDict
 
     DynamicList<word> allSiteIdNames;
 
-    forAll(idList, i)
+    forAll(idList_, i)
     {
-        const word& id(idList[i]);
+        const word& id(idList_[i]);
 
         if (!moleculePropertiesDict.found(id))
         {
@@ -105,10 +87,9 @@ void Foam::potential::potential::readPotentialDict
         }
     }
 
-    // CREATE allSiteIdNames_ DATA MEMBER
     allSiteIdNames_.transfer(allSiteIdNames.shrink());
 
-    Info<< nl << "Unique site ids found: " << allSiteIdNames << endl;
+    Info<< nl << "Unique site ids found: " << allSiteIdNames_ << endl;
 
     List<word> tetherIdList(0);
 
@@ -143,6 +124,14 @@ void Foam::potential::potential::readPotentialDict
         forAll(removalOrder_, rO)
         {
             removalOrder_[rO] = findIndex(idList_, remOrd[rO]);
+
+            if (removalOrder_[rO] == -1)
+            {
+                FatalErrorIn("potentials.C") << nl
+                    << "removalOrder entry: " << remOrd[rO]
+                    << " not found in idList."
+                    << nl << abort(FatalError);
+            }
         }
     }
 
@@ -158,7 +147,12 @@ void Foam::potential::potential::readPotentialDict
 
     const dictionary& pairDict = potentialDict.subDict("pair");
 
-    pairPotentials_.buildPotentials(idList_, pairDict, mesh_);
+    pairPotentials_.buildPotentials
+    (
+        allSiteIdNames_,
+        pairDict,
+        mesh_
+    );
 
     // *************************************************************************
     // Tether potentials
@@ -174,7 +168,12 @@ void Foam::potential::potential::readPotentialDict
 
         const dictionary& tetherDict = potentialDict.subDict("tether");
 
-        tetherPotentials_.buildPotentials(idList_, tetherDict, tetherIdList);
+        tetherPotentials_.buildPotentials
+        (
+            allSiteIdNames_,
+            tetherDict,
+            tetherIdList
+        );
     }
 
     // *************************************************************************
@@ -201,6 +200,21 @@ void Foam::potential::potential::readPotentialDict
     Info << nl << tab << "gravity = " << gravity_ << endl;
 }
 
-// * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::potential::potential(const polyMesh& mesh)
+:
+    mesh_(mesh)
+{
+    readPotentialDict();
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::potential::~potential()
+{}
+
 
 // ************************************************************************* //

src/lagrangian/molecularDynamics/potential/potential/potential.H

@@ -59,6 +59,8 @@ class potential
 
         List<word> idList_;
 
+        List<word> allSiteIdNames_;
+
         scalar potentialEnergyLimit_;
 
         labelList removalOrder_;
@@ -71,6 +73,8 @@ class potential
 
     // Private Member Functions
 
+        void readPotentialDict();
+
         //- Disallow default bitwise copy construct
         potential(const potential&);
 
@@ -98,6 +102,8 @@ public:
 
             inline const List<word>& idList() const;
 
+            inline const List<word>& allSiteIdNames() const;
+
             inline scalar potentialEnergyLimit() const;;
 
             inline label nPairPotentials() const;

src/lagrangian/molecularDynamics/potential/potential/potentialI.H

@@ -40,6 +40,12 @@ inline const Foam::List<Foam::word>& Foam::potential::idList() const
 }
 
 
+inline const Foam::List<Foam::word>& Foam::potential::allSiteIdNames() const
+{
+    return allSiteIdNames_;
+}
+
+
 inline Foam::scalar Foam::potential::potentialEnergyLimit() const
 {
     return potentialEnergyLimit_;