Loading...
Searching...
No Matches
writePointFields.H
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2020-2025 OpenCFD Ltd.
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
12
13InNamespace
14 Foam
15
16Description
17 Read point fields from disk and write as ensight data
18
19\*---------------------------------------------------------------------------*/
20
21#ifndef FoamToEnsight_writePointFields_H
22#define FoamToEnsight_writePointFields_H
23
24#include "readFields.H"
25#include "ensightMesh.H"
26#include "fvMesh.H"
27
28// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
29
30namespace Foam
31{
32
33template<class Type>
35(
36 ensightCase& ensCase,
37 const ensightMesh& ensMesh,
38 const tmp<PointField<Type>>& tfield
39)
40{
41 if (!tfield)
42 {
43 return false;
44 }
45
46 const auto& field = tfield();
47
48 // PointData = true
49 autoPtr<ensightFile> os = ensCase.newData<Type>(field.name(), true);
50
52 (
53 os.ref(),
54 field,
55 ensMesh
56 );
57
58 tfield.clear();
59 return wrote;
60}
61
62
63template<class Type>
65(
66 ensightCase& ensCase,
67 const ensightMesh& ensMesh,
68 const IOobjectList& objects
69)
70{
71 typedef PointField<Type> FieldType;
72
73 const pointMesh& pMesh = pointMesh::New(ensMesh.mesh());
74
75 label count = 0;
76
77 for (const IOobject& io : objects.csorted<FieldType>())
78 {
79 if
80 (
82 (
83 ensCase,
84 ensMesh,
86 )
87 )
88 {
89 Info<< ' ' << io.name();
90 ++count;
91 }
92 }
93
94 return count;
95}
96
97
99(
100 ensightCase& ensCase,
101 const ensightMesh& ensMesh,
102 const IOobjectList& objects
103)
104{
105 label count = 0;
106 const label total = objects.size();
107 do
108 {
109 #undef doLocalWriteCode
110 #define doLocalWriteCode(Type) \
111 { \
112 count += writePointFields<Type> \
113 ( \
114 ensCase, \
115 ensMesh, \
116 objects \
117 ); \
118 if (count >= total) break; /* early exit */ \
119 }
120
121 doLocalWriteCode(scalar);
126
127 #undef doLocalWriteCode
128 }
129 while (false);
130
131 return count;
132}
133
134} // End namespace Foam
135
136// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
137
138#endif
139
140// ************************************************************************* //
Helper routines for reading a field or fields, for foamToEnsight.
label size() const noexcept
The number of elements in table.
Definition HashTable.H:358
List of IOobjects with searching and retrieving facilities. Implemented as a HashTable,...
UPtrList< const IOobject > csorted() const
The sorted list of IOobjects with headerClassName == Type::typeName.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition IOobject.H:191
static FOAM_NO_DANGLING_REFERENCE const pointMesh & New(const polyMesh &mesh, Args &&... args)
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition autoPtr.H:65
Supports writing of ensight cases as well as providing common factory methods to open new files.
Definition ensightCase.H:64
autoPtr< ensightFile > newData(const word &varName, const bool isPointData=false) const
Open stream for new data file (on master), with current index.
Encapsulation of volume meshes for writing in ensight format. It manages cellZones,...
Definition ensightMesh.H:79
const polyMesh & mesh() const noexcept
Reference to the underlying polyMesh.
Mesh representing a set of points created from polyMesh.
Definition pointMesh.H:49
A class for managing temporary objects.
Definition tmp.H:75
rDeltaTY field()
OBJstream os(runTime.globalPath()/outputName)
const auto & io
#define doLocalWriteCode(Type)
bool writePointField(ensightOutput::floatBufferType &scratch, ensightFile &os, const GeometricField< Type, pointPatchField, pointMesh > &pf, const ensightMesh &ensMesh)
Write point field component-wise.
Namespace for OpenFOAM.
label writeAllPointFields(ensightCase &ensCase, const ensightMesh &ensMesh, const IOobjectList &objects)
messageStream Info
Information stream (stdout output on master, null elsewhere).
Tensor< scalar > tensor
Definition symmTensor.H:57
bool writePointField(ensightCase &ensCase, const ensightMesh &ensMesh, const tmp< PointField< Type > > &tfield)
GeometricField< Type, pointPatchField, pointMesh > PointField
A point field for a given type.
label writePointFields(ensightCase &ensCase, const ensightMesh &ensMesh, const IOobjectList &objects)
Vector< scalar > vector
Definition vector.H:57
SphericalTensor< scalar > sphericalTensor
SphericalTensor of scalars, i.e. SphericalTensor<scalar>.
tmp< GeoField > getField(const IOobject &io, const typename GeoField::Mesh &mesh)
Get the field or FatalError.
Definition readFields.H:51
SymmTensor< scalar > symmTensor
SymmTensor of scalars, i.e. SymmTensor<scalar>.
Definition symmTensor.H:55