openfoam/src/lagrangian/basic/particle/particleTemplates.C
Will Bainbridge e9fb8b8572 tetIndices: Removed duplicate logic
The logic for generating tetrahedra from a face base point and an offset
was duplicated in a few places. It is now confined to the tetIndices
class.
2017-06-01 09:59:38 +01:00

502 lines
13 KiB
C

/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2017 OpenFOAM Foundation
\\/ M anipulation | Copyright (C) 2016 OpenCFD Ltd.
-------------------------------------------------------------------------------
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 "IOPosition.H"
#include "cyclicPolyPatch.H"
#include "cyclicAMIPolyPatch.H"
#include "processorPolyPatch.H"
#include "symmetryPlanePolyPatch.H"
#include "symmetryPolyPatch.H"
#include "wallPolyPatch.H"
#include "wedgePolyPatch.H"
#include "meshTools.H"
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
template<class TrackData>
void Foam::particle::prepareForParallelTransfer
(
const label patchi,
TrackData& td
)
{
// Convert the face index to be local to the processor patch
facei_ = mesh_.boundaryMesh()[patchi].whichFace(facei_);
}
template<class TrackData>
void Foam::particle::correctAfterParallelTransfer
(
const label patchi,
TrackData& td
)
{
const coupledPolyPatch& ppp =
refCast<const coupledPolyPatch>(mesh_.boundaryMesh()[patchi]);
if (!ppp.parallel())
{
const tensor& T =
(
ppp.forwardT().size() == 1
? ppp.forwardT()[0]
: ppp.forwardT()[facei_]
);
transformProperties(T);
}
else if (ppp.separated())
{
const vector& s =
(
(ppp.separation().size() == 1)
? ppp.separation()[0]
: ppp.separation()[facei_]
);
transformProperties(-s);
}
// Set the topology
celli_ = ppp.faceCells()[facei_];
facei_ += ppp.start();
tetFacei_ = facei_;
// Faces either side of a coupled patch are numbered in opposite directions
// as their normals both point away from their connected cells. The tet
// point therefore counts in the opposite direction from the base point.
tetPti_ = mesh_.faces()[tetFacei_].size() - 1 - tetPti_;
// Reflect to account for the change of triangle orientation in the new cell
reflect();
// Note that the position does not need transforming explicitly. The face-
// triangle on the receive patch is the transformation of the one on the
// send patch, so whilst the barycentric coordinates remain the same, the
// change of triangle implicitly transforms the position.
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class CloudType>
void Foam::particle::readFields(CloudType& c)
{
if (!c.size())
{
return;
}
IOobject procIO(c.fieldIOobject("origProcId", IOobject::MUST_READ));
if (procIO.typeHeaderOk<IOField<label>>(true))
{
IOField<label> origProcId(procIO);
c.checkFieldIOobject(c, origProcId);
IOField<label> origId(c.fieldIOobject("origId", IOobject::MUST_READ));
c.checkFieldIOobject(c, origId);
label i = 0;
forAllIter(typename CloudType, c, iter)
{
particle& p = iter();
p.origProc_ = origProcId[i];
p.origId_ = origId[i];
i++;
}
}
}
template<class CloudType>
void Foam::particle::writeFields(const CloudType& c)
{
IOPosition<CloudType> ioP(c);
ioP.write();
label np = c.size();
IOField<label> origProc
(
c.fieldIOobject("origProcId", IOobject::NO_READ),
np
);
IOField<label> origId
(
c.fieldIOobject("origId", IOobject::NO_READ),
np
);
label i = 0;
forAllConstIter(typename CloudType, c, iter)
{
origProc[i] = iter().origProc_;
origId[i] = iter().origId_;
i++;
}
origProc.write();
origId.write();
}
template<class CloudType>
void Foam::particle::writeObjects(const CloudType& c, objectRegistry& obr)
{
label np = c.size();
IOField<vector>& position
(
cloud::createIOField<vector>("position", np, obr)
);
IOField<label>& origProc(cloud::createIOField<label>("origProc", np, obr));
IOField<label>& origId(cloud::createIOField<label>("origId", np, obr));
label i = 0;
forAllConstIter(typename CloudType, c, iter)
{
position[i] = iter().position_;
origProc[i] = iter().origProc_;
origId[i] = iter().origId_;
i++;
}
}
template<class TrackData>
void Foam::particle::trackToFace
(
const vector& displacement,
const scalar fraction,
TrackData& td
)
{
// Track
trackToFace(displacement, fraction);
// If the track is complete, return
if (!onFace())
{
return;
}
// Hit face/patch processing
typedef typename TrackData::cloudType::particleType particleType;
particleType& p = static_cast<particleType&>(*this);
p.hitFace(td);
if (onInternalFace())
{
changeCell();
}
else
{
label origFacei = facei_;
label patchi = mesh_.boundaryMesh().whichPatch(facei_);
// No action is taken for tetPti_ for tetFacei_ here. These are handled
// by the patch interaction call or later during processor transfer.
const tetIndices faceHitTetIs(celli_, tetFacei_, tetPti_);
if
(
!p.hitPatch
(
mesh_.boundaryMesh()[patchi],
td,
patchi,
stepFraction(),
faceHitTetIs
)
)
{
// Did patch interaction model switch patches?
if (facei_ != origFacei)
{
patchi = mesh_.boundaryMesh().whichPatch(facei_);
}
const polyPatch& patch = mesh_.boundaryMesh()[patchi];
if (isA<wedgePolyPatch>(patch))
{
p.hitWedgePatch
(
static_cast<const wedgePolyPatch&>(patch), td
);
}
else if (isA<symmetryPlanePolyPatch>(patch))
{
p.hitSymmetryPlanePatch
(
static_cast<const symmetryPlanePolyPatch&>(patch), td
);
}
else if (isA<symmetryPolyPatch>(patch))
{
p.hitSymmetryPatch
(
static_cast<const symmetryPolyPatch&>(patch), td
);
}
else if (isA<cyclicPolyPatch>(patch))
{
p.hitCyclicPatch
(
static_cast<const cyclicPolyPatch&>(patch), td
);
}
else if (isA<cyclicAMIPolyPatch>(patch))
{
p.hitCyclicAMIPatch
(
static_cast<const cyclicAMIPolyPatch&>(patch),
td,
displacement
);
}
else if (isA<processorPolyPatch>(patch))
{
p.hitProcessorPatch
(
static_cast<const processorPolyPatch&>(patch), td
);
}
else if (isA<wallPolyPatch>(patch))
{
p.hitWallPatch
(
static_cast<const wallPolyPatch&>(patch), td, faceHitTetIs
);
}
else
{
p.hitPatch(patch, td);
}
}
}
}
template<class TrackData>
void Foam::particle::hitFace(TrackData&)
{}
template<class TrackData>
bool Foam::particle::hitPatch
(
const polyPatch&,
TrackData&,
const label,
const scalar,
const tetIndices&
)
{
return false;
}
template<class TrackData>
void Foam::particle::hitWedgePatch
(
const wedgePolyPatch& wpp,
TrackData&
)
{
FatalErrorInFunction
<< "Hitting a wedge patch should not be possible."
<< abort(FatalError);
vector nf = normal();
nf /= mag(nf);
transformProperties(I - 2.0*nf*nf);
}
template<class TrackData>
void Foam::particle::hitSymmetryPlanePatch
(
const symmetryPlanePolyPatch& spp,
TrackData&
)
{
vector nf = normal();
nf /= mag(nf);
transformProperties(I - 2.0*nf*nf);
}
template<class TrackData>
void Foam::particle::hitSymmetryPatch
(
const symmetryPolyPatch& spp,
TrackData&
)
{
vector nf = normal();
nf /= mag(nf);
transformProperties(I - 2.0*nf*nf);
}
template<class TrackData>
void Foam::particle::hitCyclicPatch
(
const cyclicPolyPatch& cpp,
TrackData& td
)
{
const cyclicPolyPatch& receiveCpp = cpp.neighbPatch();
const label receiveFacei = receiveCpp.whichFace(facei_);
// Set the topology
facei_ = tetFacei_ = cpp.transformGlobalFace(facei_);
celli_ = mesh_.faceOwner()[facei_];
// See note in correctAfterParallelTransfer for tetPti addressing ...
tetPti_ = mesh_.faces()[tetFacei_].size() - 1 - tetPti_;
// Reflect to account for the change of triangle orientation in the new cell
reflect();
// Transform the properties
if (!receiveCpp.parallel())
{
const tensor& T =
(
receiveCpp.forwardT().size() == 1
? receiveCpp.forwardT()[0]
: receiveCpp.forwardT()[receiveFacei]
);
transformProperties(T);
}
else if (receiveCpp.separated())
{
const vector& s =
(
(receiveCpp.separation().size() == 1)
? receiveCpp.separation()[0]
: receiveCpp.separation()[receiveFacei]
);
transformProperties(-s);
}
}
template<class TrackData>
void Foam::particle::hitCyclicAMIPatch
(
const cyclicAMIPolyPatch& cpp,
TrackData& td,
const vector& direction
)
{
vector pos = position();
const cyclicAMIPolyPatch& receiveCpp = cpp.neighbPatch();
const label sendFacei = cpp.whichFace(facei_);
const label receiveFacei = cpp.pointFace(sendFacei, direction, pos);
if (receiveFacei < 0)
{
// If the patch face of the particle is not known assume that the
// particle is lost and mark it to be deleted.
td.keepParticle = false;
WarningInFunction
<< "Particle lost across " << cyclicAMIPolyPatch::typeName
<< " patches " << cpp.name() << " and " << receiveCpp.name()
<< " at position " << pos << endl;
}
// Set the topology
facei_ = tetFacei_ = receiveFacei + receiveCpp.start();
// Locate the particle on the recieving side
vector directionT = direction;
cpp.reverseTransformDirection(directionT, sendFacei);
locate
(
pos,
&directionT,
mesh_.faceOwner()[facei_],
false,
"Particle crossed between " + cyclicAMIPolyPatch::typeName +
" patches " + cpp.name() + " and " + receiveCpp.name() +
" to a location outside of the mesh."
);
// The particle must remain associated with a face for the tracking to
// register as incomplete
facei_ = tetFacei_;
// Transform the properties
if (!receiveCpp.parallel())
{
const tensor& T =
(
receiveCpp.forwardT().size() == 1
? receiveCpp.forwardT()[0]
: receiveCpp.forwardT()[receiveFacei]
);
transformProperties(T);
}
else if (receiveCpp.separated())
{
const vector& s =
(
(receiveCpp.separation().size() == 1)
? receiveCpp.separation()[0]
: receiveCpp.separation()[receiveFacei]
);
transformProperties(-s);
}
}
template<class TrackData>
void Foam::particle::hitProcessorPatch(const processorPolyPatch&, TrackData&)
{}
template<class TrackData>
void Foam::particle::hitWallPatch
(
const wallPolyPatch&,
TrackData&,
const tetIndices&
)
{}
template<class TrackData>
void Foam::particle::hitPatch(const polyPatch&, TrackData&)
{}
// ************************************************************************* //