zhyang-dev/openfoam
1
0
Fork 0
mirror of https://develop.openfoam.com/Development/openfoam.git synced 2025-11-28 03:28:01 +00:00
Code Issues Actions Packages Projects Releases Wiki Activity

CrdirectInteractionList building switchable using PP or PFEE. New solver to test rewritten moleculeCloud etc

Browse Source
This commit is contained in:
graham
2008-09-11 19:04:37 +01:00
parent 8dbfc59c5e
commit ec6d172446
12 changed files with 479 additions and 94 deletions
Download Patch File Download Diff File Expand all files Collapse all files

2
src/lagrangian/molecularDynamics/molecule/Make/files
View File

@ -8,7 +8,7 @@ directInteractionList = $(interactionLists)/directInteractionList
$(referralLists)/receivingReferralList.C
$(referralLists)/sendingReferralList.C
// $(referredCellList)/referredCellList.C
// $(referredCell)/referredCell.C
$(referredCell)/referredCell.C
$(referredMolecule)/referredMolecule.C
$(directInteractionList)/directInteractionList.C
$(interactionLists)/interactionLists.C

81
src/lagrangian/molecularDynamics/molecule/interactionLists/directInteractionList/directInteractionList.C
View File

@ -25,37 +25,29 @@ License
\*---------------------------------------------------------------------------*/
#include "directInteractionList.H"
#include "interactionLists.H"
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
void Foam::directInteractionList::buildDirectInteractionList
(
scalar rCutMax,
bool pointPointListBuild
)
{
Info<< nl << "Building list of direct interaction neighbours" << endl;
scalar rCutMaxSqr = rCutMax*rCutMax;
const polyMesh& mesh(interactionLists_.mesh());
const polyMesh& mesh(il_.mesh());
List<DynamicList<label> > directInteractionList(mesh.nCells());
if (pointPointListBuild)
{
Info<< "Point-Point direct interaction list build." << endl;
Info<< tab << "Point-Point direct interaction list build." << endl;
label pointJIndex;
forAll (mesh.points(), pointIIndex)
{
const point& ptI
(
mesh.points()[pointIIndex]
);
for
(
pointJIndex = pointIIndex;
@ -63,12 +55,7 @@ void Foam::directInteractionList::buildDirectInteractionList
++pointJIndex
)
{
const point& ptJ
(
mesh.points()[pointJIndex]
);
if (magSqr(ptI - ptJ) <= rCutMaxSqr)
if (il_.testPointPointDistance(pointIIndex, pointJIndex))
{
const labelList& ptICells
(
@ -111,13 +98,13 @@ void Foam::directInteractionList::buildDirectInteractionList
}
else
{
Info<< "Point-Face, Edge-Edge direct interaction list build." << endl;
Info<< tab << "Point-Face, Edge-Edge direct interaction list build." << endl;
forAll (mesh.points(), p)
{
forAll(mesh.faces(), f)
{
if(testPointFaceDistance(p, f))
if(il_.testPointFaceDistance(p, f))
{
const labelList& pCells(mesh.pointCells()[p]);
@ -153,11 +140,7 @@ void Foam::directInteractionList::buildDirectInteractionList
if (cellN > cellI)
{
if
(
findIndex(directInteractionList[cellI], cellN)
== -1
)
if (findIndex(directInteractionList[cellI], cellN) == -1)
{
directInteractionList[cellI].append(cellN);
}
@ -165,11 +148,7 @@ void Foam::directInteractionList::buildDirectInteractionList
if (cellI > cellN)
{
if
(
findIndex(directInteractionList[cellN], cellI)
== -1
)
if (findIndex(directInteractionList[cellN], cellI) == -1)
{
directInteractionList[cellN].append(cellI);
}
@ -195,7 +174,7 @@ void Foam::directInteractionList::buildDirectInteractionList
{
const edge& eJ(mesh.edges()[edgeJIndex]);
if (testEdgeEdgeDistance(eI, eJ))
if (il_.testEdgeEdgeDistance(eI, eJ))
{
const labelList& eICells(mesh.edgeCells()[edgeIIndex]);
@ -211,11 +190,7 @@ void Foam::directInteractionList::buildDirectInteractionList
if (cellJ > cellI)
{
if
(
findIndex(directInteractionList[cellI], cellJ)
== -1
)
if (findIndex(directInteractionList[cellI], cellJ) == -1)
{
directInteractionList[cellI].append(cellJ);
}
@ -223,11 +198,7 @@ void Foam::directInteractionList::buildDirectInteractionList
if (cellI > cellJ)
{
if
(
findIndex(directInteractionList[cellJ], cellI)
== -1
)
if (findIndex(directInteractionList[cellJ], cellI) == -1)
{
directInteractionList[cellJ].append(cellI);
}
@ -246,6 +217,13 @@ void Foam::directInteractionList::buildDirectInteractionList
directInteractionList[transDIL].shrink()
);
}
// sorting DILs
forAll((*this), dIL)
{
sort((*this)[dIL]);
}
}
@ -253,38 +231,29 @@ void Foam::directInteractionList::buildDirectInteractionList
Foam::directInteractionList::directInteractionList
(
const interactionLists& interactionLists,
scalar rCutMax,
const interactionLists& il,
bool pointPointListBuild
)
:
labelListList(interactionLists.mesh().nCells()),
interactionLists_(interactionLists)
labelListList(il.mesh().nCells()),
il_(il)
{
buildDirectInteractionList(rCutMax, pointPointListBuild);
buildDirectInteractionList(pointPointListBuild);
}
Foam::directInteractionList::directInteractionList
(
const interactionLists& interactionLists
const interactionLists& il
)
:
labelListList(interactionLists.mesh().nCells()),
interactionLists_(interactionLists)
labelListList(il.mesh().nCells()),
il_(il)
{
Info<< "Read directInteractionList from disk not implemented" << endl;
}
// * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
Foam::autoPtr<Foam::directInteractionList> Foam::directInteractionList::New()
{
return autoPtr<directInteractionList>(new directInteractionList);
}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::directInteractionList::~directInteractionList()

8
src/lagrangian/molecularDynamics/molecule/interactionLists/directInteractionList/directInteractionList.H
View File

@ -57,13 +57,12 @@ class directInteractionList
{
// Private data
const interactionLists& interactionLists_;
const interactionLists& il_;
// Private Member Functions
void buildDirectInteractionList
(
scalar rCutMax,
bool pointPointListBuild
);
@ -81,8 +80,7 @@ public:
directInteractionList
(
const interactionLists& interactionLists,
scalar rCutMax,
bool pointPointListBuild = false
bool pointPointListBuild
);
//- Construct from file
@ -100,7 +98,7 @@ public:
// Access
inline const interactionLists& interactionLists() const;
inline const interactionLists& il() const;
// Check

4
src/lagrangian/molecularDynamics/molecule/interactionLists/directInteractionList/directInteractionListI.H
View File

@ -29,9 +29,9 @@ License
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
inline const Foam::interactionLists&
Foam::interactionLists::interactionLists() const
Foam::directInteractionList::il() const
{
return interactionLists_;
return il_;
}

282
src/lagrangian/molecularDynamics/molecule/interactionLists/interactionLists.C
View File

@ -26,20 +26,24 @@ License
#include "interactionLists.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
Foam::scalar Foam::interactionLists::transTol = 1e-12;
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::interactionLists::interactionLists
interactionLists
(
const polyMesh& mesh,
scalar rCutMax,
scalar rCutMaxSqr,
bool pointPointListBuild
)
:
mesh_(mesh),
dil_(*this, rCutMax, pointPointListBuild)
rCutMaxSqr_(rCutMaxSqr),
dil_(*this, pointPointListBuild)
{}
@ -58,5 +62,277 @@ Foam::interactionLists::~interactionLists()
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
bool Foam::interactionLists::testPointPointDistance
(
const label ptI,
const label ptJ
) const
{
return (magSqr(mesh_.points()[ptI] - mesh_.points()[ptJ]) <= rCutMaxSqr_);
}
bool Foam::interactionLists::testEdgeEdgeDistance
(
const edge& eI,
const edge& eJ
) const
{
const vector& eJs(mesh_.points()[eJ.start()]);
const vector& eJe(mesh_.points()[eJ.end()]);
return testEdgeEdgeDistance(eI, eJs, eJe);
}
bool Foam::interactionLists::testPointFaceDistance
(
const label p,
const label faceNo
) const
{
const vector& pointPosition(mesh_.points()[p]);
return testPointFaceDistance(pointPosition, faceNo);
}
bool Foam::interactionLists::testPointFaceDistance
(
const label p,
const referredCell& refCell
) const
{
const vector& pointPosition(mesh_.points()[p]);
forAll (refCell.faces(), rCF)
{
if
(
testPointFaceDistance
(
pointPosition,
refCell.faces()[rCF],
refCell.vertexPositions(),
refCell.faceCentres()[rCF],
refCell.faceAreas()[rCF]
)
)
{
return true;
}
}
return false;
}
bool Foam::interactionLists::testPointFaceDistance
(
const vectorList& pointsToTest,
const label faceNo
) const
{
forAll(pointsToTest, pTT)
{
const vector& p(pointsToTest[pTT]);
// if any point in the list is in range of the face
// then the rest do not need to be tested and
// true can be returned
if (testPointFaceDistance(p, faceNo))
{
return true;
}
}
return false;
}
bool Foam::interactionLists::testPointFaceDistance
(
const vector& p,
const label faceNo
) const
{
const face& faceToTest(mesh_.faces()[faceNo]);
const vector& faceC(mesh_.faceCentres()[faceNo]);
const vector& faceA(mesh_.faceAreas()[faceNo]);
const vectorList& points(mesh_.points());
return testPointFaceDistance
(
p,
faceToTest,
points,
faceC,
faceA
);
}
bool Foam::interactionLists::testPointFaceDistance
(
const vector& p,
const labelList& faceToTest,
const vectorList& points,
const vector& faceC,
const vector& faceA
) const
{
vector faceN(faceA/mag(faceA));
scalar perpDist((p - faceC) & faceN);
if (magSqr(perpDist) > rCutMaxSqr_)
{
return false;
}
vector pointOnPlane = (p - faceN * perpDist);
if (magSqr(faceC - pointOnPlane) < rCutMaxSqr_*1e-8)
{
// If pointOnPlane is very close to the face centre
// then defining the local axes will be inaccurate
// and it is very likely that pointOnPlane will be
// inside the face, so return true if the points
// are in range to be safe
return (magSqr(pointOnPlane - p) <= rCutMaxSqr_);
}
vector xAxis = (faceC - pointOnPlane)/mag(faceC - pointOnPlane);
vector yAxis =
((faceC - pointOnPlane) ^ faceN)
/mag((faceC - pointOnPlane) ^ faceN);
List<vector2D> local2DVertices(faceToTest.size());
forAll(faceToTest, fTT)
{
const vector& V(points[faceToTest[fTT]]);
if (magSqr(V-p) <= rCutMaxSqr_)
{
return true;
}
local2DVertices[fTT] = vector2D
(
((V - pointOnPlane) & xAxis),
((V - pointOnPlane) & yAxis)
);
}
scalar localFaceCx((faceC - pointOnPlane) & xAxis);
scalar la_valid = -1;
forAll(local2DVertices, fV)
{
const vector2D& va(local2DVertices[fV]);
const vector2D& vb
(
local2DVertices[(fV + 1) % local2DVertices.size()]
);
if (mag(vb.y()-va.y()) > SMALL)
{
scalar la =
(
va.x() - va.y()*((vb.x() - va.x())/(vb.y() - va.y()))
)
/localFaceCx;
scalar lv = -va.y()/(vb.y() - va.y());
if (la >= 0 && la <= 1 && lv >= 0 && lv <= 1)
{
la_valid = la;
break;
}
}
}
if (la_valid < 0)
{
// perpendicular point inside face, nearest point is pointOnPlane;
return (magSqr(pointOnPlane-p) <= rCutMaxSqr_);
}
else
{
// perpendicular point outside face, nearest point is
// on edge that generated la_valid;
return
(
magSqr(pointOnPlane + la_valid*(faceC - pointOnPlane) - p)
<= rCutMaxSqr_
);
}
// if the algorithm hasn't returned anything by now then something has
// gone wrong.
FatalErrorIn("interactionLists.C") << nl
<< "point " << p << " to face " << faceToTest
<< " comparison did not find a nearest point"
<< " to be inside or outside face."
<< abort(FatalError);
return false;
}
bool Foam::interactionLists::testEdgeEdgeDistance
(
const edge& eI,
const vector& eJs,
const vector& eJe
) const
{
vector a(eI.vec(mesh_.points()));
vector b(eJe - eJs);
const vector& eIs(mesh_.points()[eI.start()]);
vector c(eJs - eIs);
vector crossab = a ^ b;
scalar magCrossSqr = magSqr(crossab);
if (magCrossSqr > VSMALL)
{
// If the edges are parallel then a point-face
// search will pick them up
scalar s = ((c ^ b) & crossab)/magCrossSqr;
scalar t = ((c ^ a) & crossab)/magCrossSqr;
// Check for end points outside of range 0..1
// If the closest point is outside this range
// a point-face search will have found it.
return
(
s >= 0
&& s <= 1
&& t >= 0
&& t <= 1
&& magSqr(eIs + a*s - eJs - b*t) <= rCutMaxSqr_
);
}
return false;
}
// ************************************************************************* //

79
src/lagrangian/molecularDynamics/molecule/interactionLists/interactionLists.H
View File

@ -39,6 +39,7 @@ SourceFiles
#include "polyMesh.H"
#include "directInteractionList.H"
#include "referredCell.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -55,8 +56,24 @@ class interactionLists
const polyMesh& mesh_;
scalar rCutMaxSqr_;
directInteractionList dil_;
// referredCellList ril_;
// labelList realCellsWithinRCutMaxOfAnyReferringPatch_;
// labelList realFacesWithinRCutMaxOfAnyReferringPatch_;
// labelList realEdgesWithinRCutMaxOfAnyReferringPatch_;
// labelList realPointsWithinRCutMaxOfAnyReferringPatch_;
// List<sendingReferralList> cellSendingReferralLists_;
// List<receivingReferralList> cellReceivingReferralLists_;
// Private Member Functions
//- Disallow default bitwise copy construct
@ -68,13 +85,18 @@ class interactionLists
public:
// Static data members
//- Tolerance for checking that faces on a patch segment
static scalar transTol;
// Constructors
//- Construct and create all information from the mesh
interactionLists
(
const polyMesh& mesh,
scalar rCutMax,
scalar rCutMaxSqr,
bool pointPointListBuild = false
);
@ -88,11 +110,64 @@ public:
// Member Functions
bool testPointPointDistance
(
const label ptI,
const label ptJ
) const;
bool testPointFaceDistance
(
const label p,
const label faceNo
) const;
bool testPointFaceDistance
(
const label p,
const referredCell& refCell
) const;
bool testPointFaceDistance
(
const vectorList& pointsToTest,
const label faceNo
) const;
bool testPointFaceDistance
(
const vector& p,
const label faceNo
) const;
bool testPointFaceDistance
(
const vector& p,
const labelList& faceToTest,
const vectorList& points,
const vector& faceC,
const vector& faceA
) const;
bool testEdgeEdgeDistance
(
const edge& eI,
const edge& eJ
) const;
bool testEdgeEdgeDistance
(
const edge& eI,
const vector& eJs,
const vector& eJe
) const;
// Access
inline const polyMesh& mesh() const;
inline const directInteractionList& dil();
inline const directInteractionList& dil() const;
};

29
src/lagrangian/molecularDynamics/molecule/interactionLists/interactionListsI.H
View File

@ -26,31 +26,18 @@ License
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
// * * * * * * * * * * * * * * * * Destructors * * * * * * * * * * * * * * * //
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
const Foam::polyMesh& Foam::interactionLists::mesh() const
{
return mesh_;
}
// * * * * * * * * * * * * * * * Friend Functions * * * * * * * * * * * * * //
// * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
// * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * * //
// * * * * * * * * * * * * * * * Ostream Operator * * * * * * * * * * * * * //
const Foam::directInteractionList& Foam::interactionLists::dil() const
{
return dil_;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

6
src/lagrangian/molecularDynamics/molecule/interactionLists/referredCell/referredCell.C
View File

@ -25,7 +25,7 @@ License
\*----------------------------------------------------------------------------*/
#include "referredCell.H"
#include "moleculeCloud.H"
#include "interactionLists.H"
namespace Foam
{
@ -378,7 +378,7 @@ bool referredCell::duplicate(const referredCell& refCellDupl) const
(
sourceProc_ == refCellDupl.sourceProc()
&& sourceCell_ == refCellDupl.sourceCell()
&& mag(offset_ - refCellDupl.offset()) < moleculeCloud::transTol
&& mag(offset_ - refCellDupl.offset()) < interactionLists::transTol
);
}
@ -389,7 +389,7 @@ bool referredCell::duplicate(const label procNo,const label nCells) const
(
sourceProc_ == procNo
&& sourceCell_ < nCells
&& mag(offset_) < moleculeCloud::transTol
&& mag(offset_) < interactionLists::transTol
);
}

1
src/lagrangian/molecularDynamics/molecule/mdTools/md.H
View File

@ -2,7 +2,6 @@
#define md_H
# include "moleculeCloud.H"
# include "correlationFunction.H"
# include "distribution.H"
# include "reducedUnits.H"