foamDataToFluent.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 | Copyright (C) 2011-2016 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 Application
25  foamDataToFluent
26 
27 Description
28  Translates OpenFOAM data to Fluent format.
29 
30 \*---------------------------------------------------------------------------*/
31 
32 #include "fvCFD.H"
33 #include "writeFluentFields.H"
34 #include "OFstream.H"
35 #include "IOobjectList.H"
36 
37 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
38 
39 int main(int argc, char *argv[])
40 {
41  argList::noParallel();
42  timeSelector::addOptions(false); // no constant
43 
44  #include "setRootCase.H"
45  #include "createTime.H"
46 
47  instantList timeDirs = timeSelector::select0(runTime, args);
48 
49  #include "createMesh.H"
50 
51  // make a directory called proInterface in the case
52  mkDir(runTime.rootPath()/runTime.caseName()/"fluentInterface");
53 
54  forAll(timeDirs, timeI)
55  {
56  runTime.setTime(timeDirs[timeI], timeI);
57 
58  Info<< "Time = " << runTime.timeName() << endl;
59 
60  if (mesh.readUpdate())
61  {
62  Info<< " Read new mesh" << endl;
63  }
64 
65  // make a directory called proInterface in the case
66  mkDir(runTime.rootPath()/runTime.caseName()/"fluentInterface");
67 
68  // open a file for the mesh
69  OFstream fluentDataFile
70  (
71  runTime.rootPath()/
72  runTime.caseName()/
73  "fluentInterface"/
74  runTime.caseName() + runTime.timeName() + ".dat"
75  );
76 
77  fluentDataFile
78  << "(0 \"FOAM to Fluent data File\")" << endl << endl;
79 
80  // Writing number of faces
81  label nFaces = mesh.nFaces();
82 
83  forAll(mesh.boundary(), patchi)
84  {
85  nFaces += mesh.boundary()[patchi].size();
86  }
87 
88  fluentDataFile
89  << "(33 (" << mesh.nCells() << " " << nFaces << " "
90  << mesh.nPoints() << "))" << endl;
91 
92  IOdictionary foamDataToFluentDict
93  (
94  IOobject
95  (
96  "foamDataToFluentDict",
97  runTime.system(),
98  mesh,
99  IOobject::MUST_READ_IF_MODIFIED,
100  IOobject::NO_WRITE
101  )
102  );
103 
104 
105  // Search for list of objects for this time
106  IOobjectList objects(mesh, runTime.timeName());
107 
108 
109  // Converting volScalarField
110  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
111 
112  // Search list of objects for volScalarFields
113  IOobjectList scalarFields(objects.lookupClass("volScalarField"));
114 
115  forAllIter(IOobjectList, scalarFields, iter)
116  {
117  // Read field
118  volScalarField field(*iter(), mesh);
119 
120  // lookup field from dictionary and convert field
121  label unitNumber;
122  if
123  (
124  foamDataToFluentDict.readIfPresent(field.name(), unitNumber)
125  && unitNumber > 0
126  )
127  {
128  Info<< " Converting field " << field.name() << endl;
129  writeFluentField(field, unitNumber, fluentDataFile);
130  }
131  }
132 
133 
134  // Converting volVectorField
135  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
136 
137  // Search list of objects for volVectorFields
138  IOobjectList vectorFields(objects.lookupClass("volVectorField"));
139 
140  forAllIter(IOobjectList, vectorFields, iter)
141  {
142  // Read field
143  volVectorField field(*iter(), mesh);
144 
145  // lookup field from dictionary and convert field
146  label unitNumber;
147  if
148  (
149  foamDataToFluentDict.readIfPresent(field.name(), unitNumber)
150  && unitNumber > 0
151  )
152  {
153  Info<< " Converting field " << field.name() << endl;
154  writeFluentField(field, unitNumber, fluentDataFile);
155  }
156  }
157 
158  Info<< endl;
159  }
160 
161  Info<< "End\n" << endl;
162 
163  return 0;
164 }
165 
166 
167 // ************************************************************************* //
List< instant > instantList
List of instants.
Definition: instantList.H:42
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
#define forAllIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:453
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
dynamicFvMesh & mesh
void writeFluentField(const volScalarField &phi, const label fluentFieldIdentifier, Ostream &stream)
bool mkDir(const fileName &, mode_t=0777)
Make a directory and return an error if it could not be created.
Definition: POSIX.C:295
label patchi
messageStream Info
Foam::argList args(argc, argv)