From a20caaa2affd3fa8bd0aad15e01c6605dc29c9f9 Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Fri, 22 Jul 2011 10:34:34 +0100
Subject: [PATCH] ENH:molecule: reverted to old molecule.
Removed MoleculeCloud, controllers, monoatomic, polyatomic molecules
---
.../Templates/MoleculeCloud/MoleculeCloud.C | 1265 -----------------
.../Templates/MoleculeCloud/MoleculeCloud.H | 252 ----
.../Templates/MoleculeCloud/MoleculeCloudI.H | 416 ------
.../baseClasses/moleculeCloud/moleculeCloud.C | 48 -
.../baseClasses/moleculeCloud/moleculeCloud.H | 83 --
.../derived/monoatomicCloud/monoatomicCloud.H | 52 -
.../derived/polyatomicCloud/polyatomicCloud.H | 52 -
.../basic/controllers/controllers.C | 573 --------
.../basic/controllers/controllers.H | 163 ---
.../basic/controllers/controllersI.H | 64 -
.../basic/fluxController/fluxController.C | 524 -------
.../basic/fluxController/fluxController.H | 272 ----
.../basic/stateController/stateController.C | 490 -------
.../basic/stateController/stateController.H | 270 ----
.../molecules/constPropSite/constPropSite.C | 114 --
.../molecules/constPropSite/constPropSite.H | 185 ---
.../molecules/constPropSite/constPropSiteI.H | 139 --
.../molecules/monoatomic/monoatomic.C | 199 ---
.../molecules/monoatomic/monoatomic.H | 418 ------
.../molecules/monoatomic/monoatomicI.H | 328 -----
.../molecules/monoatomic/monoatomicIO.C | 305 ----
.../molecules/polyatomic/polyatomic.C | 320 -----
.../molecules/polyatomic/polyatomic.H | 483 -------
.../molecules/polyatomic/polyatomicI.H | 653 ---------
.../molecules/polyatomic/polyatomicIO.C | 414 ------
25 files changed, 8082 deletions(-)
delete mode 100644 src/lagrangian/molecularDynamics/molecule/clouds/Templates/MoleculeCloud/MoleculeCloud.C
delete mode 100644 src/lagrangian/molecularDynamics/molecule/clouds/Templates/MoleculeCloud/MoleculeCloud.H
delete mode 100644 src/lagrangian/molecularDynamics/molecule/clouds/Templates/MoleculeCloud/MoleculeCloudI.H
delete mode 100644 src/lagrangian/molecularDynamics/molecule/clouds/baseClasses/moleculeCloud/moleculeCloud.C
delete mode 100644 src/lagrangian/molecularDynamics/molecule/clouds/baseClasses/moleculeCloud/moleculeCloud.H
delete mode 100644 src/lagrangian/molecularDynamics/molecule/clouds/derived/monoatomicCloud/monoatomicCloud.H
delete mode 100644 src/lagrangian/molecularDynamics/molecule/clouds/derived/polyatomicCloud/polyatomicCloud.H
delete mode 100644 src/lagrangian/molecularDynamics/molecule/controllers/basic/controllers/controllers.C
delete mode 100644 src/lagrangian/molecularDynamics/molecule/controllers/basic/controllers/controllers.H
delete mode 100644 src/lagrangian/molecularDynamics/molecule/controllers/basic/controllers/controllersI.H
delete mode 100644 src/lagrangian/molecularDynamics/molecule/controllers/basic/fluxController/fluxController.C
delete mode 100644 src/lagrangian/molecularDynamics/molecule/controllers/basic/fluxController/fluxController.H
delete mode 100644 src/lagrangian/molecularDynamics/molecule/controllers/basic/stateController/stateController.C
delete mode 100644 src/lagrangian/molecularDynamics/molecule/controllers/basic/stateController/stateController.H
delete mode 100644 src/lagrangian/molecularDynamics/molecule/molecules/constPropSite/constPropSite.C
delete mode 100644 src/lagrangian/molecularDynamics/molecule/molecules/constPropSite/constPropSite.H
delete mode 100644 src/lagrangian/molecularDynamics/molecule/molecules/constPropSite/constPropSiteI.H
delete mode 100644 src/lagrangian/molecularDynamics/molecule/molecules/monoatomic/monoatomic.C
delete mode 100644 src/lagrangian/molecularDynamics/molecule/molecules/monoatomic/monoatomic.H
delete mode 100644 src/lagrangian/molecularDynamics/molecule/molecules/monoatomic/monoatomicI.H
delete mode 100644 src/lagrangian/molecularDynamics/molecule/molecules/monoatomic/monoatomicIO.C
delete mode 100644 src/lagrangian/molecularDynamics/molecule/molecules/polyatomic/polyatomic.C
delete mode 100644 src/lagrangian/molecularDynamics/molecule/molecules/polyatomic/polyatomic.H
delete mode 100644 src/lagrangian/molecularDynamics/molecule/molecules/polyatomic/polyatomicI.H
delete mode 100644 src/lagrangian/molecularDynamics/molecule/molecules/polyatomic/polyatomicIO.C
diff --git a/src/lagrangian/molecularDynamics/molecule/clouds/Templates/MoleculeCloud/MoleculeCloud.C b/src/lagrangian/molecularDynamics/molecule/clouds/Templates/MoleculeCloud/MoleculeCloud.C
deleted file mode 100644
index c6d54036f7e..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/clouds/Templates/MoleculeCloud/MoleculeCloud.C
+++ /dev/null
@@ -1,1265 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2008-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-\*----------------------------------------------------------------------------*/
-
-#include "MoleculeCloud.H"
-#include "fvMesh.H"
-#include "mathematicalConstants.H"
-
-using namespace Foam::constant::mathematical;
-
-
-// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
-
-template<class MoleculeType>
-void Foam::MoleculeCloud<MoleculeType>::buildConstProps()
-{
- const List<word>& idList(pot_.idList());
-
- constPropList_.setSize(idList.size());
-
- const List<word>& siteIdList(pot_.siteIdList());
-
- IOdictionary molPropertiesDict
- (
- IOobject
- (
- this->name() + "Properties",
- mesh_.time().constant(),
- mesh_,
- IOobject::MUST_READ_IF_MODIFIED,
- IOobject::NO_WRITE,
- false
- )
- );
-
- Info<< nl << "Reading " << molPropertiesDict.name() << endl;
-
- forAll(idList, i)
- {
- const word& id = idList[i];
-
- const dictionary& molDict(molPropertiesDict.subDict(id));
-
- List<word> siteIdNames = molDict.lookup("siteIds");
-
- List<label> siteIds(siteIdNames.size());
-
- forAll(siteIdNames, sI)
- {
- const word& siteId = siteIdNames[sI];
-
- siteIds[sI] = findIndex(siteIdList, siteId);
-
- if (siteIds[sI] == -1)
- {
- FatalErrorIn
- (
- "void Foam::MoleculeCloud<MoleculeType>::buildConstProps()"
- )
- << siteId << " site not found."
- << nl << abort(FatalError);
- }
- }
-
- constPropList_[i] = typename MoleculeType::constantProperties
- (
- molDict,
- siteIds
- );
- }
-}
-
-
-template<class MoleculeType>
-void Foam::MoleculeCloud<MoleculeType>::setSiteSizesAndPositions()
-{
- forAllIter(typename MoleculeCloud<MoleculeType>, *this, mol)
- {
- const typename MoleculeType::constantProperties& cP
- (
- constProps(mol().id())
- );
-
- mol().setSiteSizes(cP.nSites());
-
- mol().setSitePositions(cP);
- }
-}
-
-
-template<class MoleculeType>
-void Foam::MoleculeCloud<MoleculeType>::buildCellOccupancy()
-{
- forAll(cellOccupancy_, cO)
- {
- cellOccupancy_[cO].clear();
- }
-
- forAllIter(typename MoleculeCloud<MoleculeType>, *this, mol)
- {
- cellOccupancy_[mol().cell()].append(&mol());
- }
-
- forAll(cellOccupancy_, cO)
- {
- cellOccupancy_[cO].shrink();
- }
-}
-
-
-template<class MoleculeType>
-void Foam::MoleculeCloud<MoleculeType>::calculatePairForce()
-{
- PstreamBuffers pBufs(Pstream::nonBlocking);
-
- // Start sending referred data
- il_.sendReferredData(cellOccupancy(), pBufs);
-
- MoleculeType* molI = NULL;
- MoleculeType* molJ = NULL;
-
- {
- // Real-Real interactions
-
- const labelListList& dil = il_.dil();
-
- forAll(dil, d)
- {
- forAll(cellOccupancy_[d],cellIMols)
- {
- molI = cellOccupancy_[d][cellIMols];
-
- forAll(dil[d], interactingCells)
- {
- List<MoleculeType*> cellJ =
- cellOccupancy_[dil[d][interactingCells]];
-
- forAll(cellJ, cellJMols)
- {
- molJ = cellJ[cellJMols];
-
- evaluatePair(*molI, *molJ);
- }
- }
-
- forAll(cellOccupancy_[d], cellIOtherMols)
- {
- molJ = cellOccupancy_[d][cellIOtherMols];
-
- if (molJ > molI)
- {
- evaluatePair(*molI, *molJ);
- }
- }
- }
- }
- }
-
- // Receive referred data
- il_.receiveReferredData(pBufs);
-
- {
- // Real-Referred interactions
-
- const labelListList& ril = il_.ril();
-
- List<IDLList<MoleculeType> >& referredMols = il_.referredParticles();
-
- forAll(ril, r)
- {
- const List<label>& realCells = ril[r];
-
- IDLList<MoleculeType>& refMols = referredMols[r];
-
- forAllIter
- (
- typename IDLList<MoleculeType>,
- refMols,
- refMol
- )
- {
- forAll(realCells, rC)
- {
- List<MoleculeType*> cellI = cellOccupancy_[realCells[rC]];
-
- forAll(cellI, cellIMols)
- {
- molI = cellI[cellIMols];
-
- evaluatePair(*molI, refMol());
- }
- }
- }
- }
- }
-}
-
-
-template<class MoleculeType>
-void Foam::MoleculeCloud<MoleculeType>::calculateTetherForce()
-{
- const tetherPotentialList& tetherPot(pot_.tetherPotentials());
-
- forAllIter(typename MoleculeCloud<MoleculeType>, *this, mol)
- {
- if (mol().tethered())
- {
- vector rIT = mol().position() - mol().specialPosition();
-
- label idI = mol().id();
-
- scalar massI = constProps(idI).mass();
-
- vector fIT = tetherPot.force(idI, rIT);
-
- mol().a() += fIT/massI;
-
- mol().potentialEnergy() += tetherPot.energy(idI, rIT);
-
- // What to do here?
- // mol().rf() += rIT*fIT;
- }
- }
-}
-
-
-template<class MoleculeType>
-void Foam::MoleculeCloud<MoleculeType>::calculateExternalForce()
-{
- forAllIter(typename MoleculeCloud<MoleculeType>, *this, mol)
- {
- mol().a() += pot_.gravity();
- }
-}
-
-
-template<class MoleculeType>
-void Foam::MoleculeCloud<MoleculeType>::removeHighEnergyOverlaps()
-{
- Info<< nl << "Removing high energy overlaps, limit = "
- << pot_.potentialEnergyLimit()
- << nl << "Removal order:";
-
- forAll(pot_.removalOrder(), rO)
- {
- Info<< ' ' << pot_.idList()[pot_.removalOrder()[rO]];
- }
-
- Info<< nl ;
-
- label initialSize = this->size();
-
- buildCellOccupancy();
-
- // Real-Real interaction
-
- MoleculeType* molI = NULL;
- MoleculeType* molJ = NULL;
-
- {
- DynamicList<MoleculeType*> molsToDelete;
-
- const labelListList& dil(il_.dil());
-
- forAll(dil, d)
- {
- forAll(cellOccupancy_[d], cellIMols)
- {
- molI = cellOccupancy_[d][cellIMols];
-
- forAll(dil[d], interactingCells)
- {
- List<MoleculeType*> cellJ =
- cellOccupancy_[dil[d][interactingCells]];
-
- forAll(cellJ, cellJMols)
- {
- molJ = cellJ[cellJMols];
-
- if (evaluatePotentialLimit(*molI, *molJ))
- {
- label idI = molI->id();
-
- label idJ = molJ->id();
-
- if
- (
- idI == idJ
- || findIndex(pot_.removalOrder(), idJ)
- < findIndex(pot_.removalOrder(), idI)
- )
- {
- if (findIndex(molsToDelete, molJ) == -1)
- {
- molsToDelete.append(molJ);
- }
- }
- else if (findIndex(molsToDelete, molI) == -1)
- {
- molsToDelete.append(molI);
- }
- }
- }
- }
-
- forAll(cellOccupancy_[d], cellIOtherMols)
- {
- molJ = cellOccupancy_[d][cellIOtherMols];
-
- if (molJ > molI)
- {
- if (evaluatePotentialLimit(*molI, *molJ))
- {
- label idI = molI->id();
-
- label idJ = molJ->id();
-
- if
- (
- idI == idJ
- || findIndex(pot_.removalOrder(), idJ)
- < findIndex(pot_.removalOrder(), idI)
- )
- {
- if (findIndex(molsToDelete, molJ) == -1)
- {
- molsToDelete.append(molJ);
- }
- }
- else if (findIndex(molsToDelete, molI) == -1)
- {
- molsToDelete.append(molI);
- }
- }
- }
- }
- }
- }
-
- forAll(molsToDelete, mTD)
- {
- deleteParticle(*(molsToDelete[mTD]));
- }
- }
-
- buildCellOccupancy();
-
- PstreamBuffers pBufs(Pstream::nonBlocking);
-
- // Start sending referred data
- il_.sendReferredData(cellOccupancy(), pBufs);
-
- // Receive referred data
- il_.receiveReferredData(pBufs);
-
- // Real-Referred interaction
-
- {
- DynamicList<MoleculeType*> molsToDelete;
-
- const labelListList& ril(il_.ril());
-
- List<IDLList<MoleculeType> >& referredMols = il_.referredParticles();
-
- forAll(ril, r)
- {
- IDLList<MoleculeType>& refMols = referredMols[r];
-
- forAllIter
- (
- typename IDLList<MoleculeType>,
- refMols,
- refMol
- )
- {
- molJ = &refMol();
-
- const List<label>& realCells = ril[r];
-
- forAll(realCells, rC)
- {
- label cellI = realCells[rC];
-
- List<MoleculeType*> cellIMols = cellOccupancy_[cellI];
-
- forAll(cellIMols, cIM)
- {
- molI = cellIMols[cIM];
-
- if (evaluatePotentialLimit(*molI, *molJ))
- {
- label idI = molI->id();
-
- label idJ = molJ->id();
-
- if
- (
- findIndex(pot_.removalOrder(), idI)
- < findIndex(pot_.removalOrder(), idJ)
- )
- {
- if (findIndex(molsToDelete, molI) == -1)
- {
- molsToDelete.append(molI);
- }
- }
- else if
- (
- findIndex(pot_.removalOrder(), idI)
- == findIndex(pot_.removalOrder(), idJ)
- )
- {
- // Remove one of the molecules
- // arbitrarily, assuring that a
- // consistent decision is made for
- // both real-referred pairs.
-
- if (molI->origId() > molJ->origId())
- {
- if (findIndex(molsToDelete, molI) == -1)
- {
- molsToDelete.append(molI);
- }
- }
- }
- }
- }
- }
- }
- }
-
- forAll(molsToDelete, mTD)
- {
- deleteParticle(*(molsToDelete[mTD]));
- }
- }
-
- buildCellOccupancy();
-
- // Start sending referred data
- il_.sendReferredData(cellOccupancy(), pBufs);
-
- // Receive referred data
- il_.receiveReferredData(pBufs);
-
- label molsRemoved = initialSize - this->size();
-
- if (Pstream::parRun())
- {
- reduce(molsRemoved, sumOp<label>());
- }
-
- Info<< tab << molsRemoved << " molecules removed" << endl;
-}
-
-
-template<class MoleculeType>
-void Foam::MoleculeCloud<MoleculeType>::initialiseMolecules
-(
- const dictionary& mdInitialiseDict
-)
-{
- Info<< nl
- << "Initialising molecules in each zone specified in "
- << mdInitialiseDict.name()
- << endl;
-
- const cellZoneMesh& cellZones = mesh_.cellZones();
-
- if (!cellZones.size())
- {
- FatalErrorIn
- (
- "void Foam::MoleculeCloud<MoleculeType>::initialiseMolecules"
- )
- << "No cellZones found in the mesh."
- << abort(FatalError);
- }
-
- forAll(cellZones, z)
- {
- const cellZone& zone(cellZones[z]);
-
- if (zone.size())
- {
- if (!mdInitialiseDict.found(zone.name()))
- {
- Info<< "No specification subDictionary for zone "
- << zone.name() << " found, skipping." << endl;
- }
- else
- {
- const dictionary& zoneDict =
- mdInitialiseDict.subDict(zone.name());
-
- const scalar temperature
- (
- readScalar(zoneDict.lookup("temperature"))
- );
-
- const vector bulkVelocity(zoneDict.lookup("bulkVelocity"));
-
- List<word> latticeIds
- (
- zoneDict.lookup("latticeIds")
- );
-
- List<vector> latticePositions
- (
- zoneDict.lookup("latticePositions")
- );
-
- if (latticeIds.size() != latticePositions.size())
- {
- FatalErrorIn
- (
- "Foam::MoleculeCloud<MoleculeType>::initialiseMolecules"
- )
- << "latticeIds and latticePositions must be the same "
- << " size." << nl
- << abort(FatalError);
- }
-
- diagTensor latticeCellShape
- (
- zoneDict.lookup("latticeCellShape")
- );
-
- scalar latticeCellScale = 0.0;
-
- if (zoneDict.found("numberDensity"))
- {
- scalar numberDensity = readScalar
- (
- zoneDict.lookup("numberDensity")
- );
-
- if (numberDensity < VSMALL)
- {
- WarningIn
- (
- "MoleculeCloud<MoleculeType>::initialiseMolecules"
- )
- << "numberDensity too small, not filling zone "
- << zone.name() << endl;
-
- continue;
- }
-
- latticeCellScale = pow
- (
- latticeIds.size()/(det(latticeCellShape)*numberDensity),
- (1.0/3.0)
- );
- }
- else if (zoneDict.found("massDensity"))
- {
- scalar unitCellMass = 0.0;
-
- forAll(latticeIds, i)
- {
- label id = findIndex(pot_.idList(), latticeIds[i]);
-
- const typename MoleculeType::constantProperties& cP
- (
- constProps(id)
- );
-
- unitCellMass += cP.mass();
- }
-
- scalar massDensity = readScalar
- (
- zoneDict.lookup("massDensity")
- );
-
- if (massDensity < VSMALL)
- {
- WarningIn
- (
- "MoleculeCloud<MoleculeType>::initialiseMolecules"
- )
- << "massDensity too small, not filling zone "
- << zone.name() << endl;
-
- continue;
- }
-
- latticeCellScale = pow
- (
- unitCellMass/(det(latticeCellShape)*massDensity),
- (1.0/3.0)
- );
- }
- else
- {
- FatalErrorIn
- (
- "Foam::MoleculeCloud<MoleculeType>::initialiseMolecules"
- )
- << "massDensity or numberDensity not specified " << nl
- << abort(FatalError);
- }
-
- latticeCellShape *= latticeCellScale;
-
- vector anchor(zoneDict.lookup("anchor"));
-
- bool tethered = false;
-
- if (zoneDict.found("tetherSiteIds"))
- {
- tethered = bool
- (
- List<word>(zoneDict.lookup("tetherSiteIds")).size()
- );
- }
-
- const vector orientationAngles
- (
- zoneDict.lookup("orientationAngles")
- );
-
- scalar phi(orientationAngles.x()*pi/180.0);
-
- scalar theta(orientationAngles.y()*pi/180.0);
-
- scalar psi(orientationAngles.z()*pi/180.0);
-
- const tensor R
- (
- cos(psi)*cos(phi) - cos(theta)*sin(phi)*sin(psi),
- cos(psi)*sin(phi) + cos(theta)*cos(phi)*sin(psi),
- sin(psi)*sin(theta),
- - sin(psi)*cos(phi) - cos(theta)*sin(phi)*cos(psi),
- - sin(psi)*sin(phi) + cos(theta)*cos(phi)*cos(psi),
- cos(psi)*sin(theta),
- sin(theta)*sin(phi),
- - sin(theta)*cos(phi),
- cos(theta)
- );
-
- // Find the optimal anchor position. Finding the approximate
- // mid-point of the zone of cells and snapping to the nearest
- // lattice location.
-
- vector zoneMin = VGREAT*vector::one;
-
- vector zoneMax = -VGREAT*vector::one;
-
- forAll(zone, cell)
- {
- const point cellCentre = mesh_.cellCentres()[zone[cell]];
-
- if (cellCentre.x() > zoneMax.x())
- {
- zoneMax.x() = cellCentre.x();
- }
- if (cellCentre.x() < zoneMin.x())
- {
- zoneMin.x() = cellCentre.x();
- }
- if (cellCentre.y() > zoneMax.y())
- {
- zoneMax.y() = cellCentre.y();
- }
- if (cellCentre.y() < zoneMin.y())
- {
- zoneMin.y() = cellCentre.y();
- }
- if (cellCentre.z() > zoneMax.z())
- {
- zoneMax.z() = cellCentre.z();
- }
- if (cellCentre.z() < zoneMin.z())
- {
- zoneMin.z() = cellCentre.z();
- }
- }
-
- point zoneMid = 0.5*(zoneMin + zoneMax);
-
- point latticeMid = inv(latticeCellShape) & (R.T()
- & (zoneMid - anchor));
-
- point latticeAnchor
- (
- label(latticeMid.x() + 0.5*sign(latticeMid.x())),
- label(latticeMid.y() + 0.5*sign(latticeMid.y())),
- label(latticeMid.z() + 0.5*sign(latticeMid.z()))
- );
-
- anchor += (R & (latticeCellShape & latticeAnchor));
-
- // Continue trying to place molecule as long as at
- // least one molecule is placed in each iteration.
- // The "|| totalZoneMols == 0" condition means that the
- // algorithm will continue if the origin is outside the
- // zone.
-
- label n = 0;
-
- label totalZoneMols = 0;
-
- label molsPlacedThisIteration = 0;
-
- while
- (
- molsPlacedThisIteration != 0
- || totalZoneMols == 0
- )
- {
- label sizeBeforeIteration = this->size();
-
- bool partOfLayerInBounds = false;
-
- if (n == 0)
- {
- // Special treatment is required for the first position,
- // i.e. iteration zero.
-
- labelVector unitCellLatticePosition(0,0,0);
-
- forAll(latticePositions, p)
- {
- label id = findIndex(pot_.idList(), latticeIds[p]);
-
- const vector& latticePosition =
- vector
- (
- unitCellLatticePosition.x(),
- unitCellLatticePosition.y(),
- unitCellLatticePosition.z()
- )
- + latticePositions[p];
-
- point globalPosition =
- anchor
- + (R & (latticeCellShape & latticePosition));
-
- partOfLayerInBounds = mesh_.bounds().contains
- (
- globalPosition
- );
-
- label cell = -1;
- label tetFace = -1;
- label tetPt = -1;
-
- mesh_.findCellFacePt
- (
- globalPosition,
- cell,
- tetFace,
- tetPt
- );
-
- if (findIndex(zone, cell) != -1)
- {
- createMolecule
- (
- globalPosition,
- cell,
- tetFace,
- tetPt,
- id,
- tethered,
- temperature,
- bulkVelocity
- );
- }
- }
- }
- else
- {
- // Place top and bottom caps.
-
- labelVector unitCellLatticePosition(0,0,0);
-
- for
- (
- unitCellLatticePosition.z() = -n;
- unitCellLatticePosition.z() <= n;
- unitCellLatticePosition.z() += 2*n
- )
- {
- for
- (
- unitCellLatticePosition.y() = -n;
- unitCellLatticePosition.y() <= n;
- unitCellLatticePosition.y()++
- )
- {
- for
- (
- unitCellLatticePosition.x() = -n;
- unitCellLatticePosition.x() <= n;
- unitCellLatticePosition.x()++
- )
- {
- forAll(latticePositions, p)
- {
- label id = findIndex
- (
- pot_.idList(),
- latticeIds[p]
- );
-
- const vector& latticePosition =
- vector
- (
- unitCellLatticePosition.x(),
- unitCellLatticePosition.y(),
- unitCellLatticePosition.z()
- )
- + latticePositions[p];
-
- point globalPosition =
- anchor
- + (
- R
- & (
- latticeCellShape
- & latticePosition
- )
- );
-
- partOfLayerInBounds =
- mesh_.bounds().contains
- (
- globalPosition
- );
-
- label cell = -1;
- label tetFace = -1;
- label tetPt = -1;
-
- mesh_.findCellFacePt
- (
- globalPosition,
- cell,
- tetFace,
- tetPt
- );
-
- if (findIndex(zone, cell) != -1)
- {
- createMolecule
- (
- globalPosition,
- cell,
- tetFace,
- tetPt,
- id,
- tethered,
- temperature,
- bulkVelocity
- );
- }
- }
- }
- }
- }
-
- for
- (
- unitCellLatticePosition.z() = -(n-1);
- unitCellLatticePosition.z() <= (n-1);
- unitCellLatticePosition.z()++
- )
- {
- for (label iR = 0; iR <= 2*n -1; iR++)
- {
- unitCellLatticePosition.x() = n;
-
- unitCellLatticePosition.y() = -n + (iR + 1);
-
- for (label iK = 0; iK < 4; iK++)
- {
- forAll(latticePositions, p)
- {
- label id = findIndex
- (
- pot_.idList(),
- latticeIds[p]
- );
-
- const vector& latticePosition =
- vector
- (
- unitCellLatticePosition.x(),
- unitCellLatticePosition.y(),
- unitCellLatticePosition.z()
- )
- + latticePositions[p];
-
- point globalPosition =
- anchor
- + (
- R
- & (
- latticeCellShape
- & latticePosition
- )
- );
-
- partOfLayerInBounds =
- mesh_.bounds().contains
- (
- globalPosition
- );
-
- label cell = -1;
- label tetFace = -1;
- label tetPt = -1;
-
- mesh_.findCellFacePt
- (
- globalPosition,
- cell,
- tetFace,
- tetPt
- );
-
- if (findIndex(zone, cell) != -1)
- {
- createMolecule
- (
- globalPosition,
- cell,
- tetFace,
- tetPt,
- id,
- tethered,
- temperature,
- bulkVelocity
- );
- }
- }
-
- unitCellLatticePosition =
- labelVector
- (
- - unitCellLatticePosition.y(),
- unitCellLatticePosition.x(),
- unitCellLatticePosition.z()
- );
- }
- }
- }
- }
-
- if
- (
- totalZoneMols == 0
- && !partOfLayerInBounds
- )
- {
- WarningIn
- (
- "Foam::MoleculeCloud<MoleculeType>::"
- "initialiseMolecules()"
- )
- << "A whole layer of unit cells was placed "
- << "outside the bounds of the mesh, but no "
- << "molecules have been placed in zone '"
- << zone.name()
- << "'. This is likely to be because the zone "
- << "has few cells ("
- << zone.size()
- << " in this case) and no lattice position "
- << "fell inside them. "
- << "Aborting filling this zone."
- << endl;
-
- break;
- }
-
- molsPlacedThisIteration =
- this->size() - sizeBeforeIteration;
-
- totalZoneMols += molsPlacedThisIteration;
-
- n++;
- }
- }
- }
- }
-}
-
-
-template<class MoleculeType>
-void Foam::MoleculeCloud<MoleculeType>::createMolecule
-(
- const point& position,
- label cell,
- label tetFace,
- label tetPt,
- label id,
- bool tethered,
- scalar temperature,
- const vector& bulkVelocity
-)
-{
- if (cell == -1)
- {
- mesh_.findCellFacePt(position, cell, tetFace, tetPt);
- }
-
- if (cell == -1)
- {
- FatalErrorIn("Foam::MoleculeCloud<MoleculeType>::createMolecule")
- << "Position specified does not correspond to a mesh cell." << nl
- << abort(FatalError);
- }
-
- point specialPosition(vector::zero);
-
- label special = 0;
-
- if (tethered)
- {
- specialPosition = position;
-
- special = MoleculeType::SPECIAL_TETHERED;
- }
-
- const typename MoleculeType::constantProperties& cP(constProps(id));
-
- typename MoleculeType::trackingData td
- (
- *this,
- MoleculeType::trackingData::tpAccess
- );
-
- addParticle
- (
- new MoleculeType
- (
- mesh_,
- position,
- cell,
- tetFace,
- tetPt,
- temperature,
- bulkVelocity,
- specialPosition,
- cP,
- td,
- special,
- id
- )
- );
-}
-
-
-template<class MoleculeType>
-Foam::label Foam::MoleculeCloud<MoleculeType>::nSites() const
-{
- label n = 0;
-
- forAllConstIter(typename MoleculeCloud<MoleculeType>, *this, mol)
- {
- n += constProps(mol().id()).nSites();
- }
-
- return n;
-}
-
-
-// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
-
-template<class MoleculeType>
-Foam::MoleculeCloud<MoleculeType>::MoleculeCloud
-(
- const word& cloudName,
- const polyMesh& mesh,
- const potential& pot,
- bool readFields
-)
-:
- Cloud<MoleculeType>(mesh, cloudName, false),
- moleculeCloud(),
- mesh_(mesh),
- pot_(pot),
- cellOccupancy_(mesh_.nCells()),
- il_(mesh_, pot_.pairPotentials().rCutMax(), false),
- constPropList_(),
- rndGen_(label(971501) + 1526*Pstream::myProcNo())
-{
- if (readFields)
- {
- MoleculeType::readFields(*this);
- }
-
- buildConstProps();
-
- setSiteSizesAndPositions();
-
- removeHighEnergyOverlaps();
-
- calculateForce();
-}
-
-
-template<class MoleculeType>
-Foam::MoleculeCloud<MoleculeType>::MoleculeCloud
-(
- const word& cloudName,
- const polyMesh& mesh,
- const potential& pot,
- const dictionary& mdInitialiseDict,
- bool readFields
-)
-:
- Cloud<MoleculeType>(mesh, cloudName, false),
- moleculeCloud(),
- mesh_(mesh),
- pot_(pot),
- il_(mesh_),
- constPropList_(),
- rndGen_(label(971501) + 1526*Pstream::myProcNo())
-{
- if (readFields)
- {
- MoleculeType::readFields(*this);
- }
-
- this->clear();
-
- buildConstProps();
-
- initialiseMolecules(mdInitialiseDict);
-}
-
-
-// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
-
-template<class MoleculeType>
-void Foam::MoleculeCloud<MoleculeType>::evolve()
-{
- typename MoleculeType::trackingData td0
- (
- *this,
- MoleculeType::trackingData::tpFirstVelocityHalfStep
- );
- Cloud<MoleculeType>::move(td0, mesh_.time().deltaTValue());
-
- typename MoleculeType::trackingData td1
- (
- *this,
- MoleculeType::trackingData::tpLinearTrack
- );
- Cloud<MoleculeType>::move(td1, mesh_.time().deltaTValue());
-
- typename MoleculeType::trackingData td2
- (
- *this,
- MoleculeType::trackingData::tpRotationalTrack
- );
- Cloud<MoleculeType>::move(td2, mesh_.time().deltaTValue());
-
- calculateForce();
-
- typename MoleculeType::trackingData td3
- (
- *this,
- MoleculeType::trackingData::tpSecondVelocityHalfStep
- );
- Cloud<MoleculeType>::move(td3, mesh_.time().deltaTValue());
-
- info();
-}
-
-
-template<class MoleculeType>
-void Foam::MoleculeCloud<MoleculeType>::calculateForce()
-{
- buildCellOccupancy();
-
- // Set accumulated quantities to zero
- forAllIter(typename MoleculeCloud<MoleculeType>, *this, mol)
- {
- mol().siteForces() = vector::zero;
-
- mol().potentialEnergy() = 0.0;
-
- mol().rf() = tensor::zero;
- }
-
- calculatePairForce();
-
- calculateTetherForce();
-
- calculateExternalForce();
-}
-
-
-template<class MoleculeType>
-void Foam::MoleculeCloud<MoleculeType>::info()
-{
- // Calculates and prints the mean momentum and energy in the system
- // and the number of molecules.
-
- typename MoleculeType::trackingData td
- (
- *this,
- MoleculeType::trackingData::tpAccess
- );
-
- MoleculeType::info(td);
-}
-
-
-template<class MoleculeType>
-void Foam::MoleculeCloud<MoleculeType>::writeXYZ(const fileName& fName) const
-{
- OFstream os(fName);
-
- os << nSites() << nl
- << "MoleculeCloud<MoleculeType> site positions in angstroms" << nl;
-
- forAllConstIter(typename MoleculeCloud<MoleculeType>, *this, mol)
- {
- const typename MoleculeType::constantProperties& cP
- (
- constProps(mol().id())
- );
-
- forAll(mol().sitePositions(), i)
- {
- const point& sP = mol().sitePositions()[i];
-
- os << pot_.siteIdList()[cP.sites()[i].siteId()]
- << ' ' << sP.x()*1e10
- << ' ' << sP.y()*1e10
- << ' ' << sP.z()*1e10
- << nl;
- }
- }
-}
-
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/clouds/Templates/MoleculeCloud/MoleculeCloud.H b/src/lagrangian/molecularDynamics/molecule/clouds/Templates/MoleculeCloud/MoleculeCloud.H
deleted file mode 100644
index 38c912a17da..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/clouds/Templates/MoleculeCloud/MoleculeCloud.H
+++ /dev/null
@@ -1,252 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2008-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-Class
- Foam::MoleculeCloud
-
-Description
-
-SourceFiles
- MoleculeCloudI.H
- MoleculeCloud.C
-
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef MoleculeCloud_H
-#define MoleculeCloud_H
-
-#include "Cloud.H"
-#include "moleculeCloud.H"
-#include "IOdictionary.H"
-#include "potential.H"
-#include "InteractionLists.H"
-#include "labelVector.H"
-#include "Random.H"
-#include "fileName.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
-/*---------------------------------------------------------------------------*\
- Class MoleculeCloud Declaration
-\*---------------------------------------------------------------------------*/
-
-template<class MoleculeType>
-class MoleculeCloud
-:
- public Cloud<MoleculeType>,
- public moleculeCloud
-{
-
-private:
-
- // Private data
-
- //-
- const polyMesh& mesh_;
-
- //-
- const potential& pot_;
-
- //-
- List<DynamicList<MoleculeType*> > cellOccupancy_;
-
- //-
- InteractionLists<MoleculeType> il_;
-
- //-
- List<typename MoleculeType::constantProperties> constPropList_;
-
- //-
- Random rndGen_;
-
-
- // Private Member Functions
-
- //-
- void buildConstProps();
-
- //-
- void setSiteSizesAndPositions();
-
- //- Determine which molecules are in which cells
- void buildCellOccupancy();
-
- //-
- void calculatePairForce();
-
- //-
- inline void evaluatePair
- (
- MoleculeType& molI,
- MoleculeType& molJ
- );
-
- //-
- inline bool evaluatePotentialLimit
- (
- MoleculeType& molI,
- MoleculeType& molJ
- ) const;
-
- //-
- void calculateTetherForce();
-
- //-
- void calculateExternalForce();
-
- //-
- void removeHighEnergyOverlaps();
-
- //-
- void initialiseMolecules(const dictionary& mdInitialiseDict);
-
- //-
- void createMolecule
- (
- const point& position,
- label cell,
- label tetFace,
- label tetPt,
- label id,
- bool tethered,
- scalar temperature,
- const vector& bulkVelocity
- );
-
- //-
- label nSites() const;
-
- //- Disallow default bitwise copy construct
- MoleculeCloud(const MoleculeCloud&);
-
- //- Disallow default bitwise assignment
- void operator=(const MoleculeCloud&);
-
-
-public:
-
- // Constructors
-
- //- Construct given mesh and potential references
- MoleculeCloud
- (
- const word& cloudName,
- const polyMesh& mesh,
- const potential& pot,
- bool readFields = true
- );
-
- //- Construct given mesh, potential and mdInitialiseDict
- MoleculeCloud
- (
- const word& cloudName,
- const polyMesh& mesh,
- const potential& pot,
- const dictionary& mdInitialiseDict,
- bool readFields = true
- );
-
-
- // Member Functions
-
- //- Evolve the molecules (move, calculate forces, control state etc)
- void evolve();
-
- //-
- void calculateForce();
-
- //- Print cloud information
- void info();
-
-
- // Access
-
- //-
- inline const polyMesh& mesh() const;
-
- //-
- inline const potential& pot() const;
-
- //-
- inline const List<DynamicList<MoleculeType*> >&
- cellOccupancy() const;
-
- //-
- inline const InteractionLists<MoleculeType>& il() const;
-
- //-
- inline const List<typename MoleculeType::constantProperties>
- constProps() const;
-
- //-
- inline const typename MoleculeType::constantProperties&
- constProps(label id) const;
-
- //-
- inline Random& rndGen();
-
- //-
- inline vector equipartitionLinearVelocity
- (
- scalar temperature,
- scalar mass
- );
-
- //-
- inline vector equipartitionAngularMomentum
- (
- scalar temperature,
- const typename MoleculeType::constantProperties& cP
- );
-
-
- // Member Operators
-
- //- Write molecule sites in XYZ format
- void writeXYZ(const fileName& fName) const;
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#include "MoleculeCloudI.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#ifdef NoRepository
-# include "MoleculeCloud.C"
-#endif
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/clouds/Templates/MoleculeCloud/MoleculeCloudI.H b/src/lagrangian/molecularDynamics/molecule/clouds/Templates/MoleculeCloud/MoleculeCloudI.H
deleted file mode 100644
index 42c506024f7..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/clouds/Templates/MoleculeCloud/MoleculeCloudI.H
+++ /dev/null
@@ -1,416 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2008-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "constants.H"
-
-using namespace Foam::constant;
-
-// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
-
-template<class MoleculeType>
-inline void Foam::MoleculeCloud<MoleculeType>::evaluatePair
-(
- MoleculeType& molI,
- MoleculeType& molJ
-)
-{
- const pairPotentialList& pairPot = pot_.pairPotentials();
-
- const pairPotential& electrostatic = pairPot.electrostatic();
-
- label idI = molI.id();
-
- label idJ = molJ.id();
-
- const typename MoleculeType::constantProperties& constPropI
- (
- constProps(idI)
- );
-
- const typename MoleculeType::constantProperties& constPropJ
- (
- constProps(idJ)
- );
-
- forAll(constPropI.pairPotSites(), pI)
- {
- label sI = constPropI.pairPotSites()[pI];
-
- label idsI = constPropI.sites()[sI].siteId();
-
- forAll(constPropJ.pairPotSites(), pJ)
- {
- label sJ = constPropJ.pairPotSites()[pJ];
-
- label idsJ = constPropJ.sites()[sJ].siteId();
-
- vector rsIsJ =
- molI.sitePositions()[sI] - molJ.sitePositions()[sJ];
-
- scalar rsIsJMagSq = magSqr(rsIsJ);
-
- if (pairPot.rCutSqr(idsI, idsJ, rsIsJMagSq))
- {
- scalar rsIsJMag = mag(rsIsJ);
-
- vector fsIsJ =
- (rsIsJ/rsIsJMag)
- *pairPot.force(idsI, idsJ, rsIsJMag);
-
- molI.siteForces()[sI] += fsIsJ;
-
- molJ.siteForces()[sJ] += -fsIsJ;
-
- scalar potentialEnergy
- (
- pairPot.energy(idsI, idsJ, rsIsJMag)
- );
-
- molI.potentialEnergy() += 0.5*potentialEnergy;
-
- molJ.potentialEnergy() += 0.5*potentialEnergy;
-
- vector rIJ = molI.position() - molJ.position();
-
- tensor virialContribution =
- (rsIsJ*fsIsJ)*(rsIsJ & rIJ)/rsIsJMagSq;
-
- molI.rf() += virialContribution;
-
- molJ.rf() += virialContribution;
- }
- }
- }
-
-
- forAll(constPropI.electrostaticSites(), pI)
- {
- label sI = constPropI.electrostaticSites()[pI];
-
- forAll(constPropJ.electrostaticSites(), pJ)
- {
- label sJ = constPropJ.electrostaticSites()[pJ];
-
- vector rsIsJ =
- molI.sitePositions()[sI] - molJ.sitePositions()[sJ];
-
- scalar rsIsJMagSq = magSqr(rsIsJ);
-
- if (rsIsJMagSq <= electrostatic.rCutSqr())
- {
- scalar rsIsJMag = mag(rsIsJ);
-
- scalar chargeI = constPropI.sites()[sI].siteCharge();
-
- scalar chargeJ = constPropJ.sites()[sJ].siteCharge();
-
- vector fsIsJ =
- (rsIsJ/rsIsJMag)
- *chargeI*chargeJ*electrostatic.force(rsIsJMag);
-
- molI.siteForces()[sI] += fsIsJ;
-
- molJ.siteForces()[sJ] += -fsIsJ;
-
- scalar potentialEnergy =
- chargeI*chargeJ
- *electrostatic.energy(rsIsJMag);
-
- molI.potentialEnergy() += 0.5*potentialEnergy;
-
- molJ.potentialEnergy() += 0.5*potentialEnergy;
-
- vector rIJ = molI.position() - molJ.position();
-
- tensor virialContribution =
- (rsIsJ*fsIsJ)*(rsIsJ & rIJ)/rsIsJMagSq;
-
- molI.rf() += virialContribution;
-
- molJ.rf() += virialContribution;
- }
- }
- }
-}
-
-
-template<class MoleculeType>
-inline bool Foam::MoleculeCloud<MoleculeType>::evaluatePotentialLimit
-(
- MoleculeType& molI,
- MoleculeType& molJ
-) const
-{
- const pairPotentialList& pairPot = pot_.pairPotentials();
-
- const pairPotential& electrostatic = pairPot.electrostatic();
-
- label idI = molI.id();
-
- label idJ = molJ.id();
-
- const typename MoleculeType::constantProperties& constPropI
- (
- constProps(idI)
- );
-
- const typename MoleculeType::constantProperties& constPropJ
- (
- constProps(idJ)
- );
-
- forAll(constPropI.pairPotSites(), pI)
- {
- label sI = constPropI.pairPotSites()[pI];
-
- label idsI = constPropI.sites()[sI].siteId();
-
- forAll(constPropJ.pairPotSites(), pJ)
- {
- label sJ = constPropJ.pairPotSites()[pJ];
-
- label idsJ = constPropJ.sites()[sJ].siteId();
-
- vector rsIsJ =
- molI.sitePositions()[sI] - molJ.sitePositions()[sJ];
-
- scalar rsIsJMagSq = magSqr(rsIsJ);
-
- if (pairPot.rCutSqr(idsI, idsJ, rsIsJMagSq))
- {
- scalar rsIsJMag = mag(rsIsJ);
-
- // Guard against pairPot.energy being evaluated
- // if rIJMag < SMALL. A floating point exception will
- // happen otherwise.
-
- if (rsIsJMag < SMALL)
- {
- WarningIn
- (
- "MoleculeCloud<MoleculeType>::"
- "removeHighEnergyOverlaps()"
- )
- << "Molecule site pair closer than "
- << SMALL
- << ": mag separation = " << rsIsJMag
- << ". These may have been placed on top of each"
- << " other by a rounding error in mdInitialise in"
- << " parallel or a block filled with moleculess"
- << " twice. Removing one of the molecules."
- << endl;
-
- return true;
- }
-
- // Guard against pairPot.energy being evaluated if rIJMag <
- // rMin. A tabulation lookup error will occur otherwise.
-
- if (rsIsJMag < pairPot.rMin(idsI, idsJ))
- {
- return true;
- }
-
- if
- (
- mag(pairPot.energy(idsI, idsJ, rsIsJMag))
- > pot_.potentialEnergyLimit()
- )
- {
- return true;
- };
- }
- }
- }
-
- forAll(constPropI.electrostaticSites(), pI)
- {
- label sI = constPropI.electrostaticSites()[pI];
-
- forAll(constPropJ.electrostaticSites(), pJ)
- {
- label sJ = constPropJ.electrostaticSites()[pJ];
-
- vector rsIsJ =
- molI.sitePositions()[sI] - molJ.sitePositions()[sJ];
-
- scalar rsIsJMagSq = magSqr(rsIsJ);
-
- if (pairPot.rCutMaxSqr(rsIsJMagSq))
- {
- scalar rsIsJMag = mag(rsIsJ);
-
- // Guard against pairPot.energy being evaluated
- // if rIJMag < SMALL. A floating point exception will
- // happen otherwise.
-
- if (rsIsJMag < SMALL)
- {
- WarningIn
- (
- "MoleculeCloud<MoleculeType>::"
- "removeHighEnergyOverlaps()"
- )
- << "Molecule site pair closer than "
- << SMALL
- << ": mag separation = " << rsIsJMag
- << ". These may have been placed on top of each"
- << " other by a rounding error in mdInitialise in"
- << " parallel or a block filled with molecules"
- << " twice. Removing one of the molecules."
- << endl;
-
- return true;
- }
-
- if (rsIsJMag < electrostatic.rMin())
- {
- return true;
- }
-
- scalar chargeI = constPropI.sites()[sI].siteCharge();
-
- scalar chargeJ = constPropJ.sites()[sJ].siteCharge();
-
- if
- (
- mag(chargeI*chargeJ*electrostatic.energy(rsIsJMag))
- > pot_.potentialEnergyLimit()
- )
- {
- return true;
- };
- }
- }
- }
-
- return false;
-}
-
-
-// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
-
-template<class MoleculeType>
-inline const Foam::polyMesh& Foam::MoleculeCloud<MoleculeType>::mesh() const
-{
- return mesh_;
-}
-
-
-template<class MoleculeType>
-inline const Foam::potential& Foam::MoleculeCloud<MoleculeType>::pot() const
-{
- return pot_;
-}
-
-
-template<class MoleculeType>
-inline const Foam::List<Foam::DynamicList<MoleculeType*> >&
-Foam::MoleculeCloud<MoleculeType>::cellOccupancy() const
-{
- return cellOccupancy_;
-}
-
-
-template<class MoleculeType>
-inline const Foam::InteractionLists<MoleculeType>&
-Foam::MoleculeCloud<MoleculeType>::il() const
-{
- return il_;
-}
-
-
-template<class MoleculeType>
-inline const Foam::List<typename MoleculeType::constantProperties>
-Foam::MoleculeCloud<MoleculeType>::constProps() const
-{
- return constPropList_;
-}
-
-
-template<class MoleculeType>
-inline const typename MoleculeType::constantProperties&
-Foam::MoleculeCloud<MoleculeType>::constProps(label id) const
-{
- return constPropList_[id];
-}
-
-
-template<class MoleculeType>
-inline Foam::Random& Foam::MoleculeCloud<MoleculeType>::rndGen()
-{
- return rndGen_;
-}
-
-
-template<class MoleculeType>
-inline Foam::vector
-Foam::MoleculeCloud<MoleculeType>::equipartitionLinearVelocity
-(
- scalar temperature,
- scalar mass
-)
-{
- return sqrt(physicoChemical::k.value()*temperature/mass)*vector
- (
- rndGen_.GaussNormal(),
- rndGen_.GaussNormal(),
- rndGen_.GaussNormal()
- );
-}
-
-
-template<class MoleculeType>
-inline Foam::vector
-Foam::MoleculeCloud<MoleculeType>::equipartitionAngularMomentum
-(
- scalar temperature,
- const typename MoleculeType::constantProperties& cP
-)
-{
- scalar sqrtKbT = sqrt(physicoChemical::k.value()*temperature);
-
- if (cP.linearMolecule())
- {
- return sqrtKbT*vector
- (
- 0.0,
- sqrt(cP.momentOfInertia().yy())*rndGen_.GaussNormal(),
- sqrt(cP.momentOfInertia().zz())*rndGen_.GaussNormal()
- );
- }
- else
- {
- return sqrtKbT*vector
- (
- sqrt(cP.momentOfInertia().xx())*rndGen_.GaussNormal(),
- sqrt(cP.momentOfInertia().yy())*rndGen_.GaussNormal(),
- sqrt(cP.momentOfInertia().zz())*rndGen_.GaussNormal()
- );
- }
-}
-
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/clouds/baseClasses/moleculeCloud/moleculeCloud.C b/src/lagrangian/molecularDynamics/molecule/clouds/baseClasses/moleculeCloud/moleculeCloud.C
deleted file mode 100644
index 4fda04df7d5..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/clouds/baseClasses/moleculeCloud/moleculeCloud.C
+++ /dev/null
@@ -1,48 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2011-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "moleculeCloud.H"
-
-// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
-
-namespace Foam
-{
- defineTypeNameAndDebug(moleculeCloud, 0);
-}
-
-
-// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
-
-Foam::moleculeCloud::moleculeCloud()
-{}
-
-
-// * * * * * * * * * * * * * * * * Destructors * * * * * * * * * * * * * * //
-
-Foam::moleculeCloud::~moleculeCloud()
-{}
-
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/clouds/baseClasses/moleculeCloud/moleculeCloud.H b/src/lagrangian/molecularDynamics/molecule/clouds/baseClasses/moleculeCloud/moleculeCloud.H
deleted file mode 100644
index c6ba2096d36..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/clouds/baseClasses/moleculeCloud/moleculeCloud.H
+++ /dev/null
@@ -1,83 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2011-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-Class
- Foam::moleculeCloud
-
-Description
- Virtual abstract base class for templated moleculeCloud
-
-SourceFiles
- moleculeCloud.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef moleculeCloud_H
-#define moleculeCloud_H
-
-#include "volFields.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
-/*---------------------------------------------------------------------------*\
- Class moleculeCloud Declaration
-\*---------------------------------------------------------------------------*/
-
-class moleculeCloud
-{
- // Private Member Functions
-
- //- Disallow default bitwise copy construct
- moleculeCloud(const moleculeCloud&);
-
- //- Disallow default bitwise assignment
- void operator=(const moleculeCloud&);
-
-
-public:
-
- //- Runtime type information
- TypeName("moleculeCloud");
-
- // Constructors
-
- //- Null constructor
- moleculeCloud();
-
- //- Destructor
- virtual ~moleculeCloud();
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/clouds/derived/monoatomicCloud/monoatomicCloud.H b/src/lagrangian/molecularDynamics/molecule/clouds/derived/monoatomicCloud/monoatomicCloud.H
deleted file mode 100644
index 727bfe226e8..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/clouds/derived/monoatomicCloud/monoatomicCloud.H
+++ /dev/null
@@ -1,52 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2011-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-Class
- Foam::monoatomicCloud
-
-Description
- Cloud class to simulate monoatomic molecules
-
-SourceFiles
- monoatomicCloud.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef monoatomicCloud_H
-#define monoatomicCloud_H
-
-#include "MoleculeCloud.H"
-#include "monoatomic.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
- typedef MoleculeCloud<monoatomic> monoatomicCloud;
-}
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/clouds/derived/polyatomicCloud/polyatomicCloud.H b/src/lagrangian/molecularDynamics/molecule/clouds/derived/polyatomicCloud/polyatomicCloud.H
deleted file mode 100644
index bb902fe3a4d..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/clouds/derived/polyatomicCloud/polyatomicCloud.H
+++ /dev/null
@@ -1,52 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2011-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-Class
- Foam::polyatomicCloud
-
-Description
- Cloud class to simulate polyatomic molecules
-
-SourceFiles
- polyatomicCloud.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef polyatomicCloud_H
-#define polyatomicCloud_H
-
-#include "MoleculeCloud.H"
-#include "polyatomic.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
- typedef MoleculeCloud<polyatomic> polyatomicCloud;
-}
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/controllers/basic/controllers/controllers.C b/src/lagrangian/molecularDynamics/molecule/controllers/basic/controllers/controllers.C
deleted file mode 100644
index f0bad336ee1..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/controllers/basic/controllers/controllers.C
+++ /dev/null
@@ -1,573 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2011-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "controllers.H"
-
-// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
-
-Foam::controllers::controllers
-(
- const polyMesh& mesh
-)
-:
- time_(mesh.time()),
- controllersDict_
- (
- IOobject
- (
- "controllersDict",
- time_.system(),
- mesh,
- IOobject::NO_READ,
- IOobject::NO_WRITE
- )
- ),
- stateControllersList_(),
- sCNames_(),
- sCIds_(),
- sCFixedPathNames_(),
- stateControllers_(),
- fluxControllersList_(),
- fCNames_(),
- fCIds_(),
- fCFixedPathNames_(),
- fluxControllers_()
-{}
-
-
-Foam::controllers::controllers
-(
- const polyMesh& mesh,
- polyatomicCloud& cloud
-)
-:
- time_(mesh.time()),
- controllersDict_
- (
- IOobject
- (
- "controllersDict",
- time_.system(),
- mesh,
- IOobject::MUST_READ,
- IOobject::NO_WRITE
- )
- ),
- stateControllersList_(controllersDict_.lookup("stateControllers")),
- sCNames_(stateControllersList_.size()),
- sCIds_(stateControllersList_.size()),
- sCFixedPathNames_(stateControllersList_.size()),
- stateControllers_(stateControllersList_.size()),
- fluxControllersList_(controllersDict_.lookup("fluxControllers")),
- fCNames_(fluxControllersList_.size()),
- fCIds_(fluxControllersList_.size()),
- fCFixedPathNames_(fluxControllersList_.size()),
- fluxControllers_(fluxControllersList_.size())
-{
-
- Info << nl << "Creating controllers" << nl << endl;
-
- // state controllers
-
- if (!stateControllers_.empty())
- {
- forAll(stateControllers_, sC)
- {
- const entry& controllersI = stateControllersList_[sC];
-
- const dictionary& controllersIDict = controllersI.dict();
-
- stateControllers_[sC] = autoPtr<stateController>
- (
- stateController::New(time_, cloud, controllersIDict)
- );
-
- sCNames_[sC] = stateControllers_[sC]->type();
-
- sCIds_[sC] = sC;
- }
- }
-
- //- flux controllers
-
- if (!fluxControllers_.empty())
- {
- forAll(fluxControllers_, fC)
- {
- const entry& controllersI = fluxControllersList_[fC];
-
- const dictionary& controllersIDict = controllersI.dict();
-
- fluxControllers_[fC] = autoPtr<fluxController>
- (
- fluxController::New(time_, cloud, controllersIDict)
- );
-
- fCNames_[fC] = fluxControllers_[fC]->type();
- fCIds_[fC] = fC;
- }
- }
-
- // creating directories for state controllers
- if (!nStateControllers_.empty())
- {
- // case/controllers
- fileName controllersPath(time_.path()/"controllers");
-
- if (!isDir(controllersPath))
- {
- mkDir(controllersPath);
- }
-
- // case/controllers/<cloudName>
- fileName controllersPath(controllersPath/cloud.name());
-
- if (!isDir(controllersPath))
- {
- mkDir(controllersPath);
- }
-
- // case/controllers/<cloudName>/stateControllers
- fileName stateControllersPath(controllersPath/"stateControllers");
-
- if (!isDir(stateControllersPath))
- {
- mkDir(stateControllersPath);
- }
-
- forAll(stateControllers_, sC)
- {
- if (stateControllers_[sC]->writeInCase())
- {
- // case/controllers/<cloudName>/
- // stateControllers/<stateControllerModel>
- fileName stateControllerPath(stateControllersPath/sCNames_[sC]);
-
- if (!isDir(stateControllerPath))
- {
- mkDir(stateControllerPath);
- }
-
- const word& regionName = stateControllers_[sC]->regionName();
-
- // case/controllers/<cloudName>/
- // stateControllers/<stateControllerModel>/<cellZoneName>
- fileName zonePath(stateControllerPath/regionName);
-
- if (!isDir(zonePath))
- {
- mkDir(zonePath);
- }
-
- sCFixedPathNames_[sC] = zonePath;
- }
- }
- }
-
- // creating directories for flux controllers
- if (nFluxControllers_ > 0)
- {
- // case/controllers
- fileName controllersPath(time_.path()/"controllers");
-
- if ( !isDir(controllersPath) )
- {
- mkDir(controllersPath);
- }
-
- // case/controllers/<cloudName>
- fileName controllersPath(time_.path()/cloud.name());
-
- if ( !isDir(controllersPath) )
- {
- mkDir(controllersPath);
- }
-
- // case/controllers/<cloudName>/fluxControllers
- fileName fluxControllersPath(controllersPath/"fluxControllers");
-
- if (!isDir(fluxControllersPath))
- {
- mkDir(fluxControllersPath);
- }
-
- forAll(fluxControllers_, fC)
- {
- if (fluxControllers_[fC]->writeInCase())
- {
- // case/controllers/<cloudName>/
- // fluxControllers/<fluxControllerModel>
- fileName fluxControllerPath(fluxControllersPath/fCNames_[fC]);
-
- if (!isDir(fluxControllerPath))
- {
- mkDir(fluxControllerPath);
- }
-
- const word& regionName = fluxControllers_[fC]->regionName();
-
- // case/controllers/<cloudName>/
- // fluxControllers/<fluxControllerModel>/<faceZoneName>
- fileName zonePath(fluxControllerPath/regionName);
-
- if (!isDir(zonePath))
- {
- mkDir(zonePath);
- }
-
- fCFixedPathNames_[fC] = zonePath;
- }
- }
- }
-}
-
-
-// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
-
-controllers::~controllers()
-{}
-
-
-// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
-
-void controllers::initialConfig()
-{
- forAll(stateControllers_, sC)
- {
- stateControllers_[sC]->initialConfiguration();
- }
-
- forAll(fluxControllers_, fC)
- {
- fluxControllers_[fC]->initialConfiguration();
- }
-}
-
-
-void controllers::updateTimeInfo()
-{
- forAll(stateControllers_, sC)
- {
- stateControllers_[sC]->updateTime();
- }
-
- forAll(fluxControllers_, fC)
- {
- fluxControllers_[fC]->updateTime();
- }
-}
-
-
-void controllers::controlState()
-{
- forAll(stateControllers_, sC)
- {
- stateControllers_[sC]->controlMols();
- }
-}
-
-
-void controllers::controlVelocitiesI()
-{
- forAll(stateControllers_, sC)
- {
- stateControllers_[sC]->controlMolsBeg();
- }
-}
-
-void controllers::controlVelocitiesII()
-{
- forAll(stateControllers_, sC)
- {
- stateControllers_[sC]->controlMolsEnd();
- }
-}
-
-
-void controllers::controlPriorToForces()
-{
- forAll(stateControllers_, sC)
- {
- stateControllers_[sC]->controlBeforeForces();
- }
-}
-
-
-void controllers::calculateStateProps()
-{
- forAll(stateControllers_, sC)
- {
- stateControllers_[sC]->calculateProperties();
- }
-
- forAll(fluxControllers_, fC)
- {
- fluxControllers_[fC]->calculateProperties();
- }
-}
-
-
-void controllers::outputStateResults()
-{
- const Time& runTime = time_;
-
- if (runTime.outputTime())
- {
- // creating a set of directories in the current time directory
- {
- List<fileName> timePathNames(sCFixedPathNames_.size());
-
- if (nStateControllers_ > 0)
- {
- if (Pstream::master())
- {
- // case/<timeDir>/uniform
- fileName uniformTimePath
- (
- runTime.path()/runTime.timeName()/"uniform"
- );
-
- if (!isDir(uniformTimePath))
- {
- mkDir(uniformTimePath);
- }
-
- if (!stateControllers_.empty())
- {
- // case/<timeDir>/uniform/controllers
- fileName controllersTimePath
- (
- uniformTimePath/"controllers"
- );
-
- if (!isDir(controllersTimePath))
- {
- mkDir(controllersTimePath);
- }
-
- // case/<timeDir>/uniform/controllers/<cloudName>
- fileName cloudTimePath
- (
- controllersTimePath/cloud.name()
- );
-
- if (!isDir(cloudTimePath))
- {
- mkDir(cloudTimePath);
- }
-
- // case/<timeDir>/uniform/controllers/<cloudName>/
- fileName stateControllersTimePath
- (
- cloudTimePath/"stateControllers"
- );
-
- if (!isDir(stateControllersTimePath))
- {
- mkDir(stateControllersTimePath);
- }
-
- forAll(stateControllers_, sC)
- {
- if (stateControllers_[sC]->writeInTimeDir())
- {
- // case/<timeDir>/uniform/controllers/
- // <cloudName>/<stateControllerModel>
- fileName sCTimePath
- (
- stateControllersTimePath/sCNames_[sC]
- );
-
- if (!isDir(sCTimePath))
- {
- mkDir(sCTimePath);
- }
-
- // Creating directory for different zones but
- // of the same model
- const word& regionName =
- stateControllers_[sC]->regionName();
-
- // case/<timeDir>/uniform/controllers/
- // <cloudName>/<stateControllerModel>/
- // <cellZoneName>
- fileName zoneTimePath(sCTimePath/regionName);
-
- if (!isDir(zoneTimePath))
- {
- mkDir(zoneTimePath);
- }
-
- timePathNames[sC] = zoneTimePath;
- }
- }
- }
- }
- }
-
- // write out data
- forAll(stateControllers_, sC)
- {
- stateControllers_[sC]->output
- (
- sCFixedPathNames_[sC],
- timePathNames[sC]
- );
- }
- }
-
- {
- List<fileName> timePathNames(fCFixedPathNames_.size());
-
- if (nFluxControllers_ > 0)
- {
- if (Pstream::master())
- {
- // case/<timeDir>/uniform
- fileName uniformTimePath
- (
- runTime.path()/runTime.timeName()/"uniform"
- );
-
- if (!isDir(uniformTimePath))
- {
- mkDir(uniformTimePath);
- }
-
- if (!fluxControllers_.empty())
- {
- // case/<timeDir>/uniform/controllers
- fileName controllersTimePath
- (
- uniformTimePath/"controllers"
- );
-
- if (!isDir(controllersTimePath))
- {
- mkDir(controllersTimePath);
- }
-
- // case/<timeDir>/uniform/controllers/<cloudName>
- fileName cloudTimePath
- (
- controllersTimePath/cloud.name()
- );
-
- if (!isDir(cloudTimePath))
- {
- mkDir(cloudTimePath);
- }
-
- // case/<timeDir>/uniform/fluxControllers
- fileName controllersTimePath
- (
- cloudTimePath/"fluxControllers"
- );
-
- if (!isDir(controllersTimePath))
- {
- mkDir(controllersTimePath);
- }
-
- forAll(fluxControllers_, fC)
- {
- if (stateControllers_[fC]->writeInTimeDir())
- {
- // case/<timeDir>/uniform/controllers/
- // <cloudName>/<fluxControllerModel>
- fileName fCTimePath
- (
- controllersTimePath/fCNames_[fC]
- );
-
- if (!isDir(fCTimePath))
- {
- mkDir(fCTimePath);
- }
-
- const word& regionName =
- fluxControllers_[fC]->regionName();
-
- // case/<timeDir>/uniform/controllers/
- // <cloudName>/<fluxControllerModel>/
- // <faceZoneName>
- fileName zoneTimePath(fCTimePath/regionName);
-
- if (!isDir(zoneTimePath))
- {
- mkDir(zoneTimePath);
- }
-
- timePathNames[fC] = zoneTimePath;
- }
- }
- }
- }
- }
-
- // write out data
- forAll(fluxControllers_, fC)
- {
- fluxControllers_[fC]->output
- (
- fCFixedPathNames_[fC],
- timePathNames[fC]
- );
- }
- }
-
- // Re-read dictionaries for modified properties (run-time selection)
- {
- stateControllersList_.clear();
-
- stateControllersList_ = controllersDict_.lookup("stateControllers");
-
- forAll(stateControllers_, sC)
- {
- const entry& controllersI = stateControllersList_[sC];
- const dictionary& controllersIDict = controllersI.dict();
-
- stateControllers_[sC]->updateProperties(controllersIDict);
- }
- }
-
- {
- fluxControllersList_.clear();
-
- fluxControllersList_ = controllersDict_.lookup("fluxControllers");
-
- forAll(fluxControllers_, fC)
- {
- const entry& controllersI = fluxControllersList_[fC];
- const dictionary& controllersIDict = controllersI.dict();
-
- fluxControllers_[fC]->updateProperties(controllersIDict);
- }
- }
- }
-}
-
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/controllers/basic/controllers/controllers.H b/src/lagrangian/molecularDynamics/molecule/controllers/basic/controllers/controllers.H
deleted file mode 100644
index 22e9ca7cca1..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/controllers/basic/controllers/controllers.H
+++ /dev/null
@@ -1,163 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2008-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-Class
-
- controllers
-
-Description
-
- Stores all the information for the controllers models defined within
- the controllersDict, and selects & builds the models automatically.
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef controllers_H
-#define controllers_H
-
-#include "List.H"
-#include "IOdictionary.H"
-#include "autoPtr.H"
-#include "polyMesh.H"
-#include "timeData.H"
-#include "stateController.H"
-#include "fluxController.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
-/*---------------------------------------------------------------------------*\
- Class controllers Declaration
-\*---------------------------------------------------------------------------*/
-
-class controllers
-{
- // Private data
-
- Time& time_;
-
- //- The entire dictionary (containing multiple subDictionaries)
- IOdictionary controllersDict_;
-
- //- state controllers
- PtrList<entry> stateControllersList_;
- List<word> sCNames_;
- List<label> sCIds_;
- List<fileName> sCFixedPathNames_;
- List< autoPtr<stateController> > stateControllers_;
-
- //- flux controllers
- PtrList<entry> fluxControllersList_;
- List<word> fCNames_;
- List<label> fCIds_;
- List<fileName> fCFixedPathNames_;
- List< autoPtr<fluxController> > fluxControllers_;
-
-
-public:
-
- // Constructors
-
- //- Null Constructor
- controllers
- (
- const polyMesh& mesh
- );
-
- //- Constructor for with cloud
- controllers
- (
- const polyMesh& mesh,
- polyatomicCloud& cloud
- );
-
- //- Destructor
- ~controllers();
-
- // Member Functions
-
- //- Initial configuration call this function after the polyatomicCloud
- // is completely initialised
- void initialConfig();
-
- //- this function is to be called at the beginning of the MD time-step.
- // since we have placed a non-referenced time-data class in the
- // state-controller class.
- void updateTimeInfo();
-
- //- control molecular state -- call this after the intermolecular force
- // calulation
- void controlState();
-
- //-
- void controlVelocitiesI();
-
- //-
- void controlVelocitiesII();
-
- //-
- void controlPriorToForces();
-
- //- calculate properties -- call this at the end of the MD time-step.
- void calculateStateProps();
-
- //- output -- call this function at the end of the MD time-step
- void outputStateResults();
-
- // Access
-
- //-
- inline List< autoPtr<stateController> >& stateControllers();
-
- //-
- inline const List< autoPtr<stateController> >&
- stateControllers() const;
-
- //-
- inline List< autoPtr<fluxController> >& fluxControllers();
-
- //-
- inline const List< autoPtr<fluxController> >&
- fluxControllers() const;
-
- //-
- inline const List<word>& stateControllersNames() const;
-
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#include "controllersI.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/controllers/basic/controllers/controllersI.H b/src/lagrangian/molecularDynamics/molecule/controllers/basic/controllers/controllersI.H
deleted file mode 100644
index 66c10f07840..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/controllers/basic/controllers/controllersI.H
+++ /dev/null
@@ -1,64 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2011-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-\*---------------------------------------------------------------------------*/
-
-// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
-
-Foam::List<Foam::autoPtr<Foam::stateController> >&
-Foam::controllers::stateControllers()
-{
- return stateControllers_;
-}
-
-
-const Foam::List<Foam::autoPtr<Foam::stateController> >&
-Foam::controllers::stateControllers() const
-{
- return stateControllers_;
-}
-
-
-Foam::List<Foam::autoPtr<Foam::fluxController> >&
-Foam::controllers::fluxControllers()
-{
- return fluxControllers_;
-}
-
-
-const Foam::List< autoPtr<fluxController> >&
-Foam::controllers::fluxControllers() const
-{
- return fluxControllers_;
-}
-
-
-const Foam::List<Foam::word>& Foam::controllers::stateControllersNames() const
-{
- return sCNames_;
-}
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/controllers/basic/fluxController/fluxController.C b/src/lagrangian/molecularDynamics/molecule/controllers/basic/fluxController/fluxController.C
deleted file mode 100644
index 711b6617588..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/controllers/basic/fluxController/fluxController.C
+++ /dev/null
@@ -1,524 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2011-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "waterFluxController.H"
-#include "IFstream.H"
-#include "graph.H"
-#include "polyatomicCloud.H"
-
-namespace Foam
-{
-
-// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
-
-defineTypeNameAndDebug(waterFluxController, 0);
-
-defineRunTimeSelectionTable(waterFluxController, dictionary);
-
-// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
-
-// Construct from components
-waterFluxController::waterFluxController
-(
- Time& t,
- polyatomicCloud& cloud,
- const dictionary& dict
-)
-:
- mesh_(refCast<const fvMesh>(cloud.mesh())),
- cloud_(cloud),
- rndGen_(clock::getTime()),
- controllerDict_(dict.subDict("controllerProperties")),
- timeDict_(controllerDict_.subDict("timeProperties")),
- time_(t, timeDict_),
- regionName_(controllerDict_.lookup("zoneName")),
- regionId_(-1),
- zoneSurfaceArea_(0.0),
- internalFaces_(),
- processorFaces_(),
- control_(true),
- readStateFromFile_(true),
- singleValueController_(false),
- density_(0.0),
- velocity_(vector::zero),
- temperature_(0.0),
- pressure_(0.0),
- strainRate_(tensor::zero),
- tempGradient_(vector::zero),
- fieldController_(false),
- densities_(),
- velocities_(),
- temperatures_(),
- pressures_(),
- writeInTimeDir_(true),
- writeInCase_(true)
-{
- const faceZoneMesh& faceZones = mesh_.faceZones();
- regionId_ = faceZones.findZoneID(regionName_);
-
- if (regionId_ == -1)
- {
- FatalErrorIn("waterFluxController::waterFluxController()")
- << "Cannot find region (faceZone): " << regionName_ << nl << "in: "
- << t.time().system()/"controllersDict"
- << exit(FatalError);
- }
-
- control_ = Switch(controllerDict_.lookup("controlSwitch"));
- readStateFromFile_ = Switch(controllerDict_.lookup("readStateFromFile"));
-
- setFacesInfo();
-}
-
-
-// * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
-
-autoPtr<waterFluxController> waterFluxController::New
-(
- Time& t,
- polyatomicCloud& cloud,
- const dictionary& dict
-)
-{
- word waterFluxControllerName
- (
- dict.lookup("fluxControllerModel")
- );
-
- Info<< "Selecting fluxController "
- << waterFluxControllerName << endl;
-
- dictionaryConstructorTable::iterator cstrIter =
- dictionaryConstructorTablePtr_->find(waterFluxControllerName);
-
- if (cstrIter == dictionaryConstructorTablePtr_->end())
- {
- FatalError
- << "waterFluxController::New(const dictionary&) : " << endl
- << " unknown waterFluxController type "
- << waterFluxControllerName
- << ", constructor not in hash table" << endl << endl
- << " Valid injector types are :" << endl;
- Info<< dictionaryConstructorTablePtr_->toc() << abort(FatalError);
- }
-
- return autoPtr<waterFluxController>
- (
- cstrIter()(t, cloud, dict)
- );
-}
-
-
-// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
-
-waterFluxController::~waterFluxController()
-{}
-
-// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
-
-// void waterFluxController::updateTime()
-// {
-// time_++;
-//
-// const scalar& t = time_.time().timeOutputValue();
-//
-// if ((t - initialTime_) < timePeriod_)
-// {
-// time_.controlTimeInterval().endTime() = false;
-// // control_ = false;
-// }
-// else
-// {
-// // control_ = true;
-// }
-// }
-
-
-void waterFluxController::setFacesInfo()
-{
- const labelList& faces = controlZone();
-
- if (Pstream::parRun())
- {
- DynamicList<label> processorFaces(0);
-
- forAll(mesh_.boundaryMesh(), patchI)
- {
- const polyPatch& patch = mesh_.boundaryMesh()[patchI];
-
- if (isA<processorPolyPatch>(patch))
- {
- for (label p = 0; p < patch.size(); p++)
- {
- label patchFaceI = p + patch.start();
- label faceId = findIndex (faces, patchFaceI);
-
- if (faceId != -1)
- {
- processorFaces.append(patchFaceI);
- }
- }
- }
- }
-
- processorFaces.shrink();
-
- processorFaces_.setSize(processorFaces.size(), -1);
-
- forAll(processorFaces, f)
- {
- processorFaces_[f] = processorFaces[f];
- }
-
- label nInternalFaces = faces.size() - processorFaces.size();
- internalFaces_.setSize(nInternalFaces, -1);
-
- label counter = 0;
-
- forAll(faces, f)
- {
- const label& faceI = faces[f];
-
- if (findIndex(processorFaces, faceI) == -1)
- {
- internalFaces_[counter] = faceI;
- counter++;
- }
- }
-
-// Pout << "processorFaces: " << processorFaces_ << endl;
-// Pout << "internalFaces: " << internalFaces_ << endl;
-
- forAll(internalFaces_, f)
- {
- const label& faceI = internalFaces_[f];
- zoneSurfaceArea_ += mag(mesh_.faceAreas()[faceI]);
- }
-
- // faces on a zone located on a processor cut belong to both processors
- // (hence the 0.5)
-
- forAll(processorFaces_, f)
- {
- const label& faceI = processorFaces_[f];
- zoneSurfaceArea_ += 0.5*mag(mesh_.faceAreas()[faceI]);
- }
-
-
- if (Pstream::parRun())
- {
- for (int p = 0; p < Pstream::nProcs(); p++)
- {
- if (p != Pstream::myProcNo())
- {
- const int proc = p;
- {
- OPstream toNeighbour(Pstream::blocking, proc);
- toNeighbour << zoneSurfaceArea_;
- }
- }
- }
-
- //- receiving
- for (int p = 0; p < Pstream::nProcs(); p++)
- {
- if (p != Pstream::myProcNo())
- {
- scalar zoneSurfaceAreaProc;
-
- const int proc = p;
- {
- IPstream fromNeighbour(Pstream::blocking, proc);
- fromNeighbour >> zoneSurfaceAreaProc;
- }
-
- zoneSurfaceArea_ += zoneSurfaceAreaProc;
- }
- }
- }
- }
- else
- {
- forAll(faces, f)
- {
- const label& faceI = faces[f];
-
- zoneSurfaceArea_ += mag(mesh_.faceAreas()[faceI]);
- }
- }
-}
-
-
-
-void waterFluxController::updateTime()
-{
- time_++;
-
-// const scalar& t = time_.time().timeOutputValue();
-//
-// if ((t - initialTime_) < timePeriod_)
-// {
-// time_.controlTimeInterval().endTime() = false;
-// // control_ = false;
-// }
-// else
-// {
-// // control_ = true;
-// }
-}
-
-void waterFluxController::updateFluxControllerProperties
-(
- const dictionary& newDict
-)
-{
- controllerDict_ = newDict.subDict("controllerProperties");
-
- //- you can reset the controlling zone from here. This essentially
- // means that the coupling zone can infact move arbitrarily. To make
- // this happen we probably need to devise a technique for automatically
- // changing the cellZone else where, and then calling this function to
- // reset the controlling zone in which the controller operates in.
-
- if (controllerDict_.found("controlSwitch"))
- {
- control_ = Switch(controllerDict_.lookup("controlSwitch"));
- }
-
- if (controllerDict_.found("readStateFromFile"))
- {
- readStateFromFile_ = Switch
- (
- controllerDict_.lookup("readStateFromFile")
- );
- }
-}
-
-const labelList& waterFluxController::controlZone() const
-{
- return mesh_.faceZones()[regionId_];
-}
-
-label waterFluxController::isFaceOnControlZone(const label& faceI)
-{
- const label f = findIndex(controlZone(), faceI);
-
- return f;
-}
-
-const word& waterFluxController::regionName() const
-{
- return regionName_;
-}
-
-const scalar& waterFluxController::density() const
-{
- return density_;
-}
-
-scalar& waterFluxController::density()
-{
- return density_;
-}
-
-const vector& waterFluxController::velocity() const
-{
- return velocity_;
-}
-
-vector& waterFluxController::velocity()
-{
- return velocity_;
-}
-
-const scalar& waterFluxController::temperature() const
-{
- return temperature_;
-}
-
-scalar& waterFluxController::temperature()
-{
- return temperature_;
-}
-
-const scalar& waterFluxController::pressure() const
-{
- return pressure_;
-}
-
-scalar& waterFluxController::pressure()
-{
- return pressure_;
-}
-
-const tensor& waterFluxController::strainRate() const
-{
- return strainRate_;
-}
-
-tensor& waterFluxController::strainRate()
-{
- return strainRate_;
-}
-
-const vector& waterFluxController::tempGradient() const
-{
- return tempGradient_;
-}
-
-vector& waterFluxController::tempGradient()
-{
- return tempGradient_;
-}
-
-
-const scalarField& waterFluxController::densityField() const
-{
- return densities_;
-}
-
-scalarField& waterFluxController::densityField()
-{
- return densities_;
-}
-
-const vectorField& waterFluxController::velocityField() const
-{
- return velocities_;
-}
-vectorField& waterFluxController::velocityField()
-{
- return velocities_;
-}
-
-const scalarField& waterFluxController::temperatureField() const
-{
- return temperatures_;
-}
-
-scalarField& waterFluxController::temperatureField()
-{
- return temperatures_;
-}
-
-const scalarField& waterFluxController::pressureField() const
-{
- return pressures_;
-}
-
-scalarField& waterFluxController::pressureField()
-{
- return pressures_;
-}
-
-
-const bool& waterFluxController::singleValueController() const
-{
- return singleValueController_;
-}
-
-bool& waterFluxController::singleValueController()
-{
- return singleValueController_;
-}
-
-const bool& waterFluxController::fieldController() const
-{
- return fieldController_;
-}
-
-bool& waterFluxController::fieldController()
-{
- return fieldController_;
-}
-
-
-const bool& waterFluxController::writeInTimeDir() const
-{
- return writeInTimeDir_;
-}
-
-const bool& waterFluxController::writeInCase() const
-{
- return writeInCase_;
-}
-
-
-// const scalar waterFluxController::avReqDensity() const
-// {
-// scalar totalDensity = 0.0;
-//
-// forAll(densities_, c)
-// {
-// totalDensity += densities_[c];
-// }
-//
-// if (cells_.size() > 0)
-// {
-// totalDensity /= scalar(cells_.size());
-// }
-//
-// return totalDensity;
-// }
-//
-// const vector waterFluxController::avReqVelocity() const
-// {
-// vector totalVel = vector::zero;
-//
-// forAll(velocities_, c)
-// {
-// totalVel += velocities_[c];
-// }
-//
-// if (cells_.size() > 0)
-// {
-// totalVel /= scalar(cells_.size());
-// }
-//
-// return totalVel;
-// }
-//
-// const scalar waterFluxController::avReqTemperature() const
-// {
-// scalar totalTemp = 0.0;
-//
-// forAll(densities_, c)
-// {
-// totalTemp += temperatures_[c];
-// }
-//
-// if (cells_.size() > 0)
-// {
-// totalTemp /= scalar(cells_.size());
-// }
-//
-// return totalTemp;
-// }
-
-
-
-} // End namespace Foam
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/controllers/basic/fluxController/fluxController.H b/src/lagrangian/molecularDynamics/molecule/controllers/basic/fluxController/fluxController.H
deleted file mode 100644
index 9b8a3752772..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/controllers/basic/fluxController/fluxController.H
+++ /dev/null
@@ -1,272 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2008-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-Class
- waterFluxController
-
-Description
-
-SourceFiles
- waterFluxControllerI.H
- waterFluxController.C
- waterFluxControllerIO.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef waterFluxController_H
-#define waterFluxController_H
-
-#include "IOdictionary.H"
-#include "Time.H"
-#include "autoPtr.H"
-#include "runTimeSelectionTables.H"
-#include "vector.H"
-#include "volFields.H"
-#include "Random.H"
-#include "polyatomic.H"
-#include "timeData.H"
-#include "writeTimeData.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
-/*---------------------------------------------------------------------------*\
- Class waterFluxController Declaration
-\*---------------------------------------------------------------------------*/
-
-class waterFluxController
-{
-
-protected:
-
- // Protected data
-
-// Time& time_;
-
- const fvMesh& mesh_;
-
- polyatomicCloud& cloud_;
-
- Random rndGen_;
-
- //- subDictionary containing the properties
- dictionary controllerDict_;
-
-
- dictionary timeDict_;
-
- timeData time_;
-
- //- name of face zone
- word regionName_;
- label regionId_;
-// labelList faces_;
-
- scalar zoneSurfaceArea_;
-
- labelList internalFaces_;
- labelList processorFaces_;
-
- bool control_;
- bool readStateFromFile_;
-
- //- set all the properties below from model if required
-
- bool singleValueController_;
-
- // target values
- scalar density_;
- vector velocity_;
- scalar temperature_;
- scalar pressure_;
-
- tensor strainRate_;
- vector tempGradient_;
-
- bool fieldController_;
-
- //- targeted fields
- scalarField densities_;
- vectorField velocities_;
- scalarField temperatures_;
- scalarField pressures_;
-
- bool writeInTimeDir_;
- bool writeInCase_;
-
-
-
-
- // Private Member Functions
-
- void setFacesInfo();
-
-
-public:
-
- //- Runtime type information
- TypeName("waterFluxController");
-
- // Declare runtime constructor selection table
- declareRunTimeSelectionTable
- (
- autoPtr,
- waterFluxController,
- dictionary,
- (
- Time& t,
- polyatomicCloud& cloud,
- const dictionary& dict
- ),
- (t, cloud, dict)
- );
-
- // Constructors
-
- //- Construct from components
- waterFluxController
- (
- Time& t,
- polyatomicCloud& cloud,
- const dictionary& dict
- );
-
-
- // Selectors
-
- static autoPtr<waterFluxController> New
- (
- Time& t,
- polyatomicCloud& cloud,
- const dictionary& dict
- );
-
-
- // Destructor
-
- virtual ~waterFluxController();
-
-
- // Member Functions
-
- void updateTime();
-
- //- create an initial configuration
- virtual void initialConfiguration() = 0;
-
-
- //- calculate any required properties
- virtual void calculateProperties() = 0;
-
-
- //- control the polyatomic from the tracking function
- virtual void controlMol
- (
- polyatomic& mol,
- polyatomic::trackData& td
- ) = 0;
-
- //- output data
- virtual void output
- (
- const fileName& fixedPathName,
- const fileName& timePath
- ) = 0;
-
-
- //- E. update properties from a modified dictionary
- virtual void updateProperties(const dictionary&) = 0;
-
-
- void updateFluxControllerProperties(const dictionary&);
-
- // Access
-
- //- return the control zone cells
- const labelList& controlZone() const;
-
-
- label isFaceOnControlZone(const label& faceI);
-
- //- return the control zone name
- const word& regionName() const;
-
- //- return the targeted values
- const scalar& density() const;
- scalar& density();
-
- const vector& velocity() const;
- vector& velocity();
-
- const scalar& temperature() const;
- scalar& temperature();
-
- const scalar& pressure() const;
- scalar& pressure();
-
- const tensor& strainRate() const;
- tensor& strainRate();
-
- const vector& tempGradient() const;
- vector& tempGradient();
-
- //- return the targeted fields
- const scalarField& densityField() const;
- scalarField& densityField();
-
- const vectorField& velocityField() const;
- vectorField& velocityField();
-
- const scalarField& temperatureField() const;
- scalarField& temperatureField();
-
- const scalarField& pressureField() const;
- scalarField& pressureField();
-
-
- const bool& singleValueController() const;
- bool& singleValueController();
-
- const bool& fieldController() const;
- bool& fieldController();
-
- const bool& writeInTimeDir() const;
- const bool& writeInCase() const;
-
-// const scalar avReqDensity() const;
-// const vector avReqVelocity() const;
-// const scalar avReqTemperature() const;
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/controllers/basic/stateController/stateController.C b/src/lagrangian/molecularDynamics/molecule/controllers/basic/stateController/stateController.C
deleted file mode 100644
index b98aa126a3e..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/controllers/basic/stateController/stateController.C
+++ /dev/null
@@ -1,490 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2011-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "stateController.H"
-#include "IFstream.H"
-#include "polyatomicCloud.H"
-
-// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
-defineTypeNameAndDebug(stateController, 0);
-
-defineRunTimeSelectionTable(stateController, dictionary);
-
-}
-
-
-// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
-
-Foam::stateController::stateController
-(
- polyatomicCloud& cloud,
- const dictionary& dict
-)
-:
- mesh_(refCast<const fvMesh>(cloud.mesh())),
- cloud_(cloud),
- rndGen_(clock::getTime()),
- controllerDict_(dict.subDict("controllerProperties")),
- timeDict_(controllerDict_.subDict("timeProperties")),
- time_(mesh_.time(), timeDict_),
- timePeriod_(readScalar(timeDict_.lookup("initialTimePeriod"))), //temp
- initialTime_(time_.time().startTime().value()),
- regionName_(controllerDict_.lookup("zoneName")),
- regionId_(-1),
- control_(true),
- readStateFromFile_(true),
- singleValueController_(false),
- density_(0.0),
- velocity_(vector::zero),
- temperature_(0.0),
- pressure_(0.0),
- strainRate_(tensor::zero),
- tempGradient_(vector::zero),
- fieldController_(false),
- densities_(),
- velocities_(),
- temperatures_(),
- pressures_(),
- writeInTimeDir_(true),
- writeInCase_(true)
-{
- const cellZoneMesh& cellZones = mesh_.cellZones();
-
- regionId_ = cellZones.findZoneID(regionName_);
-
- if (regionId_ == -1)
- {
- FatalErrorIn("stateController::stateController()")
- << "Cannot find region: " << regionName_ << nl << "in: "
- << time_.time().system()/"controllersDict"
- << exit(FatalError);
- }
-
- control_ = Switch(controllerDict_.lookup("controlSwitch"));
-
- readStateFromFile_ = Switch(controllerDict_.lookup("readStateFromFile"));
-
- const scalar& avTimeInterval = time_.averageTimeInterval().deltaT();
-
- if ((timePeriod_ < avTimeInterval) && (timePeriod_ > 0.0))
- {
- timePeriod_ = avTimeInterval;
- }
-}
-
-
-// * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
-
-Foam::autoPtr<Foam::stateController> Foam::stateController::New
-(
- polyatomicCloud& cloud,
- const dictionary& dict
-)
-{
- word stateControllerName
- (
- dict.lookup("stateControllerModel")
- );
-
- Info<< "Selecting stateController "
- << stateControllerName << endl;
-
- dictionaryConstructorTable::iterator cstrIter =
- dictionaryConstructorTablePtr_->find(stateControllerName);
-
- if (cstrIter == dictionaryConstructorTablePtr_->end())
- {
- FatalError
- << "stateController::New(const dictionary&) : " << endl
- << " unknown stateController type "
- << stateControllerName
- << ", constructor not in hash table" << endl << endl
- << " Valid types are :" << endl;
- Info<< dictionaryConstructorTablePtr_->toc() << abort(FatalError);
- }
-
- return autoPtr<stateController>
- (
- cstrIter()(cloud, dict)
- );
-}
-
-
-// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
-
-Foam::stateController::~stateController()
-{}
-
-
-// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
-
-void Foam::stateController::updateTime()
-{
- time_++;
-
- const scalar& t = time_.time().timeOutputValue();
-
- if ((t - initialTime_) < timePeriod_)
- {
- time_.controlTimeInterval().endTime() = false;
- }
-}
-
-
-void Foam::stateController::updateStateControllerProperties
-(
- const dictionary& newDict
-)
-{
- controllerDict_ = newDict.subDict("controllerProperties");
-
- if (controllerDict_.found("controlSwitch"))
- {
- control_ = Switch(controllerDict_.lookup("controlSwitch"));
- }
-
- if (controllerDict_.found("readStateFromFile"))
- {
- readStateFromFile_ = Switch
- (
- controllerDict_.lookup("readStateFromFile")
- );
- }
-
- timeDict_ = controllerDict_.subDict("timeProperties");
-
- if (timeDict_.found("resetAtOutput"))
- {
- time_.resetFieldsAtOutput() = Switch(timeDict_.lookup("resetAtOutput"));
- }
-}
-
-
-const Foam::labelList& Foam::stateController::controlZone() const
-{
- return mesh_.cellZones()[regionId_];
-}
-
-const Foam::word& Foam::stateController::regionName() const
-{
- return regionName_;
-}
-
-
-Foam::scalar Foam::stateController::density() const
-{
- return density_;
-}
-
-
-Foam::scalar& Foam::stateController::density()
-{
- return density_;
-}
-
-
-const Foam::vector& Foam::stateController::velocity() const
-{
- return velocity_;
-}
-
-
-Foam::vector& Foam::stateController::velocity()
-{
- return velocity_;
-}
-
-
-Foam::scalar Foam::stateController::temperature() const
-{
- return temperature_;
-}
-
-
-Foam::scalar& Foam::stateController::temperature()
-{
- return temperature_;
-}
-
-const Foam::scalar& Foam::stateController::pressure() const
-{
- return pressure_;
-}
-
-
-Foam::scalar& Foam::stateController::pressure()
-{
- return pressure_;
-}
-
-
-const Foam::tensor& Foam::stateController::strainRate() const
-{
- return strainRate_;
-}
-
-
-Foam::tensor& Foam::stateController::strainRate()
-{
- return strainRate_;
-}
-
-
-const Foam::vector& Foam::stateController::tempGradient() const
-{
- return tempGradient_;
-}
-
-
-Foam::vector& Foam::stateController::tempGradient()
-{
- return tempGradient_;
-}
-
-
-const Foam::scalarField& Foam::stateController::densityField() const
-{
- return densities_;
-}
-
-Foam::scalarField& Foam::stateController::densityField()
-{
- return densities_;
-}
-
-
-const Foam::vectorField& Foam::stateController::velocityField() const
-{
- return velocities_;
-}
-
-
-Foam::vectorField& Foam::stateController::velocityField()
-{
- return velocities_;
-}
-
-
-const Foam::scalarField& Foam::stateController::temperatureField() const
-{
- return temperatures_;
-}
-
-
-Foam::scalarField& Foam::stateController::temperatureField()
-{
- return temperatures_;
-}
-
-
-const Foam::scalarField& Foam::stateController::pressureField() const
-{
- return pressures_;
-}
-
-
-Foam::scalarField& Foam::stateController::pressureField()
-{
- return pressures_;
-}
-
-
-bool Foam::stateController::singleValueController() const
-{
- return singleValueController_;
-}
-
-
-bool& Foam::stateController::singleValueController()
-{
- return singleValueController_;
-}
-
-
-bool Foam::stateController::fieldController() const
-{
- return fieldController_;
-}
-
-
-bool& Foam::stateController::fieldController()
-{
- return fieldController_;
-}
-
-
-bool Foam::stateController::writeInTimeDir() const
-{
- return writeInTimeDir_;
-}
-
-
-bool Foam::stateController::writeInCase() const
-{
- return writeInCase_;
-}
-
-
-Foam::scalar Foam::stateController::avReqDensity()
-{
- scalar totalDensity = 0.0;
-
- if (singleValueController_)
- {
- totalDensity = density_;
- }
- else if (fieldController_)
- {
- label controlCells = controlZone().size();
-
- forAll(densities_, c)
- {
- totalDensity += densities_[c];
- }
-
- if (Pstream::parRun())
- {
- reduce(totalDensity, sumOp<scalar>());
-
- reduce(controlCells, sumOp<label>());
- }
-
- if (controlCells > 0)
- {
- totalDensity /= scalar(controlCells);
- }
- }
-
- return totalDensity;
-}
-
-
-Foam::vector Foam::stateController::avReqVelocity()
-{
- vector totalVel = vector::zero;
-
- if (singleValueController_)
- {
- totalVel = velocity_;
- }
- else if (fieldController_)
- {
- label controlCells = controlZone().size();
-
- forAll(velocities_, c)
- {
- totalVel += velocities_[c];
- }
-
- if (Pstream::parRun())
- {
- reduce(totalVel, sumOp<vector>());
-
- reduce(controlCells, sumOp<label>());
- }
-
- if (controlCells > 0)
- {
- totalVel /= scalar(controlCells);
- }
- }
-
- return totalVel;
-}
-
-
-Foam::scalar Foam::stateController::avReqTemperature()
-{
- scalar totalTemp = 0.0;
-
- if (singleValueController_)
- {
- totalTemp = temperature_;
- }
- else if (fieldController_)
- {
- label controlCells = controlZone().size();
-
- forAll(temperatures_, c)
- {
- totalTemp += temperatures_[c];
- }
-
- if (Pstream::parRun())
- {
- reduce(totalTemp, sumOp<scalar>());
-
- reduce(controlCells, sumOp<label>());
- }
-
- if (controlCells > 0)
- {
- totalTemp /= scalar(controlCells);
- }
- }
-
- return totalTemp;
-}
-
-
-Foam::scalar Foam::stateController::avReqPressure()
-{
- scalar totalPressure = 0.0;
-
- if (singleValueController_)
- {
- totalPressure = pressure_;
- }
- else if (fieldController_)
- {
- label controlCells = controlZone().size();
-
- forAll(pressures_, c)
- {
- totalPressure += pressures_[c];
- }
-
- if (Pstream::parRun())
- {
- reduce(totalPressure, sumOp<scalar>());
-
- reduce(controlCells, sumOp<label>());
- }
-
- if (controlCells > 0)
- {
- totalPressure /= scalar(controlCells);
- }
- }
-
- return totalPressure;
-}
-
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/controllers/basic/stateController/stateController.H b/src/lagrangian/molecularDynamics/molecule/controllers/basic/stateController/stateController.H
deleted file mode 100644
index 3e142c3962e..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/controllers/basic/stateController/stateController.H
+++ /dev/null
@@ -1,270 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2011-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-Class
- stateController
-
-Description
-
- Basic/abstract class of a state controller
-
-SourceFiles
- stateControllerI.H
- stateController.C
- stateControllerIO.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef stateController_H
-#define stateController_H
-
-#include "IOdictionary.H"
-#include "autoPtr.H"
-#include "runTimeSelectionTables.H"
-#include "vector.H"
-#include "volFields.H"
-#include "Random.H"
-#include "polyatomic.H"
-#include "timeData.H"
-#include "writeTimeData.H"
-#include "selectIds.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
-/*---------------------------------------------------------------------------*\
- Class stateController Declaration
-\*---------------------------------------------------------------------------*/
-
-class stateController
-{
-
-protected:
-
- // Protected data
-
- //-
- const fvMesh& mesh_;
-
- //-
- polyatomicCloud& cloud_;
-
- //-
- Random rndGen_;
-
- //- subDictionary containing the properties
- dictionary controllerDict_;
-
- //-
- dictionary timeDict_;
-
- //-
- timeData time_;
-
- //-
- scalar timePeriod_;
-
- //-
- scalar initialTime_;
-
- //- name of control zone
- word regionName_;
-
- //-
- label regionId_;
-
- //-
- bool control_;
-
- //-
- bool readStateFromFile_;
-
- //- set all the properties below from model if required
-
- //-
- bool singleValueController_;
-
- //- target values
- scalar density_;
- vector velocity_;
- scalar temperature_;
- scalar pressure_;
- tensor strainRate_;
- vector tempGradient_;
-
- //- set this in model
- bool fieldController_;
-
- //- targeted fields
- scalarField densities_;
- vectorField velocities_;
- scalarField temperatures_;
- scalarField pressures_;
-
- // set these in model
- bool writeInTimeDir_;
- bool writeInCase_;
-
-
-public:
-
- //- Runtime type information
- TypeName("stateController");
-
- //- Declare runtime constructor selection table
- declareRunTimeSelectionTable
- (
- autoPtr,
- stateController,
- dictionary,
- (
- polyatomicCloud& cloud,
- const dictionary& dict
- ),
- (t, cloud, dict)
- );
-
- // Constructors
-
- //- Construct from components
- stateController
- (
- polyatomicCloud& cloud,
- const dictionary& dict
- );
-
-
- // Selectors
-
- static autoPtr<stateController> New
- (
- polyatomicCloud& cloud,
- const dictionary& dict
- );
-
-
- // Destructor
-
- virtual ~stateController();
-
-
- // Member Functions
-
-
- void updateTime();
-
- //- create an initial configuration
- virtual void initialConfiguration() = 0;
-
- //- calculate any required properties
- virtual void calculateProperties() = 0;
-
- //- control molecules at different stages of the integration time-step
- virtual void controlMolsBeg() = 0;
-
- virtual void controlBeforeForces() = 0;
-
- virtual void controlMols() = 0;
-
- virtual void controlMolsEnd() = 0;
-
-
- //- output data
- virtual void output
- (
- const fileName& fixedPathName,
- const fileName& timePath
- ) = 0;
-
- //- update properties from a modified dictionary
- virtual void updateProperties(const dictionary&) = 0;
-
- void updateStateControllerProperties(const dictionary&);
-
- // Access
-
- //- return the control zone cells
- const labelList& controlZone() const;
-
- //- return the control zone name
- const word& regionName() const;
-
- //- return the targeted fields
- scalar density() const;
- scalar& density();
-
- const vector& velocity() const;
- vector& velocity();
-
- scalar temperature() const;
- scalar& temperature();
-
- scalar pressure() const;
- scalar& pressure();
-
- const tensor& strainRate() const;
- tensor& strainRate();
-
- const vector& tempGradient() const;
- vector& tempGradient();
-
- //- return the targeted fields
- const scalarField& densityField() const;
- scalarField& densityField();
-
- const vectorField& velocityField() const;
- vectorField& velocityField();
-
- const scalarField& temperatureField() const;
- scalarField& temperatureField();
-
- const scalarField& pressureField() const;
- scalarField& pressureField();
-
- bool singleValueController() const;
- bool& singleValueController();
-
- bool fieldController() const;
- bool& fieldController();
-
- bool writeInTimeDir() const;
- bool writeInCase() const;
-
- scalar avReqDensity();
- vector avReqVelocity();
- scalar avReqTemperature();
- scalar avReqPressure();
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/molecules/constPropSite/constPropSite.C b/src/lagrangian/molecularDynamics/molecule/molecules/constPropSite/constPropSite.C
deleted file mode 100644
index 08ea263dd6b..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/molecules/constPropSite/constPropSite.C
+++ /dev/null
@@ -1,114 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2011-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "constPropSite.H"
-
-// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
-
-Foam::constPropSite::constPropSite()
-:
- siteReferencePosition_(vector::zero),
- siteMass_(0.0),
- siteCharge_(0.0),
- siteId_(0),
- name_(),
- pairPotentialSite_(false),
- electrostaticSite_(false)
-{}
-
-
-Foam::constPropSite::constPropSite
-(
- const vector& siteReferencePosition,
- const scalar& siteMass,
- const scalar& siteCharge,
- const label& siteId,
- const word& name,
- const bool& pairPotentialSite,
- const bool& electrostaticSite
-)
-:
- siteReferencePosition_(siteReferencePosition),
- siteMass_(siteMass),
- siteCharge_(siteCharge),
- siteId_(siteId),
- name_(name),
- pairPotentialSite_(pairPotentialSite),
- electrostaticSite_(electrostaticSite)
-{}
-
-
-// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
-
-Foam::constPropSite::~constPropSite()
-{}
-
-
-// * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
-
-Foam::Istream& Foam::operator>>(Istream& is, constPropSite& cPS)
-{
- is >> cPS.siteReferencePosition_
- >> cPS.siteMass_
- >> cPS.siteCharge_
- >> cPS.siteId_
- >> cPS.name_
- >> cPS.pairPotentialSite_
- >> cPS.electrostaticSite_;
-
- // Check state of Istream
- is.check
- (
- "Foam::Istream& Foam::operator>>"
- "(Foam::Istream&, Foam::constPropSite&)"
- );
-
- return is;
-}
-
-
-Foam::Ostream& Foam::operator<<(Ostream& os, const constPropSite& cPS)
-{
-
- os << token::SPACE << cPS.siteReferencePosition()
- << token::SPACE << cPS.siteMass()
- << token::SPACE << cPS.siteCharge()
- << token::SPACE << cPS.siteId()
- << token::SPACE << cPS.name()
- << token::SPACE << cPS.pairPotentialSite()
- << token::SPACE << cPS.electrostaticSite();
-
- // Check state of Ostream
- os.check
- (
- "Foam::Ostream& Foam::operator<<(Foam::Ostream&, "
- "const Foam::constPropSite&)"
- );
-
- return os;
-}
-
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/molecules/constPropSite/constPropSite.H b/src/lagrangian/molecularDynamics/molecule/molecules/constPropSite/constPropSite.H
deleted file mode 100644
index 910c0ccaef4..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/molecules/constPropSite/constPropSite.H
+++ /dev/null
@@ -1,185 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2011-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-Class
-
-
-Description
-
-
-SourceFiles
- constPropSiteI.H
- constPropSite.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef constPropSite_H
-#define constPropSite_H
-
-#include "vector.H"
-#include "IOstreams.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
-// Forward declaration of classes
-class Istream;
-class Ostream;
-
-// Forward declaration of friend functions and operators
-class constPropSite;
-Istream& operator>>(Istream&, constPropSite&);
-Ostream& operator<<(Ostream&, const constPropSite&);
-
-
-/*---------------------------------------------------------------------------*\
- Class constPropSite Declaration
-\*---------------------------------------------------------------------------*/
-
-class constPropSite
-{
- // Private data
-
- //-
- vector siteReferencePosition_;
-
- //-
- scalar siteMass_;
-
- //-
- scalar siteCharge_;
-
- //-
- label siteId_;
-
- //-
- word name_;
-
- //-
- bool pairPotentialSite_;
-
- //-
- bool electrostaticSite_;
-
-
-public:
-
- // Constructors
-
- //- Construct null
- constPropSite();
-
- //- Construct from components
- constPropSite
- (
- const vector& siteReferencePosition,
- const scalar& siteMass,
- const scalar& siteCharge,
- const label& siteId,
- const word& name,
- const bool& pairPotentialSite,
- const bool& electrostaticSite
- );
-
-
- // Selectors
-
- // Destructor
-
- ~constPropSite();
-
- // Member Functions
-
- // Access
-
- //-
- inline const vector& siteReferencePosition() const;
-
- //-
- inline vector& siteReferencePosition();
-
- //-
- inline const scalar& siteMass() const;
-
- //-
- inline scalar& siteMass();
-
- //-
- inline const scalar& siteCharge() const;
-
- //-
- inline scalar& siteCharge();
-
- //-
- inline const label& siteId() const;
-
- //-
- inline label& siteId();
-
- //-
- inline const word& name() const;
-
- //-
- inline word& name();
-
- //-
- inline const bool& pairPotentialSite() const;
-
- //-
- inline bool& pairPotentialSite();
-
- //-
- inline const bool& electrostaticSite() const;
-
- //-
- inline bool& electrostaticSite();
-
-
- // Member Operators
-
- inline bool operator==(const constPropSite&) const;
- inline bool operator!=(const constPropSite&) const;
-
-
- // IOstream Operators
-
- friend Istream& operator>>(Istream&, constPropSite&);
- friend Ostream& operator<<(Ostream&, const constPropSite&);
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#include "constPropSiteI.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/molecules/constPropSite/constPropSiteI.H b/src/lagrangian/molecularDynamics/molecule/molecules/constPropSite/constPropSiteI.H
deleted file mode 100644
index 9c648c94c8b..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/molecules/constPropSite/constPropSiteI.H
+++ /dev/null
@@ -1,139 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2011-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
-
-inline const Foam::vector& Foam::constPropSite::siteReferencePosition() const
-{
- return siteReferencePosition_;
-}
-
-
-inline Foam::vector& Foam::constPropSite::siteReferencePosition()
-{
- return siteReferencePosition_;
-}
-
-
-inline const Foam::scalar& Foam::constPropSite::siteMass() const
-{
- return siteMass_;
-}
-
-
-inline Foam::scalar& Foam::constPropSite::siteMass()
-{
- return siteMass_;
-}
-
-
-inline const Foam::scalar& Foam::constPropSite::siteCharge() const
-{
- return siteCharge_;
-}
-
-
-inline Foam::scalar& Foam::constPropSite::siteCharge()
-{
- return siteCharge_;
-}
-
-
-inline const Foam::label& Foam::constPropSite::siteId() const
-{
- return siteId_;
-}
-
-
-inline Foam::label& Foam::constPropSite::siteId()
-{
- return siteId_;
-}
-
-
-inline const Foam::word& Foam::constPropSite::name() const
-{
- return name_;
-}
-
-
-inline Foam::word& Foam::constPropSite::name()
-{
- return name_;
-}
-
-
-inline const bool& Foam::constPropSite::pairPotentialSite() const
-{
- return pairPotentialSite_;
-}
-
-
-inline bool& Foam::constPropSite::pairPotentialSite()
-{
- return pairPotentialSite_;
-}
-
-
-inline const bool& Foam::constPropSite::electrostaticSite() const
-{
- return electrostaticSite_;
-}
-
-
-inline bool& Foam::constPropSite::electrostaticSite()
-{
- return electrostaticSite_;
-}
-
-
-// * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * * //
-
-bool Foam::constPropSite::operator==
-(
- const constPropSite& rhs
-) const
-{
- return
- siteReferencePosition_ == rhs.siteReferencePosition_
- && siteMass_ == rhs.siteMass_
- && siteCharge_ == rhs.siteCharge_
- && siteId_ == rhs.siteId_
- && name_ == rhs.name_
- && pairPotentialSite_ == rhs.pairPotentialSite_
- && electrostaticSite_ == rhs.electrostaticSite_;
-}
-
-
-bool Foam::constPropSite::operator!=
-(
- const constPropSite& rhs
-) const
-{
- return !(*this == rhs);
-}
-
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/molecules/monoatomic/monoatomic.C b/src/lagrangian/molecularDynamics/molecule/molecules/monoatomic/monoatomic.C
deleted file mode 100644
index f5b03164878..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/molecules/monoatomic/monoatomic.C
+++ /dev/null
@@ -1,199 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2008-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-\*----------------------------------------------------------------------------*/
-
-#include "monoatomic.H"
-#include "Random.H"
-#include "Time.H"
-
-// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
-
-namespace Foam
-{
- defineTypeNameAndDebug(monoatomic, 0);
- defineTemplateTypeNameAndDebug(Cloud<monoatomic>, 0);
-}
-
-// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
-
-bool Foam::monoatomic::move
-(
- monoatomic::trackingData& td,
- const scalar trackTime
-)
-{
- td.switchProcessor = false;
- td.keepParticle = true;
-
- if (special_ != SPECIAL_FROZEN)
- {
- return td.keepParticle;
- }
-
- const constantProperties& constProps(td.cloud().constProps(id_));
-
- if (td.part() == trackingData::tpFirstVelocityHalfStep)
- {
- // First leapfrog velocity adjust part, required before tracking+force
- // part
-
- v_ += 0.5*trackTime*a_;
- }
- else if (td.part() == trackingData::tpLinearTrack)
- {
- // Leapfrog tracking part
-
- scalar tEnd = (1.0 - stepFraction())*trackTime;
- scalar dtMax = tEnd;
-
- while (td.keepParticle && !td.switchProcessor && tEnd > ROOTVSMALL)
- {
- // set the lagrangian time-step
- scalar dt = min(dtMax, tEnd);
-
- dt *= trackToFace(position() + dt*v_, td);
-
- tEnd -= dt;
- stepFraction() = 1.0 - tEnd/trackTime;
- }
-
- setSitePositions(constProps);
- }
- else if (td.part() == trackingData::tpSecondVelocityHalfStep)
- {
- // Second leapfrog velocity adjust part, required after tracking+force
- // part
-
- a_ = siteForces_[0]/constProps.mass();
-
- v_ += 0.5*trackTime*a_;
- }
- else if (td.part() != trackingData::tpRotationalTrack)
- {
- FatalErrorIn("monoatomic::move(trackingData&, const scalar)") << nl
- << td.part() << " is an invalid part of the integration method."
- << abort(FatalError);
- }
-
- return td.keepParticle;
-}
-
-
-void Foam::monoatomic::transformProperties(const tensor& T)
-{
- particle::transformProperties(T);
-
- v_ = transform(T, v_);
-
- a_ = transform(T, a_);
-
- rf_ = transform(T, rf_);
-
- sitePositions_[0] = position_ + (T & (sitePositions_[0] - position_));
-
- siteForces_[0] = T & siteForces_[0];
-}
-
-
-void Foam::monoatomic::transformProperties(const vector& separation)
-{
- particle::transformProperties(separation);
-
- if (special_ == SPECIAL_TETHERED)
- {
- specialPosition_ += separation;
- }
-
- sitePositions_[0] += separation;
-}
-
-
-void Foam::monoatomic::setSitePositions(const constantProperties& constProps)
-{
- sitePositions_[0] = position_;
-}
-
-
-void Foam::monoatomic::setSiteSizes(label size)
-{
- // Nothing required, size controlled internally
-}
-
-
-bool Foam::monoatomic::hitPatch
-(
- const polyPatch&,
- trackingData&,
- const label,
- const scalar,
- const tetIndices&
-)
-{
- return false;
-}
-
-
-void Foam::monoatomic::hitProcessorPatch
-(
- const processorPolyPatch&,
- trackingData& td
-)
-{
- td.switchProcessor = true;
-}
-
-
-void Foam::monoatomic::hitWallPatch
-(
- const wallPolyPatch& wpp,
- trackingData& td,
- const tetIndices& tetIs
-)
-{
- // Use of the normal from tetIs is not required as
- // hasWallImpactDistance is false.
- vector nw = normal();
- nw /= mag(nw);
-
- scalar vn = v_ & nw;
-
- // Specular reflection
- if (vn > 0)
- {
- v_ -= 2*vn*nw;
- }
-}
-
-
-void Foam::monoatomic::hitPatch
-(
- const polyPatch&,
- trackingData& td
-)
-{
- td.keepParticle = false;
-}
-
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/molecules/monoatomic/monoatomic.H b/src/lagrangian/molecularDynamics/molecule/molecules/monoatomic/monoatomic.H
deleted file mode 100644
index d1e999435b9..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/molecules/monoatomic/monoatomic.H
+++ /dev/null
@@ -1,418 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2008-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-Class
- Foam::monoatomic
-
-Description
- Foam::monoatomic
-
-SourceFiles
- monoatomicI.H
- monoatomic.C
- monoatomicIO.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef monoatomic_H
-#define monoatomic_H
-
-#include "particle.H"
-#include "IOstream.H"
-#include "autoPtr.H"
-#include "diagTensor.H"
-#include "constPropSite.H"
-#include "MoleculeCloud.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
-/*---------------------------------------------------------------------------*\
- Class monoatomic Declaration
-\*---------------------------------------------------------------------------*/
-
-class monoatomic
-:
- public particle
-{
-
-public:
-
- // Values of special that are less than zero are for built-in functionality.
- // Values greater than zero are user specifiable/expandable (i.e. test
- // special_ >= SPECIAL_USER)
-
- enum specialTypes
- {
- SPECIAL_TETHERED = -1,
- SPECIAL_FROZEN = -2,
- NOT_SPECIAL = 0,
- SPECIAL_USER = 1
- };
-
- //- Class to hold monoatomic constant properties
- class constantProperties
- {
- // Private data
-
- //- Sites of mass, charge or interaction
- FixedList<constPropSite, 1> sites_;
-
- //- Which sites require electrostatic interactions
- FixedList<label, 1> electrostaticSites_;
-
- //- Which sites require pair interactions
- FixedList<label, 1> pairPotSites_;
-
- //-
- scalar mass_;
-
-
- public:
-
- //-
- inline constantProperties();
-
- //- Construct from dictionary
- inline constantProperties
- (
- const dictionary& dict,
- const List<label>& siteIds
- );
-
- // Member functions
-
- //-
- inline const FixedList<constPropSite, 1>& sites() const;
-
- //-
- inline const FixedList<label, 1>& pairPotSites() const;
-
- //-
- inline const FixedList<label, 1>& electrostaticSites() const;
-
- //-
- inline label degreesOfFreedom() const;
-
- //-
- inline scalar mass() const;
-
- //-
- inline label nSites() const;
- };
-
-
- //- Class used to pass tracking data to the trackToFace function
- class trackingData
- :
- public particle::TrackingData<MoleculeCloud<monoatomic> >
- {
- public:
-
- enum trackPart
- {
- tpFirstVelocityHalfStep,
- tpLinearTrack,
- tpRotationalTrack,
- tpSecondVelocityHalfStep,
- tpAccess
- };
-
- private:
-
- // label specifying which part of the integration algorithm is taking
- label part_;
-
-
- public:
-
- // Constructors
-
- trackingData(MoleculeCloud<monoatomic>& cloud, trackPart part)
- :
- particle::TrackingData<MoleculeCloud<monoatomic> >(cloud),
- part_(part)
- {}
-
- // Member functions
-
- inline label part() const
- {
- return part_;
- }
- };
-
-
-private:
-
- // Private data
-
- //- Linear velocity of monoatomic
- vector v_;
-
- //- Total linear acceleration of monoatomic
- vector a_;
-
- //-
- vector specialPosition_;
-
- //-
- scalar potentialEnergy_;
-
- // - r_ij f_ij, stress dyad
- tensor rf_;
-
-// // - r_ij outer product f_ij: virial contribution
-// tensor rDotf_;
-
- //-
- label special_;
-
- //-
- label id_;
-
- //-
- FixedList<vector, 1> siteForces_;
-
- //-
- FixedList<vector, 1> sitePositions_;
-
-
-public:
-
- //- Runtime type information
- TypeName("monoatomic");
-
- friend class Cloud<monoatomic>;
-
- // Constructors
-
- //- Construct with macroscopic description
- inline monoatomic
- (
- const polyMesh& mesh,
- const vector& position,
- const label cellI,
- const label tetFaceI,
- const label tetPtI,
- const scalar temperature,
- const vector& bulkVelocity,
- const vector& specialPosition,
- const constantProperties& constProps,
- trackingData& td,
- const label special,
- const label id
- );
-
- //- Construct from all components
- inline monoatomic
- (
- const polyMesh& mesh,
- const vector& position,
- const label cellI,
- const label tetFaceI,
- const label tetPtI,
- const vector& v,
- const vector& a,
- const vector& specialPosition,
- const constantProperties& constProps,
- const label special,
- const label id
- );
-
- //- Construct from Istream
- monoatomic
- (
- const polyMesh& mesh,
- Istream& is,
- bool readFields = true
- );
-
- //- Construct and return a clone
- autoPtr<particle> clone() const
- {
- return autoPtr<particle>(new monoatomic(*this));
- }
-
- //- Factory class to read-construct particles used for
- // parallel transfer
- class iNew
- {
- const polyMesh& mesh_;
-
- public:
-
- iNew(const polyMesh& mesh)
- :
- mesh_(mesh)
- {}
-
- autoPtr<monoatomic> operator()(Istream& is) const
- {
- return autoPtr<monoatomic>(new monoatomic(mesh_, is, true));
- }
- };
-
-
- // Member Functions
-
- // Tracking
-
- //-
- bool move(trackingData&, const scalar trackTime);
-
- //-
- virtual void transformProperties(const tensor& T);
-
- //-
- virtual void transformProperties(const vector& separation);
-
- //-
- void setSitePositions(const constantProperties& constProps);
-
- //-
- void setSiteSizes(label size);
-
-
- // Access
-
- //-
- inline const vector& v() const;
-
- //-
- inline vector& v();
-
- //-
- inline const vector& a() const;
-
- //-
- inline vector& a();
-
- //-
- inline const FixedList<vector, 1>& siteForces() const;
-
- //-
- inline FixedList<vector, 1>& siteForces();
-
- //-
- inline const FixedList<vector, 1>& sitePositions() const;
-
- //-
- inline FixedList<vector, 1>& sitePositions();
-
- //-
- inline const vector& specialPosition() const;
-
- //-
- inline vector& specialPosition();
-
- //-
- inline scalar potentialEnergy() const;
-
- //-
- inline scalar& potentialEnergy();
-
- //-
- inline const tensor& rf() const;
-
- //-
- inline tensor& rf();
-
- //-
- inline label special() const;
-
- //-
- inline bool tethered() const;
-
- //-
- inline label id() const;
-
-
- // Member Operators
-
- //- Overridable function to handle the particle hitting a patch
- // Executed before other patch-hitting functions
- bool hitPatch
- (
- const polyPatch&,
- trackingData& td,
- const label patchI,
- const scalar trackFraction,
- const tetIndices& tetIs
- );
-
- //- Overridable function to handle the particle hitting a processorPatch
- void hitProcessorPatch
- (
- const processorPolyPatch&,
- trackingData& td
- );
-
- //- Overridable function to handle the particle hitting a wallPatch
- void hitWallPatch
- (
- const wallPolyPatch&,
- trackingData& td,
- const tetIndices&
- );
-
- //- Overridable function to handle the particle hitting a polyPatch
- void hitPatch
- (
- const polyPatch&,
- trackingData& td
- );
-
-
- // I-O
-
- //- Read
- static void readFields(Cloud<monoatomic>& mC);
-
- //- Write
- static void writeFields(const Cloud<monoatomic>& mC);
-
- //- Show info
- static void info(trackingData& td);
-
-
- // IOstream Operators
-
- friend Ostream& operator<<(Ostream&, const monoatomic&);
-};
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#include "monoatomicI.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/molecules/monoatomic/monoatomicI.H b/src/lagrangian/molecularDynamics/molecule/molecules/monoatomic/monoatomicI.H
deleted file mode 100644
index 4485d53b8dd..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/molecules/monoatomic/monoatomicI.H
+++ /dev/null
@@ -1,328 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2008-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
-
-inline Foam::monoatomic::constantProperties::constantProperties()
-:
- sites_(),
- electrostaticSites_(-1),
- pairPotSites_(-1),
- mass_(0)
-{}
-
-
-inline Foam::monoatomic::constantProperties::constantProperties
-(
- const dictionary& dict,
- const List<label>& siteIds
-)
-:
- sites_(),
- electrostaticSites_(-1),
- pairPotSites_(-1),
- mass_(0)
-{
- if (siteIds.size() != 1)
- {
- FatalErrorIn
- (
- "Foam::monoatomic::constantProperties::constantProperties"
- "("
- "const dictionary& dict, "
- "const List<label>& siteIds"
- ")"
- )
- << "monoatomic, single site only, given: " << dict
- << nl << abort(FatalError);
- }
-
- FixedList<word, 1> siteIdNames = dict.lookup("siteIds");
-
- FixedList<vector, 1> siteReferencePositions
- (
- dict.lookup("siteReferencePositions")
- );
-
- FixedList<scalar, 1> siteMasses(dict.lookup("siteMasses"));
-
- FixedList<scalar, 1> siteCharges(dict.lookup("siteCharges"));
-
- FixedList<word, 1> pairPotentialIds(dict.lookup("pairPotentialSiteIds"));
-
- constPropSite site = sites_[0];
-
- site = constPropSite
- (
- siteReferencePositions[0],
- siteMasses[0],
- siteCharges[0],
- siteIds[0],
- siteIdNames[0],
- (findIndex(pairPotentialIds, siteIdNames[0]) != -1), // pair
- (mag(siteCharges[0]) > VSMALL) // charge
- );
-
- mass_ = site.siteMass();
-
- if (site.pairPotentialSite())
- {
- pairPotSites_[0] = 0;
- }
-
- if (site.electrostaticSite())
- {
- electrostaticSites_[0] = 0;
- }
-
- if (!site.pairPotentialSite() && !site.electrostaticSite())
- {
- WarningIn
- (
- "Foam::monoatomic::constantProperties::constantProperties"
- "("
- "const dictionary& dict"
- ")"
- )
- << siteIdNames[0] << " is a non-interacting site." << endl;
- }
-
- // Single site monoatomic - no rotational motion.
- site.siteReferencePosition() = vector::zero;
-}
-
-
-inline Foam::monoatomic::monoatomic
-(
- const polyMesh& mesh,
- const vector& position,
- const label cellI,
- const label tetFaceI,
- const label tetPtI,
- const scalar temperature,
- const vector& bulkVelocity,
- const vector& specialPosition,
- const constantProperties& cP,
- monoatomic::trackingData& td,
- const label special,
- const label id
-)
-:
- particle(mesh, position, cellI, tetFaceI, tetPtI),
- v_(bulkVelocity),
- a_(vector::zero),
- specialPosition_(specialPosition),
- potentialEnergy_(0.0),
- rf_(tensor::zero),
- special_(special),
- id_(id),
- siteForces_(vector::zero),
- sitePositions_()
-{
- setSitePositions(cP);
-
- v_ += td.cloud().equipartitionLinearVelocity(temperature, cP.mass());
-}
-
-
-inline Foam::monoatomic::monoatomic
-(
- const polyMesh& mesh,
- const vector& position,
- const label cellI,
- const label tetFaceI,
- const label tetPtI,
- const vector& v,
- const vector& a,
- const vector& specialPosition,
- const constantProperties& cP,
- const label special,
- const label id
-
-)
-:
- particle(mesh, position, cellI, tetFaceI, tetPtI),
- v_(v),
- a_(a),
- specialPosition_(specialPosition),
- potentialEnergy_(0.0),
- rf_(tensor::zero),
- special_(special),
- id_(id),
- siteForces_(vector::zero),
- sitePositions_()
-{
- setSitePositions(cP);
-}
-
-// * * * * * * * constantProperties Member Functions * * * * * * * * * * * * //
-
-inline const Foam::FixedList<Foam::constPropSite, 1>&
-Foam::monoatomic::constantProperties::sites() const
-{
- return sites_;
-}
-
-
-inline const Foam::FixedList<Foam::label, 1>&
-Foam::monoatomic::constantProperties::pairPotSites() const
-{
- return pairPotSites_;
-}
-
-
-inline const Foam::FixedList<Foam::label, 1>&
-Foam::monoatomic::constantProperties::electrostaticSites() const
-{
- return electrostaticSites_;
-}
-
-
-inline Foam::label
-Foam::monoatomic::constantProperties::degreesOfFreedom() const
-{
- return 3;
-}
-
-
-inline Foam::scalar Foam::monoatomic::constantProperties::mass() const
-{
- return mass_;
-}
-
-
-inline Foam::label Foam::monoatomic::constantProperties::nSites() const
-{
- return 1;
-}
-
-
-// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
-
-inline const Foam::vector& Foam::monoatomic::v() const
-{
- return v_;
-}
-
-
-inline Foam::vector& Foam::monoatomic::v()
-{
- return v_;
-}
-
-
-inline const Foam::vector& Foam::monoatomic::a() const
-{
- return a_;
-}
-
-
-inline Foam::vector& Foam::monoatomic::a()
-{
- return a_;
-}
-
-
-inline const Foam::FixedList<Foam::vector, 1>&
-Foam::monoatomic::siteForces() const
-{
- return siteForces_;
-}
-
-
-inline Foam::FixedList<Foam::vector, 1>& Foam::monoatomic::siteForces()
-{
- return siteForces_;
-}
-
-
-inline const Foam::FixedList<Foam::vector, 1>&
-Foam::monoatomic::sitePositions() const
-{
- return sitePositions_;
-}
-
-
-inline Foam::FixedList<Foam::vector, 1>& Foam::monoatomic::sitePositions()
-{
- return sitePositions_;
-}
-
-
-inline const Foam::vector& Foam::monoatomic::specialPosition() const
-{
- return specialPosition_;
-}
-
-
-inline Foam::vector& Foam::monoatomic::specialPosition()
-{
- return specialPosition_;
-}
-
-
-inline Foam::scalar Foam::monoatomic::potentialEnergy() const
-{
- return potentialEnergy_;
-}
-
-
-inline Foam::scalar& Foam::monoatomic::potentialEnergy()
-{
- return potentialEnergy_;
-}
-
-
-inline const Foam::tensor& Foam::monoatomic::rf() const
-{
- return rf_;
-}
-
-
-inline Foam::tensor& Foam::monoatomic::rf()
-{
- return rf_;
-}
-
-
-inline Foam::label Foam::monoatomic::special() const
-{
- return special_;
-}
-
-
-inline bool Foam::monoatomic::tethered() const
-{
- return special_ == SPECIAL_TETHERED;
-}
-
-
-inline Foam::label Foam::monoatomic::id() const
-{
- return id_;
-}
-
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/molecules/monoatomic/monoatomicIO.C b/src/lagrangian/molecularDynamics/molecule/molecules/monoatomic/monoatomicIO.C
deleted file mode 100644
index a9a27bb9929..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/molecules/monoatomic/monoatomicIO.C
+++ /dev/null
@@ -1,305 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2008-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "monoatomic.H"
-#include "IOstreams.H"
-#include "Cloud.H"
-
-// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
-
-Foam::monoatomic::monoatomic
-(
- const polyMesh& mesh,
- Istream& is,
- bool readFields
-)
-:
- particle(mesh, is, readFields),
- v_(vector::zero),
- a_(vector::zero),
- specialPosition_(vector::zero),
- potentialEnergy_(0.0),
- rf_(tensor::zero),
- special_(0),
- id_(0),
- siteForces_(),
- sitePositions_()
-{
- if (readFields)
- {
- if (is.format() == IOstream::ASCII)
- {
- is >> v_;
- is >> a_;
- is >> siteForces_;
- potentialEnergy_ = readScalar(is);
- is >> rf_;
- special_ = readLabel(is);
- id_ = readLabel(is);
- is >> sitePositions_;
- is >> specialPosition_;
- }
- else
- {
- is.read
- (
- reinterpret_cast<char*>(&v_),
- sizeof(v_)
- + sizeof(a_)
- + sizeof(specialPosition_)
- + sizeof(potentialEnergy_)
- + sizeof(rf_)
- + sizeof(special_)
- + sizeof(id_)
- );
-
- is >> siteForces_ >> sitePositions_;
- }
- }
-
- // Check state of Istream
- is.check
- (
- "Foam::monoatomic::monoatomic"
- "("
- "const polyMesh& mesh,"
- "Istream& is,"
- "bool readFields"
- ")"
- );
-}
-
-
-void Foam::monoatomic::readFields(Cloud<monoatomic>& mC)
-{
- if (!mC.size())
- {
- return;
- }
-
- particle::readFields(mC);
-
- IOField<vector> v(mC.fieldIOobject("v", IOobject::MUST_READ));
- mC.checkFieldIOobject(mC, v);
-
- IOField<vector> a(mC.fieldIOobject("a", IOobject::MUST_READ));
- mC.checkFieldIOobject(mC, a);
-
- IOField<vector> specialPosition
- (
- mC.fieldIOobject("specialPosition", IOobject::MUST_READ)
- );
- mC.checkFieldIOobject(mC, specialPosition);
-
- IOField<label> special(mC.fieldIOobject("special", IOobject::MUST_READ));
- mC.checkFieldIOobject(mC, special);
-
- IOField<label> id(mC.fieldIOobject("id", IOobject::MUST_READ));
- mC.checkFieldIOobject(mC, id);
-
- label i = 0;
-
- forAllIter(typename Cloud<monoatomic>, mC, iter)
- {
- monoatomic& mol = iter();
-
- mol.v_ = v[i];
- mol.a_ = a[i];
- mol.specialPosition_ = specialPosition[i];
- mol.special_ = special[i];
- mol.id_ = id[i];
- i++;
- }
-}
-
-
-void Foam::monoatomic::writeFields(const Cloud<monoatomic>& mC)
-{
- particle::writeFields(mC);
-
- label np = mC.size();
-
- IOField<vector> v(mC.fieldIOobject("v", IOobject::NO_READ), np);
- IOField<vector> a(mC.fieldIOobject("a", IOobject::NO_READ), np);
- IOField<vector> specialPosition
- (
- mC.fieldIOobject("specialPosition", IOobject::NO_READ),
- np
- );
- IOField<label> special(mC.fieldIOobject("special", IOobject::NO_READ), np);
- IOField<label> id(mC.fieldIOobject("id", IOobject::NO_READ), np);
-
- label i = 0;
- forAllConstIter(typename Cloud<monoatomic>, mC, iter)
- {
- const monoatomic& mol = iter();
-
- v[i] = mol.v_;
- a[i] = mol.a_;
- specialPosition[i] = mol.specialPosition_;
- special[i] = mol.special_;
- id[i] = mol.id_;
- i++;
- }
-
- v.write();
- a.write();
- specialPosition.write();
- special.write();
- id.write();
-}
-
-
-void Foam::monoatomic::info(monoatomic::trackingData& td)
-{
- vector totalLinearMomentum(vector::zero);
- scalar maxVelocityMag = 0.0;
- scalar totalMass = 0.0;
- scalar totalLinearKE = 0.0;
- scalar totalPE = 0.0;
- scalar totalrDotf = 0.0;
- label nMols = td.cloud().size();
-
- forAllConstIter(typename Cloud<monoatomic>, td.cloud(), mol)
- {
- const label molId = mol().id();
- scalar molMass(td.cloud().constProps(molId).mass());
- totalMass += molMass;
- }
-
- forAllConstIter(typename Cloud<monoatomic>, td.cloud(), mol)
- {
- const label molId = mol().id();
- const monoatomic::constantProperties cP
- (
- td.cloud().constProps(molId)
- );
- scalar molMass(cP.mass());
- const vector& molV(mol().v());
- totalLinearMomentum += molV * molMass;
- if (mag(molV) > maxVelocityMag)
- {
- maxVelocityMag = mag(molV);
- }
- totalLinearKE += 0.5*molMass*magSqr(molV);
- totalPE += mol().potentialEnergy();
- totalrDotf += tr(mol().rf());
- }
-
- scalar meshVolume = sum(td.cloud().mesh().cellVolumes());
-
- if (Pstream::parRun())
- {
- reduce(totalLinearMomentum, sumOp<vector>());
- reduce(maxVelocityMag, maxOp<scalar>());
- reduce(totalMass, sumOp<scalar>());
- reduce(totalLinearKE, sumOp<scalar>());
- reduce(totalPE, sumOp<scalar>());
- reduce(totalrDotf, sumOp<scalar>());
- reduce(nMols, sumOp<label>());
- reduce(meshVolume, sumOp<scalar>());
- }
-
- if (nMols)
- {
- Info<< nl << "Number of molecules in " << td.cloud().name() << " = "
- << nMols << nl
- << " Overall number density = "
- << nMols/meshVolume << nl
- << " Overall mass density = "
- << totalMass/meshVolume << nl
- << " Average linear momentum per molecule = "
- << totalLinearMomentum/nMols << ' '
- << mag(totalLinearMomentum)/nMols << nl
- << " maximum |velocity| = "
- << maxVelocityMag << nl
- << " Average linear KE per molecule = "
- << totalLinearKE/nMols << nl
- << " Average angular KE per molecule = "
- << totalPE/nMols << nl
- << " Average TE per molecule = "
- <<
- (
- totalLinearKE
- + totalPE
- )
- /nMols
- << nl << endl;
- }
- else
- {
- Info<< nl << "No molecules in " << td.cloud().name() << endl;
- }
-}
-
-
-// * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
-
-Foam::Ostream& Foam::operator<<(Ostream& os, const monoatomic& mol)
-{
- if (os.format() == IOstream::ASCII)
- {
- os << token::SPACE << static_cast<const particle&>(mol)
- << token::SPACE << mol.face()
- << token::SPACE << mol.stepFraction()
- << token::SPACE << mol.v_
- << token::SPACE << mol.a_
- << token::SPACE << mol.specialPosition_
- << token::SPACE << mol.potentialEnergy_
- << token::SPACE << mol.rf_
- << token::SPACE << mol.special_
- << token::SPACE << mol.id_
- << token::SPACE << mol.siteForces_
- << token::SPACE << mol.sitePositions_;
- }
- else
- {
- os << static_cast<const particle&>(mol);
- os.write
- (
- reinterpret_cast<const char*>(&mol.v_),
- sizeof(mol.v_)
- + sizeof(mol.a_)
- + sizeof(mol.specialPosition_)
- + sizeof(mol.potentialEnergy_)
- + sizeof(mol.rf_)
- + sizeof(mol.special_)
- + sizeof(mol.id_)
- );
- os << mol.siteForces_ << mol.sitePositions_;
- }
-
- // Check state of Ostream
- os.check
- (
- "Foam::Ostream& Foam::operator<<"
- "(Foam::Ostream&, const Foam::monoatomic&)"
- );
-
- return os;
-}
-
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/molecules/polyatomic/polyatomic.C b/src/lagrangian/molecularDynamics/molecule/molecules/polyatomic/polyatomic.C
deleted file mode 100644
index 18fbfbb5c9c..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/molecules/polyatomic/polyatomic.C
+++ /dev/null
@@ -1,320 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2008-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-\*----------------------------------------------------------------------------*/
-
-#include "polyatomic.H"
-#include "Random.H"
-#include "Time.H"
-
-// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
-
-namespace Foam
-{
- defineTypeNameAndDebug(polyatomic, 0);
- defineTemplateTypeNameAndDebug(Cloud<polyatomic>, 0);
-}
-
-
-// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
-
-Foam::tensor Foam::polyatomic::rotationTensorX(scalar phi) const
-{
- return tensor
- (
- 1, 0, 0,
- 0, Foam::cos(phi), -Foam::sin(phi),
- 0, Foam::sin(phi), Foam::cos(phi)
- );
-}
-
-
-Foam::tensor Foam::polyatomic::rotationTensorY(scalar phi) const
-{
- return tensor
- (
- Foam::cos(phi), 0, Foam::sin(phi),
- 0, 1, 0,
- -Foam::sin(phi), 0, Foam::cos(phi)
- );
-}
-
-
-Foam::tensor Foam::polyatomic::rotationTensorZ(scalar phi) const
-{
- return tensor
- (
- Foam::cos(phi), -Foam::sin(phi), 0,
- Foam::sin(phi), Foam::cos(phi), 0,
- 0, 0, 1
- );
-}
-
-
-// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
-
-bool Foam::polyatomic::move
-(
- polyatomic::trackingData& td,
- const scalar trackTime
-)
-{
- td.switchProcessor = false;
- td.keepParticle = true;
-
- if (special_ != SPECIAL_FROZEN)
- {
- return td.keepParticle;
- }
-
- const constantProperties& constProps(td.cloud().constProps(id_));
-
- if (td.part() == trackingData::tpFirstVelocityHalfStep)
- {
- // First leapfrog velocity adjust part, required before tracking+force
- // part
-
- v_ += 0.5*trackTime*a_;
-
- pi_ += 0.5*trackTime*tau_;
- }
- else if (td.part() == trackingData::tpLinearTrack)
- {
- // Leapfrog tracking part
-
- scalar tEnd = (1.0 - stepFraction())*trackTime;
- scalar dtMax = tEnd;
-
- while (td.keepParticle && !td.switchProcessor && tEnd > ROOTVSMALL)
- {
- // set the lagrangian time-step
- scalar dt = min(dtMax, tEnd);
-
- dt *= trackToFace(position() + dt*v_, td);
-
- tEnd -= dt;
- stepFraction() = 1.0 - tEnd/trackTime;
- }
- }
- else if (td.part() == trackingData::tpRotationalTrack)
- {
- // Leapfrog orientation adjustment, carried out before force calculation
- // but after tracking stage, i.e. rotation carried once linear motion
- // complete.
-
- if (!constProps.pointMolecule())
- {
- const diagTensor& momentOfInertia(constProps.momentOfInertia());
-
- tensor R;
-
- if (!constProps.linearMolecule())
- {
- R = rotationTensorX(0.5*trackTime*pi_.x()/momentOfInertia.xx());
- pi_ = pi_ & R;
- Q_ = Q_ & R;
- }
-
- R = rotationTensorY(0.5*trackTime*pi_.y()/momentOfInertia.yy());
- pi_ = pi_ & R;
- Q_ = Q_ & R;
-
- R = rotationTensorZ(trackTime*pi_.z()/momentOfInertia.zz());
- pi_ = pi_ & R;
- Q_ = Q_ & R;
-
- R = rotationTensorY(0.5*trackTime*pi_.y()/momentOfInertia.yy());
- pi_ = pi_ & R;
- Q_ = Q_ & R;
-
- if (!constProps.linearMolecule())
- {
- R = rotationTensorX(0.5*trackTime*pi_.x()/momentOfInertia.xx());
- pi_ = pi_ & R;
- Q_ = Q_ & R;
- }
- }
-
- setSitePositions(constProps);
- }
- else if (td.part() == trackingData::tpSecondVelocityHalfStep)
- {
- // Second leapfrog velocity adjust part, required after tracking+force
- // part
-
- scalar m = constProps.mass();
-
- a_ = vector::zero;
-
- tau_ = vector::zero;
-
- forAll(siteForces_, sI)
- {
- const vector& f = siteForces_[sI];
-
- a_ += f/m;
-
- tau_ +=
- constProps.sites()[sI].siteReferencePosition()
- ^ (Q_.T() & f);
- }
-
- v_ += 0.5*trackTime*a_;
-
- pi_ += 0.5*trackTime*tau_;
-
- if (constProps.pointMolecule())
- {
- tau_ = vector::zero;
-
- pi_ = vector::zero;
- }
-
- if (constProps.linearMolecule())
- {
- tau_.x() = 0.0;
-
- pi_.x() = 0.0;
- }
- }
- else
- {
- FatalErrorIn("polyatomic::move(trackingData&, const scalar)") << nl
- << td.part() << " is an invalid part of the integration method."
- << abort(FatalError);
- }
-
- return td.keepParticle;
-}
-
-
-void Foam::polyatomic::transformProperties(const tensor& T)
-{
- particle::transformProperties(T);
-
- Q_ = T & Q_;
-
- v_ = transform(T, v_);
-
- a_ = transform(T, a_);
-
- pi_ = Q_.T() & transform(T, Q_ & pi_);
-
- tau_ = Q_.T() & transform(T, Q_ & tau_);
-
- rf_ = transform(T, rf_);
-
- sitePositions_ = position_ + (T & (sitePositions_ - position_));
-
- siteForces_ = T & siteForces_;
-}
-
-
-void Foam::polyatomic::transformProperties(const vector& separation)
-{
- particle::transformProperties(separation);
-
- if (special_ == SPECIAL_TETHERED)
- {
- specialPosition_ += separation;
- }
-
- sitePositions_ = sitePositions_ + separation;
-}
-
-
-void Foam::polyatomic::setSitePositions(const constantProperties& constProps)
-{
- forAll(constProps.sites(), sI)
- {
- sitePositions_[sI] =
- position_
- + (Q_ & constProps.sites()[sI].siteReferencePosition());
- }
-}
-
-
-void Foam::polyatomic::setSiteSizes(label size)
-{
- sitePositions_.setSize(size);
-
- siteForces_.setSize(size);
-}
-
-
-bool Foam::polyatomic::hitPatch
-(
- const polyPatch&,
- trackingData&,
- const label,
- const scalar,
- const tetIndices&
-)
-{
- return false;
-}
-
-
-void Foam::polyatomic::hitProcessorPatch
-(
- const processorPolyPatch&,
- trackingData& td
-)
-{
- td.switchProcessor = true;
-}
-
-
-void Foam::polyatomic::hitWallPatch
-(
- const wallPolyPatch& wpp,
- trackingData& td,
- const tetIndices& tetIs
-)
-{
- // Use of the normal from tetIs is not required as
- // hasWallImpactDistance is false.
- vector nw = normal();
- nw /= mag(nw);
-
- scalar vn = v_ & nw;
-
- // Specular reflection
- if (vn > 0)
- {
- v_ -= 2*vn*nw;
- }
-}
-
-
-void Foam::polyatomic::hitPatch
-(
- const polyPatch&,
- trackingData& td
-)
-{
- td.keepParticle = false;
-}
-
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/molecules/polyatomic/polyatomic.H b/src/lagrangian/molecularDynamics/molecule/molecules/polyatomic/polyatomic.H
deleted file mode 100644
index 0a7797516b7..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/molecules/polyatomic/polyatomic.H
+++ /dev/null
@@ -1,483 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2008-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-Class
- Foam::polyatomic
-
-Description
- Foam::polyatomic
-
-SourceFiles
- polyatomicI.H
- polyatomic.C
- polyatomicIO.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef polyatomic_H
-#define polyatomic_H
-
-#include "particle.H"
-#include "IOstream.H"
-#include "autoPtr.H"
-#include "diagTensor.H"
-#include "constPropSite.H"
-#include "MoleculeCloud.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
-/*---------------------------------------------------------------------------*\
- Class polyatomic Declaration
-\*---------------------------------------------------------------------------*/
-
-class polyatomic
-:
- public particle
-{
-
-public:
-
- // Values of special that are less than zero are for built-in functionality.
- // Values greater than zero are user specifiable/expandable (i.e. test
- // special_ >= SPECIAL_USER)
-
- enum specialTypes
- {
- SPECIAL_TETHERED = -1,
- SPECIAL_FROZEN = -2,
- NOT_SPECIAL = 0,
- SPECIAL_USER = 1
- };
-
- //- Class to hold polyatomic constant properties
- class constantProperties
- {
- // Private data
-
- //- Sites of mass, charge or interaction
- List<constPropSite> sites_;
-
- //- Which sites require electrostatic interactions
- List<label> electrostaticSites_;
-
- //- Which sites require pair interactions
- List<label> pairPotSites_;
-
- //- Moment of intertia (in principal axis configiration)
- diagTensor momentOfInertia_;
-
- //-
- scalar mass_;
-
-
- // Private Member Functions
-
- //-
- bool linearMoleculeTest() const;
-
-
- public:
-
- //-
- inline constantProperties();
-
- //- Construct from dictionary
- inline constantProperties
- (
- const dictionary& dict,
- const labelList& siteIds
- );
-
- // Member functions
-
- //-
- inline const List<constPropSite>& sites() const;
-
- //-
- inline const List<label>& pairPotSites() const;
-
- //-
- inline const List<label>& electrostaticSites() const;
-
- //-
- inline const diagTensor& momentOfInertia() const;
-
- //-
- inline bool linearMolecule() const;
-
- //-
- inline bool pointMolecule() const;
-
- //-
- inline label degreesOfFreedom() const;
-
- //-
- inline scalar mass() const;
-
- //-
- inline label nSites() const;
- };
-
-
- //- Class used to pass tracking data to the trackToFace function
- class trackingData
- :
- public particle::TrackingData<MoleculeCloud<polyatomic> >
- {
- public:
-
- enum trackPart
- {
- tpFirstVelocityHalfStep,
- tpLinearTrack,
- tpRotationalTrack,
- tpSecondVelocityHalfStep,
- tpAccess
- };
-
-
- private:
-
- // label specifying which part of the integration algorithm is taking
- label part_;
-
-
- public:
-
- // Constructors
-
- trackingData(MoleculeCloud<polyatomic>& cloud, trackPart part)
- :
- particle::TrackingData<MoleculeCloud<polyatomic> >(cloud),
- part_(part)
- {}
-
- // Member functions
-
- inline label part() const
- {
- return part_;
- }
- };
-
-
-private:
-
- // Private data
-
- //- Orientation, stored as the rotation tensor to transform
- // from the polyatomic to the global reference frame, i.e.:
- // globalVector = Q_ & bodyLocalVector
- // bodyLocalVector = Q_.T() & globalVector
- tensor Q_;
-
- //- Linear velocity of polyatomic
- vector v_;
-
- //- Total linear acceleration of polyatomic
- vector a_;
-
- //- Angular momentum of polyatomic, in body local reference frame
- vector pi_;
-
- //- Total torque on polyatomic, in body local reference frame
- vector tau_;
-
- //-
- vector specialPosition_;
-
- //-
- scalar potentialEnergy_;
-
- // - r_ij f_ij, stress dyad
- tensor rf_;
-
-// // - r_ij outer product f_ij: virial contribution
-// tensor rDotf_;
-
- //-
- label special_;
-
- //-
- label id_;
-
- //-
- List<vector> siteForces_;
-
- //-
- List<vector> sitePositions_;
-
-
- // Private Member Functions
-
-
- //-
- tensor rotationTensorX(scalar deltaT) const;
-
- //-
- tensor rotationTensorY(scalar deltaT) const;
-
- //-
- tensor rotationTensorZ(scalar deltaT) const;
-
-
-public:
-
- //- Runtime type information
- TypeName("polyatomic");
-
- friend class Cloud<polyatomic>;
-
- // Constructors
-
- //- Construct with macroscopic description
- inline polyatomic
- (
- const polyMesh& mesh,
- const vector& position,
- const label cellI,
- const label tetFaceI,
- const label tetPtI,
- const scalar temperature,
- const vector& bulkVelocity,
- const vector& specialPosition,
- const constantProperties& constProps,
- polyatomic::trackingData& td,
- const label special,
- const label id
- );
-
- //- Construct from all components
- inline polyatomic
- (
- const polyMesh& mesh,
- const vector& position,
- const label cellI,
- const label tetFaceI,
- const label tetPtI,
- const tensor& Q,
- const vector& v,
- const vector& a,
- const vector& pi,
- const vector& tau,
- const vector& specialPosition,
- const constantProperties& constProps,
- const label special,
- const label id
- );
-
- //- Construct from Istream
- polyatomic
- (
- const polyMesh& mesh,
- Istream& is,
- bool readFields = true
- );
-
- //- Construct and return a clone
- autoPtr<particle> clone() const
- {
- return autoPtr<particle>(new polyatomic(*this));
- }
-
- //- Factory class to read-construct particles used for
- // parallel transfer
- class iNew
- {
- const polyMesh& mesh_;
-
- public:
-
- iNew(const polyMesh& mesh)
- :
- mesh_(mesh)
- {}
-
- autoPtr<polyatomic> operator()(Istream& is) const
- {
- return autoPtr<polyatomic>(new polyatomic(mesh_, is, true));
- }
- };
-
-
- // Member Functions
-
- // Tracking
-
- //-
- bool move(trackingData&, const scalar trackTime);
-
- //-
- virtual void transformProperties(const tensor& T);
-
- //-
- virtual void transformProperties(const vector& separation);
-
- //-
- void setSitePositions(const constantProperties& constProps);
-
- //-
- void setSiteSizes(label size);
-
-
- // Access
-
- //-
- inline const tensor& Q() const;
-
- //-
- inline tensor& Q();
-
- //-
- inline const vector& v() const;
-
- //-
- inline vector& v();
-
- //-
- inline const vector& a() const;
-
- //-
- inline vector& a();
-
- //-
- inline const vector& pi() const;
-
- //-
- inline vector& pi();
-
- //-
- inline const vector& tau() const;
-
- //-
- inline vector& tau();
-
- //-
- inline const List<vector>& siteForces() const;
-
- //-
- inline List<vector>& siteForces();
-
- //-
- inline const List<vector>& sitePositions() const;
-
- //-
- inline List<vector>& sitePositions();
-
- //-
- inline const vector& specialPosition() const;
-
- //-
- inline vector& specialPosition();
-
- //-
- inline scalar potentialEnergy() const;
-
- //-
- inline scalar& potentialEnergy();
-
- //-
- inline const tensor& rf() const;
-
- //-
- inline tensor& rf();
-
- //-
- inline label special() const;
-
- //-
- inline bool tethered() const;
-
- //-
- inline label id() const;
-
-
- // Member Operators
-
- //- Overridable function to handle the particle hitting a patch
- // Executed before other patch-hitting functions
- bool hitPatch
- (
- const polyPatch&,
- trackingData& td,
- const label patchI,
- const scalar trackFraction,
- const tetIndices& tetIs
- );
-
- //- Overridable function to handle the particle hitting a processorPatch
- void hitProcessorPatch
- (
- const processorPolyPatch&,
- trackingData& td
- );
-
- //- Overridable function to handle the particle hitting a wallPatch
- void hitWallPatch
- (
- const wallPolyPatch&,
- trackingData& td,
- const tetIndices&
- );
-
- //- Overridable function to handle the particle hitting a polyPatch
- void hitPatch
- (
- const polyPatch&,
- trackingData& td
- );
-
-
- // I-O
-
- //- Read
- static void readFields(Cloud<polyatomic>& mC);
-
- //- Write
- static void writeFields(const Cloud<polyatomic>& mC);
-
- //- Show info
- static void info(trackingData& td);
-
-
- // IOstream Operators
-
- friend Ostream& operator<<(Ostream&, const polyatomic&);
-};
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#include "polyatomicI.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/molecules/polyatomic/polyatomicI.H b/src/lagrangian/molecularDynamics/molecule/molecules/polyatomic/polyatomicI.H
deleted file mode 100644
index 57149e5fa9b..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/molecules/polyatomic/polyatomicI.H
+++ /dev/null
@@ -1,653 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2008-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
-
-inline Foam::polyatomic::constantProperties::constantProperties()
-:
- sites_(),
- electrostaticSites_(),
- pairPotSites_(),
- momentOfInertia_(diagTensor(0, 0, 0)),
- mass_(0)
-{}
-
-
-inline Foam::polyatomic::constantProperties::constantProperties
-(
- const dictionary& dict,
- const labelList& siteIds
-)
-:
- sites_(),
- electrostaticSites_(),
- pairPotSites_(),
- momentOfInertia_(),
- mass_(0)
-{
- List<word> siteIdNames = dict.lookup("siteIds");
-
- List<vector> siteReferencePositions(dict.lookup("siteReferencePositions"));
-
- List<scalar> siteMasses(dict.lookup("siteMasses"));
-
- List<scalar> siteCharges(dict.lookup("siteCharges"));
-
- List<word> pairPotentialIds(dict.lookup("pairPotentialSiteIds"));
-
- sites_.setSize(siteReferencePositions.size());
-
- if
- (
- (siteIdNames.size() != sites_.size())
- || (siteCharges.size() != sites_.size())
- )
- {
- FatalErrorIn
- (
- "Foam::polyatomic::constantProperties::constantProperties"
- "("
- "const dictionary& dict"
- ")"
- )
- << "Sizes of site id, charge and "
- << "referencePositions are not the same: " << nl
- << siteIdNames << nl
- << siteReferencePositions << nl
- << siteCharges << nl
- << abort(FatalError);
- }
-
- electrostaticSites_.setSize(sites_.size(), -1);
- pairPotSites_.setSize(sites_.size(), -1);
-
- label pairI = 0;
- label electrostaticI = 0;
-
- forAll(sites_, sI)
- {
- sites_[sI] = constPropSite
- (
- siteReferencePositions[sI],
- siteMasses[sI],
- siteCharges[sI],
- siteIds[sI],
- siteIdNames[sI],
- (findIndex(pairPotentialIds, siteIdNames[sI]) != -1), // pair
- (mag(siteCharges[sI]) > VSMALL) // charge
- );
-
- if (sites_[sI].pairPotentialSite())
- {
- pairPotSites_[pairI++] = sI;
- }
-
- if (sites_[sI].electrostaticSite())
- {
- electrostaticSites_[electrostaticI++] = sI;
- }
-
- if (!sites_[sI].pairPotentialSite() && !sites_[sI].electrostaticSite())
- {
- WarningIn
- (
- "Foam::polyatomic::constantProperties::constantProperties"
- "("
- "const dictionary& dict"
- ")"
- )
- << siteIdNames[sI] << " is a non-interacting site." << endl;
- }
- }
-
- pairPotSites_.setSize(pairI);
- electrostaticSites_.setSize(electrostaticI);
-
- // Calculate the centre of mass of the body and subtract it from each
- // position
-
- vector centreOfMass(vector::zero);
-
- forAll(sites_, sI)
- {
- mass_ += sites_[sI].siteMass();
-
- centreOfMass +=
- sites_[sI].siteReferencePosition()*sites_[sI].siteMass();
- }
-
- centreOfMass /= mass_;
-
- forAll(sites_, sI)
- {
- sites_[sI].siteReferencePosition() -= centreOfMass;
- }
-
- if (sites_.size() == 1)
- {
- // Single site polyatomic - no rotational motion.
-
- sites_[0].siteReferencePosition() = vector::zero;
-
- momentOfInertia_ = diagTensor(-1, -1, -1);
- }
- else if (linearMoleculeTest())
- {
- // Linear polyatomic.
-
- Info<< nl << "Linear polyatomic." << endl;
-
- vector dir =
- sites_[1].siteReferencePosition()
- - sites_[0].siteReferencePosition();
-
- dir /= mag(dir);
-
- tensor Q = rotationTensor(dir, vector(1,0,0));
-
- forAll(sites_, sI)
- {
- sites_[sI].siteReferencePosition() =
- (Q & sites_[sI].siteReferencePosition());
- }
-
- // The rotation was around the centre of mass but remove any
- // components that have crept in due to floating point errors
-
- centreOfMass = vector::zero;
-
- forAll(sites_, sI)
- {
- centreOfMass +=
- sites_[sI].siteReferencePosition()*sites_[sI].siteMass();
- }
-
- centreOfMass /= mass_;
-
- forAll(sites_, sI)
- {
- sites_[sI].siteReferencePosition() -= centreOfMass;
- }
-
- diagTensor momOfInertia = diagTensor::zero;
-
- forAll(sites_, sI)
- {
- const vector& p(sites_[sI].siteReferencePosition());
-
- momOfInertia +=
- sites_[sI].siteMass()*diagTensor(0, p.x()*p.x(), p.x()*p.x());
- }
-
- momentOfInertia_ = diagTensor
- (
- -1,
- momOfInertia.yy(),
- momOfInertia.zz()
- );
- }
- else
- {
- // Fully 6DOF polyatomic
-
- // Calculate the inertia tensor in the current orientation
-
- tensor momOfInertia(tensor::zero);
-
- forAll(sites_, sI)
- {
- const vector& p(sites_[sI].siteReferencePosition());
-
- momOfInertia += sites_[sI].siteMass()*tensor
- (
- p.y()*p.y() + p.z()*p.z(), -p.x()*p.y(), -p.x()*p.z(),
- -p.y()*p.x(), p.x()*p.x() + p.z()*p.z(), -p.y()*p.z(),
- -p.z()*p.x(), -p.z()*p.y(), p.x()*p.x() + p.y()*p.y()
- );
- }
-
- if (eigenValues(momOfInertia).x() < VSMALL)
- {
- FatalErrorIn("polyatomic::constantProperties::constantProperties")
- << "An eigenvalue of the inertia tensor is zero, the "
- << "polyatomic " << dict.name() << " is not a valid 6DOF shape."
- << nl << abort(FatalError);
- }
-
- // Normalise the inertia tensor magnitude to avoid SMALL numbers in the
- // components causing problems
-
- momOfInertia /= eigenValues(momOfInertia).x();
-
- tensor e = eigenVectors(momOfInertia);
-
- // Calculate the transformation between the principle axes to the
- // global axes
-
- tensor Q =
- vector(1,0,0)*e.x() + vector(0,1,0)*e.y() + vector(0,0,1)*e.z();
-
- // Transform the site positions
-
- forAll(sites_, sI)
- {
- sites_[sI].siteReferencePosition() =
- (Q & sites_[sI].siteReferencePosition());
- }
-
- // Recalculating the moment of inertia with the new site positions
-
- // The rotation was around the centre of mass but remove any
- // components that have crept in due to floating point errors
-
- centreOfMass = vector::zero;
-
- forAll(sites_, sI)
- {
- centreOfMass +=
- sites_[sI].siteReferencePosition()*sites_[sI].siteMass();
- }
-
- centreOfMass /= mass_;
-
- forAll(sites_, sI)
- {
- sites_[sI].siteReferencePosition() -= centreOfMass;
- }
-
- // Calculate the moment of inertia in the principle component
- // reference frame
-
- momOfInertia = tensor::zero;
-
- forAll(sites_, sI)
- {
- const vector& p(sites_[sI].siteReferencePosition());
-
- momOfInertia += sites_[sI].siteMass()*tensor
- (
- p.y()*p.y() + p.z()*p.z(), -p.x()*p.y(), -p.x()*p.z(),
- -p.y()*p.x(), p.x()*p.x() + p.z()*p.z(), -p.y()*p.z(),
- -p.z()*p.x(), -p.z()*p.y(), p.x()*p.x() + p.y()*p.y()
- );
- }
-
- momentOfInertia_ = diagTensor
- (
- momOfInertia.xx(),
- momOfInertia.yy(),
- momOfInertia.zz()
- );
- }
-}
-
-
-inline Foam::polyatomic::polyatomic
-(
- const polyMesh& mesh,
- const vector& position,
- const label cellI,
- const label tetFaceI,
- const label tetPtI,
- const scalar temperature,
- const vector& bulkVelocity,
- const vector& specialPosition,
- const constantProperties& cP,
- polyatomic::trackingData& td,
- const label special,
- const label id
-)
-:
- particle(mesh, position, cellI, tetFaceI, tetPtI),
- Q_(I),
- v_(bulkVelocity),
- a_(vector::zero),
- pi_(vector::zero),
- tau_(vector::zero),
- specialPosition_(specialPosition),
- potentialEnergy_(0.0),
- rf_(tensor::zero),
- special_(special),
- id_(id),
- siteForces_(cP.nSites(), vector::zero),
- sitePositions_(cP.nSites())
-{
- setSitePositions(cP);
-
- v_ += td.cloud().equipartitionLinearVelocity(temperature, cP.mass());
-
- if (!cP.pointMolecule())
- {
- pi_ = td.cloud().equipartitionAngularMomentum(temperature, cP);
-
- scalar phi(td.cloud().rndGen().scalar01()*twoPi);
-
- scalar theta(td.cloud().rndGen().scalar01()*twoPi);
-
- scalar psi(td.cloud().rndGen().scalar01()*twoPi);
-
- Q_ = tensor
- (
- cos(psi)*cos(phi) - cos(theta)*sin(phi)*sin(psi),
- cos(psi)*sin(phi) + cos(theta)*cos(phi)*sin(psi),
- sin(psi)*sin(theta),
- - sin(psi)*cos(phi) - cos(theta)*sin(phi)*cos(psi),
- - sin(psi)*sin(phi) + cos(theta)*cos(phi)*cos(psi),
- cos(psi)*sin(theta),
- sin(theta)*sin(phi),
- - sin(theta)*cos(phi),
- cos(theta)
- );
- }
-}
-
-
-Foam::polyatomic::polyatomic
-(
- const polyMesh& mesh,
- const vector& position,
- const label cellI,
- const label tetFaceI,
- const label tetPtI,
- const tensor& Q,
- const vector& v,
- const vector& a,
- const vector& pi,
- const vector& tau,
- const vector& specialPosition,
- const constantProperties& cP,
- const label special,
- const label id
-)
-:
- particle(mesh, position, cellI, tetFaceI, tetPtI),
- Q_(Q),
- v_(v),
- a_(a),
- pi_(pi),
- tau_(tau),
- specialPosition_(specialPosition),
- potentialEnergy_(0.0),
- rf_(tensor::zero),
- special_(special),
- id_(id),
- siteForces_(cP.nSites(), vector::zero),
- sitePositions_(cP.nSites())
-{
- setSitePositions(cP);
-}
-
-
-// * * * * * * * constantProperties Private Member Functions * * * * * * * * //
-
-inline bool Foam::polyatomic::constantProperties::linearMoleculeTest() const
-{
- if (sites_.size() == 2)
- {
- return true;
- }
-
- vector refDir =
- sites_[1].siteReferencePosition()
- - sites_[0].siteReferencePosition();
-
- refDir /= mag(refDir);
-
- for
- (
- label i = 2;
- i < sites_.size();
- i++
- )
- {
- vector dir =
- sites_[i].siteReferencePosition()
- - sites_[i - 1].siteReferencePosition();
-
- dir /= mag(dir);
-
- if (mag(refDir & dir) < 1 - SMALL)
- {
- return false;
- }
- }
-
- return true;
-}
-
-
-// * * * * * * * constantProperties Member Functions * * * * * * * * * * * * //
-
-inline const Foam::List<Foam::constPropSite>&
-Foam::polyatomic::constantProperties::sites() const
-{
- return sites_;
-}
-
-
-inline const Foam::List<Foam::label>&
-Foam::polyatomic::constantProperties::pairPotSites() const
-{
- return pairPotSites_;
-}
-
-
-inline const Foam::List<Foam::label>&
-Foam::polyatomic::constantProperties::electrostaticSites() const
-{
- return electrostaticSites_;
-}
-
-
-inline const Foam::diagTensor&
-Foam::polyatomic::constantProperties::momentOfInertia() const
-{
- return momentOfInertia_;
-}
-
-
-inline bool Foam::polyatomic::constantProperties::linearMolecule() const
-{
- return ((momentOfInertia_.xx() < 0) && (momentOfInertia_.yy() > 0));
-}
-
-
-inline bool Foam::polyatomic::constantProperties::pointMolecule() const
-{
- return (momentOfInertia_.zz() < 0);
-}
-
-
-inline Foam::label
-Foam::polyatomic::constantProperties::degreesOfFreedom() const
-{
- if (linearMolecule())
- {
- return 5;
- }
- else if (pointMolecule())
- {
- return 3;
- }
- else
- {
- return 6;
- }
-}
-
-
-inline Foam::scalar Foam::polyatomic::constantProperties::mass() const
-{
- return mass_;
-}
-
-
-inline Foam::label Foam::polyatomic::constantProperties::nSites() const
-{
- return sites_.size();
-}
-
-
-// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
-
-inline const Foam::tensor& Foam::polyatomic::Q() const
-{
- return Q_;
-}
-
-
-inline Foam::tensor& Foam::polyatomic::Q()
-{
- return Q_;
-}
-
-
-inline const Foam::vector& Foam::polyatomic::v() const
-{
- return v_;
-}
-
-
-inline Foam::vector& Foam::polyatomic::v()
-{
- return v_;
-}
-
-
-inline const Foam::vector& Foam::polyatomic::a() const
-{
- return a_;
-}
-
-
-inline Foam::vector& Foam::polyatomic::a()
-{
- return a_;
-}
-
-
-inline const Foam::vector& Foam::polyatomic::pi() const
-{
- return pi_;
-}
-
-
-inline Foam::vector& Foam::polyatomic::pi()
-{
- return pi_;
-}
-
-
-inline const Foam::vector& Foam::polyatomic::tau() const
-{
- return tau_;
-}
-
-
-inline Foam::vector& Foam::polyatomic::tau()
-{
- return tau_;
-}
-
-
-inline const Foam::List<Foam::vector>& Foam::polyatomic::siteForces() const
-{
- return siteForces_;
-}
-
-
-inline Foam::List<Foam::vector>& Foam::polyatomic::siteForces()
-{
- return siteForces_;
-}
-
-
-inline const Foam::List<Foam::vector>& Foam::polyatomic::sitePositions() const
-{
- return sitePositions_;
-}
-
-
-inline Foam::List<Foam::vector>& Foam::polyatomic::sitePositions()
-{
- return sitePositions_;
-}
-
-
-inline const Foam::vector& Foam::polyatomic::specialPosition() const
-{
- return specialPosition_;
-}
-
-
-inline Foam::vector& Foam::polyatomic::specialPosition()
-{
- return specialPosition_;
-}
-
-
-inline Foam::scalar Foam::polyatomic::potentialEnergy() const
-{
- return potentialEnergy_;
-}
-
-
-inline Foam::scalar& Foam::polyatomic::potentialEnergy()
-{
- return potentialEnergy_;
-}
-
-
-inline const Foam::tensor& Foam::polyatomic::rf() const
-{
- return rf_;
-}
-
-
-inline Foam::tensor& Foam::polyatomic::rf()
-{
- return rf_;
-}
-
-
-inline Foam::label Foam::polyatomic::special() const
-{
- return special_;
-}
-
-
-inline bool Foam::polyatomic::tethered() const
-{
- return special_ == SPECIAL_TETHERED;
-}
-
-
-inline Foam::label Foam::polyatomic::id() const
-{
- return id_;
-}
-
-
-// ************************************************************************* //
diff --git a/src/lagrangian/molecularDynamics/molecule/molecules/polyatomic/polyatomicIO.C b/src/lagrangian/molecularDynamics/molecule/molecules/polyatomic/polyatomicIO.C
deleted file mode 100644
index 86dcbfe49f9..00000000000
--- a/src/lagrangian/molecularDynamics/molecule/molecules/polyatomic/polyatomicIO.C
+++ /dev/null
@@ -1,414 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2008-2011 OpenCFD Ltd.
- \\/ M anipulation |
--------------------------------------------------------------------------------
-License
- This file is part of OpenFOAM.
-
- OpenFOAM is free software: you can redistribute it and/or modify it
- under the terms of the GNU General Public License as published by
- the Free Software Foundation, either version 3 of the License, or
- (at your option) any later version.
-
- OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
- FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
- for more details.
-
- You should have received a copy of the GNU General Public License
- along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "polyatomic.H"
-#include "IOstreams.H"
-#include "Cloud.H"
-
-// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
-
-Foam::polyatomic::polyatomic
-(
- const polyMesh& mesh,
- Istream& is,
- bool readFields
-)
-:
- particle(mesh, is, readFields),
- Q_(tensor::zero),
- v_(vector::zero),
- a_(vector::zero),
- pi_(vector::zero),
- tau_(vector::zero),
- specialPosition_(vector::zero),
- potentialEnergy_(0.0),
- rf_(tensor::zero),
- special_(0),
- id_(0),
- siteForces_(0),
- sitePositions_(0)
-{
- if (readFields)
- {
- if (is.format() == IOstream::ASCII)
- {
- is >> Q_;
- is >> v_;
- is >> a_;
- is >> pi_;
- is >> tau_;
- is >> siteForces_;
- potentialEnergy_ = readScalar(is);
- is >> rf_;
- special_ = readLabel(is);
- id_ = readLabel(is);
- is >> sitePositions_;
- is >> specialPosition_;
- }
- else
- {
- is.read
- (
- reinterpret_cast<char*>(&Q_),
- sizeof(Q_)
- + sizeof(v_)
- + sizeof(a_)
- + sizeof(pi_)
- + sizeof(tau_)
- + sizeof(specialPosition_)
- + sizeof(potentialEnergy_)
- + sizeof(rf_)
- + sizeof(special_)
- + sizeof(id_)
- );
-
- is >> siteForces_ >> sitePositions_;
- }
- }
-
- // Check state of Istream
- is.check
- (
- "Foam::polyatomic::polyatomic"
- "("
- "const polyMesh& mesh,"
- "Istream& is,"
- "bool readFields"
- ")"
- );
-}
-
-
-void Foam::polyatomic::readFields(Cloud<polyatomic>& mC)
-{
- if (!mC.size())
- {
- return;
- }
-
- particle::readFields(mC);
-
- IOField<tensor> Q(mC.fieldIOobject("Q", IOobject::MUST_READ));
- mC.checkFieldIOobject(mC, Q);
-
- IOField<vector> v(mC.fieldIOobject("v", IOobject::MUST_READ));
- mC.checkFieldIOobject(mC, v);
-
- IOField<vector> a(mC.fieldIOobject("a", IOobject::MUST_READ));
- mC.checkFieldIOobject(mC, a);
-
- IOField<vector> pi(mC.fieldIOobject("pi", IOobject::MUST_READ));
- mC.checkFieldIOobject(mC, pi);
-
- IOField<vector> tau(mC.fieldIOobject("tau", IOobject::MUST_READ));
- mC.checkFieldIOobject(mC, tau);
-
- IOField<vector> specialPosition
- (
- mC.fieldIOobject("specialPosition", IOobject::MUST_READ)
- );
- mC.checkFieldIOobject(mC, specialPosition);
-
- IOField<label> special(mC.fieldIOobject("special", IOobject::MUST_READ));
- mC.checkFieldIOobject(mC, special);
-
- IOField<label> id(mC.fieldIOobject("id", IOobject::MUST_READ));
- mC.checkFieldIOobject(mC, id);
-
- label i = 0;
-
- forAllIter(typename Cloud<polyatomic>, mC, iter)
- {
- polyatomic& mol = iter();
-
- mol.Q_ = Q[i];
- mol.v_ = v[i];
- mol.a_ = a[i];
- mol.pi_ = pi[i];
- mol.tau_ = tau[i];
- mol.specialPosition_ = specialPosition[i];
- mol.special_ = special[i];
- mol.id_ = id[i];
- i++;
- }
-}
-
-
-void Foam::polyatomic::writeFields(const Cloud<polyatomic>& mC)
-{
- particle::writeFields(mC);
-
- label np = mC.size();
-
- IOField<tensor> Q(mC.fieldIOobject("Q", IOobject::NO_READ), np);
- IOField<vector> v(mC.fieldIOobject("v", IOobject::NO_READ), np);
- IOField<vector> a(mC.fieldIOobject("a", IOobject::NO_READ), np);
- IOField<vector> pi(mC.fieldIOobject("pi", IOobject::NO_READ), np);
- IOField<vector> tau(mC.fieldIOobject("tau", IOobject::NO_READ), np);
- IOField<vector> specialPosition
- (
- mC.fieldIOobject("specialPosition", IOobject::NO_READ),
- np
- );
- IOField<label> special(mC.fieldIOobject("special", IOobject::NO_READ), np);
- IOField<label> id(mC.fieldIOobject("id", IOobject::NO_READ), np);
-
- // Post processing fields
-
- IOField<vector> piGlobal
- (
- mC.fieldIOobject("piGlobal", IOobject::NO_READ),
- np
- );
-
- IOField<vector> tauGlobal
- (
- mC.fieldIOobject("tauGlobal", IOobject::NO_READ),
- np
- );
-
- IOField<vector> orientation1
- (
- mC.fieldIOobject("orientation1", IOobject::NO_READ),
- np
- );
-
- IOField<vector> orientation2
- (
- mC.fieldIOobject("orientation2", IOobject::NO_READ),
- np
- );
-
- IOField<vector> orientation3
- (
- mC.fieldIOobject("orientation3", IOobject::NO_READ),
- np
- );
-
- label i = 0;
- forAllConstIter(typename Cloud<polyatomic>, mC, iter)
- {
- const polyatomic& mol = iter();
-
- Q[i] = mol.Q_;
- v[i] = mol.v_;
- a[i] = mol.a_;
- pi[i] = mol.pi_;
- tau[i] = mol.tau_;
- specialPosition[i] = mol.specialPosition_;
- special[i] = mol.special_;
- id[i] = mol.id_;
-
- piGlobal[i] = mol.Q_ & mol.pi_;
- tauGlobal[i] = mol.Q_ & mol.tau_;
-
- orientation1[i] = mol.Q_ & vector(1,0,0);
- orientation2[i] = mol.Q_ & vector(0,1,0);
- orientation3[i] = mol.Q_ & vector(0,0,1);
-
- i++;
- }
-
- Q.write();
- v.write();
- a.write();
- pi.write();
- tau.write();
- specialPosition.write();
- special.write();
- id.write();
-
- piGlobal.write();
- tauGlobal.write();
-
- orientation1.write();
- orientation2.write();
- orientation3.write();
-}
-
-
-void Foam::polyatomic::info(polyatomic::trackingData& td)
-{
- vector totalLinearMomentum(vector::zero);
- vector totalAngularMomentum(vector::zero);
- scalar maxVelocityMag = 0.0;
- scalar totalMass = 0.0;
- scalar totalLinearKE = 0.0;
- scalar totalAngularKE = 0.0;
- scalar totalPE = 0.0;
- scalar totalrDotf = 0.0;
- //vector CentreOfMass(vector::zero);
- label nMols = td.cloud().size();
- label dofs = 0;
-
- forAllConstIter(typename Cloud<polyatomic>, td.cloud(), mol)
- {
- const label molId = mol().id();
- scalar molMass(td.cloud().constProps(molId).mass());
- totalMass += molMass;
- //CentreOfMass += mol().position()*molMass;
- }
-
- // if (nMols)
- // {
- // CentreOfMass /= totalMass;
- // }
-
- forAllConstIter(typename Cloud<polyatomic>, td.cloud(), mol)
- {
- const label molId = mol().id();
- const polyatomic::constantProperties cP
- (
- td.cloud().constProps(molId)
- );
- scalar molMass(cP.mass());
- const diagTensor& molMoI(cP.momentOfInertia());
- const vector& molV(mol().v());
- const vector& molOmega(inv(molMoI) & mol().pi());
- vector molPiGlobal = mol().Q() & mol().pi();
- totalLinearMomentum += molV * molMass;
- totalAngularMomentum += molPiGlobal;
- //+((mol().position() - CentreOfMass) ^ (molV * molMass));
- if (mag(molV) > maxVelocityMag)
- {
- maxVelocityMag = mag(molV);
- }
- totalLinearKE += 0.5*molMass*magSqr(molV);
- totalAngularKE += 0.5*(molOmega & molMoI & molOmega);
- totalPE += mol().potentialEnergy();
- totalrDotf += tr(mol().rf());
- dofs += cP.degreesOfFreedom();
- }
-
- scalar meshVolume = sum(td.cloud().mesh().cellVolumes());
-
- if (Pstream::parRun())
- {
- reduce(totalLinearMomentum, sumOp<vector>());
- reduce(totalAngularMomentum, sumOp<vector>());
- reduce(maxVelocityMag, maxOp<scalar>());
- reduce(totalMass, sumOp<scalar>());
- reduce(totalLinearKE, sumOp<scalar>());
- reduce(totalAngularKE, sumOp<scalar>());
- reduce(totalPE, sumOp<scalar>());
- reduce(totalrDotf, sumOp<scalar>());
- reduce(nMols, sumOp<label>());
- reduce(dofs, sumOp<label>());
- reduce(meshVolume, sumOp<scalar>());
- }
-
- if (nMols)
- {
- Info<< nl << "Number of molecules in " << td.cloud().name() << " = "
- << nMols << nl
- << " Overall number density = "
- << nMols/meshVolume << nl
- << " Overall mass density = "
- << totalMass/meshVolume << nl
- << " Average linear momentum per molecule = "
- << totalLinearMomentum/nMols << ' '
- << mag(totalLinearMomentum)/nMols << nl
- << " Average angular momentum per molecule = "
- << totalAngularMomentum << ' '
- << mag(totalAngularMomentum)/nMols << nl
- << " maximum |velocity| = "
- << maxVelocityMag << nl
- << " Average linear KE per molecule = "
- << totalLinearKE/nMols << nl
- << " Average angular KE per molecule = "
- << totalAngularKE/nMols << nl
- << " Average PE per molecule = "
- << totalPE/nMols << nl
- << " Average TE per molecule = "
- <<
- (
- totalLinearKE
- + totalAngularKE
- + totalPE
- )
- /nMols
- << nl << endl;
- }
- else
- {
- Info<< nl << "No molecules in " << td.cloud().name() << endl;
- }
-}
-
-
-// * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
-
-Foam::Ostream& Foam::operator<<(Ostream& os, const polyatomic& mol)
-{
- if (os.format() == IOstream::ASCII)
- {
- os << token::SPACE << static_cast<const particle&>(mol)
- << token::SPACE << mol.face()
- << token::SPACE << mol.stepFraction()
- << token::SPACE << mol.Q_
- << token::SPACE << mol.v_
- << token::SPACE << mol.a_
- << token::SPACE << mol.pi_
- << token::SPACE << mol.tau_
- << token::SPACE << mol.specialPosition_
- << token::SPACE << mol.potentialEnergy_
- << token::SPACE << mol.rf_
- << token::SPACE << mol.special_
- << token::SPACE << mol.id_
- << token::SPACE << mol.siteForces_
- << token::SPACE << mol.sitePositions_;
- }
- else
- {
- os << static_cast<const particle&>(mol);
- os.write
- (
- reinterpret_cast<const char*>(&mol.Q_),
- sizeof(mol.Q_)
- + sizeof(mol.v_)
- + sizeof(mol.a_)
- + sizeof(mol.pi_)
- + sizeof(mol.tau_)
- + sizeof(mol.specialPosition_)
- + sizeof(mol.potentialEnergy_)
- + sizeof(mol.rf_)
- + sizeof(mol.special_)
- + sizeof(mol.id_)
- );
- os << mol.siteForces_ << mol.sitePositions_;
- }
-
- // Check state of Ostream
- os.check
- (
- "Foam::Ostream& Foam::operator<<"
- "(Foam::Ostream&, const Foam::polyatomic&)"
- );
-
- return os;
-}
-
-
-// ************************************************************************* //
--
GitLab