50 template<
class ReadGeoField,
class MappedGeoField>
57 const typename MappedGeoField::value_type& nullValue,
61 typedef typename MappedGeoField::value_type Type;
69 tetFields.
resize(fieldObjects.size());
75 << ReadGeoField::typeName <<
' ' <<
io.
name() <<
endl;
77 ReadGeoField readField(
io,
mesh);
92 readField.registerObject()
102 fld.setSize(map.
size(), nullValue);
105 label index = map[pointi];
109 label celli = index-1;
110 fld[pointi] = readField[celli];
114 const label facei = -index-1;
122 fld[pointi] = readField.boundaryField()[patchi][localFacei];
140 tetFields[i].correctBoundaryConditions();
147 int main(
int argc,
char *argv[])
155 "Convert polyMesh results to tetDualMesh" 170 Info<<
"Create tetDualMesh for time = " 188 "pointDualAddressing",
196 if (pointDualAddressing.size() != tetDualMesh.
nPoints())
199 <<
"Size " << pointDualAddressing.size()
200 <<
" of addressing map " << pointDualAddressing.objectPath()
201 <<
" differs from number of points in mesh " 209 label nPatchFaces = 0;
211 forAll(pointDualAddressing, pointi)
213 label index = pointDualAddressing[pointi];
225 label facei = -index-1;
229 <<
"Face " << facei <<
" from index " << index
230 <<
" is not a boundary face." 244 Info<<
"tetDualMesh points : " << tetDualMesh.
nPoints()
245 <<
" of which mapped to" <<
nl 246 <<
" cells : " << nCells <<
nl 247 <<
" patch faces : " << nPatchFaces <<
nl 248 <<
" not mapped : " << nUnmapped <<
nl 257 ReadAndMapFields<volScalarField, pointScalarField>
268 ReadAndMapFields<volVectorField, pointVectorField>
279 ReadAndMapFields<volSphericalTensorField, pointSphericalTensorField>
290 ReadAndMapFields<volSymmTensorField, pointSymmTensorField>
301 ReadAndMapFields<volTensorField, pointTensorField>
311 tetDualMesh.objectRegistry::write();
static void noFunctionObjects(bool addWithOption=false)
Remove '-noFunctionObjects' option and ignore any occurrences.
static void addNote(const string ¬e)
Add extra notes for the usage information.
void size(const label n)
Older name for setAddressableSize.
List of IOobjects with searching and retrieving facilities. Implemented as a HashTable, so the various sorted methods should be used if traversing in parallel.
errorManipArg< error, int > exit(error &err, const int errNo=1)
const fileName & facesInstance() const
Return the current instance directory for faces.
label nPoints() const noexcept
Number of mesh points.
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const word & name() const noexcept
Return the object name.
static const pointMesh & New(const polyMesh &mesh, Args &&... args)
Get existing or create MeshObject registered with typeName.
constexpr char nl
The newline '\n' character (0x0a)
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
Ostream & endl(Ostream &os)
Add newline and flush stream.
const labelList & patchID() const
Per boundary face label the patch index.
static bool setTimeIfPresent(Time &runTime, const argList &args, const bool forceInitial=false)
Set the runTime based on -constant (if present), -time (value), or -latestTime.
Generic dimensioned Type class.
#define forAll(list, i)
Loop across all elements in list.
Generic templated field type.
const polyBoundaryMesh & boundaryMesh() const noexcept
Return boundary mesh.
static void addOptions_singleTime()
Add single-time timeSelector options to argList::validOptions()
label nInternalFaces() const noexcept
Number of internal faces.
A list of pointers to objects of type <T>, without allocation/deallocation management of the pointers...
const T * set(const label i) const
Return const pointer to element (can be nullptr), or nullptr for out-of-range access (ie...
static word timeName(const scalar t, const int precision=precision_)
Return a time name for the given scalar time value formatted with the given precision.
void resize(const label newLen)
Adjust size of PtrList.
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< ' ';}gmvFile<< nl;for(const word &name :lagrangianScalarNames){ IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers...
Mesh data needed to do the Finite Volume discretisation.
Automatically write from objectRegistry::writeObject()
messageStream Info
Information stream (stdout output on master, null elsewhere)
IOobject io("surfaceFilmProperties", mesh.time().constant(), mesh, IOobject::READ_IF_PRESENT, IOobject::NO_WRITE, IOobject::NO_REGISTER)
Foam::argList args(argc, argv)
Defines the attributes of an object for which implicit objectRegistry management is supported...
UPtrList< const IOobject > csorted() const
The sorted list of IOobjects with headerClassName == Type::typeName.
void reduce(T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce inplace (cf. MPI Allreduce) using linear/tree communication schedule.
static constexpr const zero Zero
Global zero (0)