Go to the documentation of this file.
64 Info<<
"Updating vertex markup. curLeft: "
65 << curLeft <<
" curRight: " << curRight <<
endl;
67 tmp<scalarField> tvertexMarkup(
new scalarField(
p.size()));
72 if (
p[pI].
x() < curLeft - SMALL)
74 vertexMarkup[pI] = -1;
76 else if (
p[pI].
x() > curRight + SMALL)
90 void Foam::movingConeTopoFvMesh::addZonesAndModifiers()
99 || topoChanger_.size()
103 <<
"Zones and modifiers already present. Skipping."
109 Info<<
"Time = " << time().timeName() <<
endl
110 <<
"Adding zones and modifiers to the mesh" <<
endl;
116 boolList flipZone1(fc.size(),
false);
117 label nZoneFaces1 = 0;
120 boolList flipZone2(fc.size(),
false);
121 label nZoneFaces2 = 0;
127 fc[facei].
x() > -0.003501
128 && fc[facei].
x() < -0.003499
131 if ((fa[facei] &
vector(1, 0, 0)) < 0)
133 flipZone1[nZoneFaces1] =
true;
136 zone1[nZoneFaces1] = facei;
137 Info<<
"face " << facei <<
" for zone 1. Flip: "
138 << flipZone1[nZoneFaces1] <<
endl;
143 fc[facei].
x() > -0.00701
144 && fc[facei].
x() < -0.00699
147 zone2[nZoneFaces2] = facei;
149 if ((fa[facei] &
vector(1, 0, 0)) > 0)
151 flipZone2[nZoneFaces2] =
true;
154 Info<<
"face " << facei <<
" for zone 2. Flip: "
155 << flipZone2[nZoneFaces2] <<
endl;
160 zone1.setSize(nZoneFaces1);
161 flipZone1.setSize(nZoneFaces1);
163 zone2.setSize(nZoneFaces2);
164 flipZone2.setSize(nZoneFaces2);
169 List<pointZone*> pz(0);
170 List<faceZone*> fz(2);
171 List<cellZone*> cz(0);
178 "rightExtrusionFaces",
180 std::move(flipZone1),
189 "leftExtrusionFaces",
191 std::move(flipZone2),
199 Info<<
"Adding mesh zones." <<
endl;
200 addZones(pz, fz, cz);
205 List<polyMeshModifier*> tm(2);
209 new layerAdditionRemoval
214 "rightExtrusionFaces",
215 motionDict_.subDict(
"right").get<scalar>(
"minThickness"),
216 motionDict_.subDict(
"right").get<scalar>(
"maxThickness")
220 tm[nMods] =
new layerAdditionRemoval
225 "leftExtrusionFaces",
226 motionDict_.subDict(
"left").get<scalar>(
"minThickness"),
227 motionDict_.subDict(
"left").get<scalar>(
"maxThickness")
232 Info<<
"Adding " << nMods <<
" mesh modifiers" <<
endl;
233 topoChanger_.addTopologyModifiers(tm);
241 Foam::movingConeTopoFvMesh::movingConeTopoFvMesh(
const IOobject& io)
257 ).optionalSubDict(typeName +
"Coeffs")
259 motionVelAmplitude_(motionDict_.get<
vector>(
"motionVelAmplitude")),
260 motionVelPeriod_(motionDict_.get<scalar>(
"motionVelPeriod")),
263 motionVelAmplitude_*
sin(time().value()*
pi/motionVelPeriod_)
265 leftEdge_(motionDict_.get<scalar>(
"leftEdge")),
266 curLeft_(motionDict_.get<scalar>(
"leftObstacleEdge")),
267 curRight_(motionDict_.get<scalar>(
"rightObstacleEdge"))
270 <<
" Initial curMotionVel_:" << curMotionVel_
273 addZonesAndModifiers();
279 faceZones().findZoneID(
"leftExtrusionFaces")
287 faceZones().findZoneID(
"rightExtrusionFaces")
291 motionMask_ = vertexMarkup
318 motionVelAmplitude_*
sin(time().value()*
pi/motionVelPeriod_);
320 Pout<<
"time:" << time().value() <<
" curMotionVel_:" << curMotionVel_
321 <<
" curLeft:" << curLeft_ <<
" curRight:" << curRight_
324 if (topoChangeMap.
valid())
326 Info<<
"Topology change. Calculating motion points" <<
endl;
328 if (topoChangeMap().hasMotionPoints())
330 Info<<
"Topology change. Has premotion points" <<
endl;
335 topoChangeMap().preMotionPoints(),
342 topoChangeMap().preMotionPoints()
345 )*curMotionVel_*time().deltaT().value();
349 Info<<
"Topology change. Already set mesh points" <<
endl;
364 )*curMotionVel_*time().deltaT().value();
369 Info<<
"No topology change" <<
endl;
375 )*curMotionVel_*time().deltaT().value();
379 Info <<
"Executing mesh motion" <<
endl;
380 movePoints(newPoints);
388 faceZones().findZoneID(
"leftExtrusionFaces")
396 faceZones().findZoneID(
"rightExtrusionFaces")
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
List< label > labelList
A List of labels.
vectorField pointField
pointField is a vectorField.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
virtual const pointField & points() const
Return raw points.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
#define InfoInFunction
Report an information message using Foam::Info.
virtual ~movingConeTopoFvMesh()
Destructor.
A class for managing temporary objects.
dimensionedScalar sin(const dimensionedScalar &ds)
List< bool > boolList
A List of bools.
bool valid() const noexcept
True if the managed pointer is non-null.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const Type & value() const
Return const reference to value.
dimensionedScalar pos0(const dimensionedScalar &ds)
prefixOSstream Pout
An Ostream wrapper for parallel output to std::cout.
virtual bool update()
Update the mesh for both mesh motion and topology change.
#define forAll(list, i)
Loop across all elements in list.
Field< vector > vectorField
Specialisation of Field<T> for vector.
const faceZoneMesh & faceZones() const
Return face zone mesh.
messageStream Info
Information stream (uses stdout - output is on the master only)
Macros for easy insertion into run-time selection tables.
Vector< scalar > vector
A scalar version of the templated Vector.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
constexpr scalar pi(M_PI)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Abstract base class for a topology changing fvMesh.
void write(vtk::formatter &fmt, const Type &val, const label n=1)
Component-wise write of a value (N times)
const Time & time() const
Return the top-level database.
constant condensation/saturation model.
defineTypeNameAndDebug(combustionModel, 0)
dimensioned< Type > average(const DimensionedField< Type, GeoMesh > &df)