singleCellMesh.C
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) 2011-2015 OpenFOAM Foundation
9 Copyright (C) 2016-2021 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27Application
28 singleCellMesh
29
30Group
31 grpMeshManipulationUtilities
32
33Description
34 Reads all fields and maps them to a mesh with all internal faces removed
35 (singleCellFvMesh) which gets written to region "singleCell".
36
37 Used to generate mesh and fields that can be used for boundary-only data.
38 Might easily result in illegal mesh though so only look at boundaries
39 in paraview.
40
41\*---------------------------------------------------------------------------*/
42
43#include "argList.H"
44#include "fvMesh.H"
45#include "volFields.H"
46#include "Time.H"
47#include "ReadFields.H"
48#include "singleCellFvMesh.H"
49#include "timeSelector.H"
50
51using namespace Foam;
52
53// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
54
55// Name of region to create
56const string singleCellName = "singleCell";
57
58
59template<class GeoField>
60void interpolateFields
61(
62 const singleCellFvMesh& scMesh,
63 const PtrList<GeoField>& flds
64)
65{
66 forAll(flds, i)
67 {
68 tmp<GeoField> scFld = scMesh.interpolate(flds[i]);
69 GeoField* scFldPtr = scFld.ptr();
70 scFldPtr->writeOpt(IOobject::AUTO_WRITE);
71 scFldPtr->store();
72 }
73}
74
75
76
77int main(int argc, char *argv[])
78{
79 argList::addNote
80 (
81 "Map fields to a mesh with all internal faces removed"
82 " (singleCellFvMesh) which gets written to region 'singleCell'"
83 );
84
85 timeSelector::addOptions(true, false); // constant(true), zero(false)
86
87 #include "setRootCase.H"
88 #include "createTime.H"
89
90 instantList timeDirs = timeSelector::select0(runTime, args);
91
92 #include "createNamedMesh.H"
93
94 if (regionName == singleCellName)
95 {
97 << "Cannot convert region " << singleCellName
98 << " since result would overwrite it. Please rename your region."
99 << exit(FatalError);
100 }
101
102 // Create the mesh
103 Info<< "Creating singleCell mesh" << nl << endl;
105 (
107 (
109 (
110 singleCellName,
111 mesh.pointsInstance(),
112 runTime,
113 IOobject::NO_READ,
114 IOobject::AUTO_WRITE
115 ),
116 mesh
117 )
118 );
119 scMesh().setInstance(mesh.pointsInstance());
120
121 // For convenience create any fv* files
122 if (!exists(scMesh().fvSolution::objectPath()))
123 {
124 mkDir(scMesh().fvSolution::path());
125 ln("../fvSolution", scMesh().fvSolution::objectPath());
126 }
127 if (!exists(scMesh().fvSchemes::objectPath()))
128 {
129 mkDir(scMesh().fvSolution::path());
130 ln("../fvSchemes", scMesh().fvSchemes::objectPath());
131 }
132
133
134 forAll(timeDirs, timeI)
135 {
136 runTime.setTime(timeDirs[timeI], timeI);
137
138 Info<< nl << "Time = " << runTime.timeName() << endl;
139
140
141 // Check for new mesh
142 if (mesh.readUpdate() != polyMesh::UNCHANGED)
143 {
144 Info<< "Detected changed mesh. Recreating singleCell mesh." << endl;
145 scMesh.clear(); // remove any registered objects
146 scMesh.reset
147 (
149 (
151 (
152 singleCellName,
153 mesh.pointsInstance(),
154 runTime,
155 IOobject::NO_READ,
156 IOobject::AUTO_WRITE
157 ),
158 mesh
159 )
160 );
161 }
162
163
164 // Read objects in time directory
165 IOobjectList objects(mesh, runTime.timeName());
166
167 // Read vol fields.
169 ReadFields(mesh, objects, vsFlds);
170
172 ReadFields(mesh, objects, vvFlds);
173
175 ReadFields(mesh, objects, vstFlds);
176
178 ReadFields(mesh, objects, vsymtFlds);
179
181 ReadFields(mesh, objects, vtFlds);
182
183 // Map and store the fields on the scMesh.
184 interpolateFields(scMesh(), vsFlds);
185 interpolateFields(scMesh(), vvFlds);
186 interpolateFields(scMesh(), vstFlds);
187 interpolateFields(scMesh(), vsymtFlds);
188 interpolateFields(scMesh(), vtFlds);
189
190
191 // Write
192 Info<< "Writing mesh to time " << runTime.timeName() << endl;
193 scMesh().write();
194 }
195
196
197 Info<< "End\n" << endl;
198
199 return 0;
200}
201
202
203// ************************************************************************* //
Field reading functions for post-processing utilities.
List of IOobjects with searching and retrieving facilities.
Definition: IOobjectList.H:59
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: PtrList.H:73
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: autoPtr.H:66
virtual bool write(const bool valid=true) const
Write mesh using IO settings from time.
Definition: fvMesh.C:1079
void clear()
Clear all entries from the registry.
void setInstance(const fileName &instance, const IOobject::writeOption wOpt=IOobject::AUTO_WRITE)
Set the instance for mesh files.
Definition: polyMeshIO.C:36
void reset(const label nPoints, const label nInternalFaces, const label nFaces, const label nCells)
Reset this primitiveMesh given the primitive array sizes.
fvMesh as subset of other mesh. Consists of one cell and all original bounday faces....
tmp< GeometricField< Type, fvPatchField, volMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &) const
A class for managing temporary objects.
Definition: tmp.H:65
T * ptr() const
Return managed pointer for reuse, or clone() the object reference.
Definition: tmpI.H:255
dynamicFvMesh & mesh
engineTime & runTime
Foam::word regionName(Foam::polyMesh::defaultRegion)
Required Variables.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Namespace for OpenFOAM.
wordList ReadFields(const typename GeoMesh::Mesh &mesh, const IOobjectList &objects, PtrList< GeometricField< Type, PatchField, GeoMesh > > &fields, const bool syncPar=true, const bool readOldTime=false)
Read Geometric fields of templated type.
bool exists(const fileName &name, const bool checkGzip=true, const bool followLink=true)
Does the name exist (as DIRECTORY or FILE) in the file system?
Definition: MSwindows.C:633
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
bool ln(const fileName &src, const fileName &dst)
Create a softlink. dst should not exist. Returns true if successful.
Definition: MSwindows.C:933
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
Foam::argList args(argc, argv)
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
mkDir(pdfPath)