datToFoam.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-2016 OpenFOAM Foundation
9 Copyright (C) 2020-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 datToFoam
29
30Group
31 grpMeshConversionUtilities
32
33Description
34 Reads in a datToFoam mesh file and outputs a points file.
35 Used in conjunction with blockMesh.
36
37\*---------------------------------------------------------------------------*/
38
39#include "argList.H"
40#include "Time.H"
41#include "Fstream.H"
42#include "polyMesh.H"
43#include "unitConversion.H"
44
45using namespace Foam;
46
47// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
48
49int main(int argc, char *argv[])
50{
51 argList::addNote
52 (
53 "Reads in a datToFoam mesh file and outputs a points file.\n"
54 "Used in conjunction with blockMesh."
55 );
56
57 argList::noParallel();
58 argList::addArgument("dat file");
59
60 argList args(argc, argv);
61
62 if (!args.check())
63 {
65 }
66
67 #include "createTime.H"
68
69 std::ifstream plot3dFile(args.get<fileName>(1));
70
71 string line;
72 std::getline(plot3dFile, line);
73 std::getline(plot3dFile, line);
74
75 IStringStream Istring(line);
76 word block;
77 string zoneName;
78 token punctuation;
79 label iPoints, jPoints;
80
81 Istring >> block;
82 Istring >> block;
83 Istring >> zoneName;
84 Istring >> punctuation;
85 Istring >> block;
86 Istring >> iPoints;
87 Istring >> block;
88 Istring >> jPoints;
89
90 Info<< "Number of vertices in i direction = " << iPoints << nl
91 << "Number of vertices in j direction = " << jPoints << endl;
92
93 // We ignore the first layer of points in i and j the biconic meshes
94 const label nPointsij = (iPoints - 1)*(jPoints - 1);
95
96 pointField points(nPointsij, Zero);
97
98 for (direction comp = 0; comp < 2; ++comp)
99 {
100 label p = 0;
101
102 for (label j = 0; j < jPoints; ++j)
103 {
104 for (label i = 0; i < iPoints; ++i)
105 {
106 scalar coord;
107 plot3dFile >> coord;
108
109 // if statement ignores the first layer in i and j
110 if (i>0 && j>0)
111 {
112 points[p++][comp] = coord;
113 }
114 }
115 }
116 }
117
118 // correct error in biconic meshes
119 for (point& p : points)
120 {
121 if (p.y() < 1e-07)
122 {
123 p.y() = 0;
124 }
125 }
126
127 pointField pointsWedge(nPointsij*2, Zero);
128
129 const scalar a = degToRad(0.1);
130 const tensor rotateZ
131 (
132 1.0, 0.0, 0.0,
133 0.0, 1.0, 0.0,
134 0.0, -::sin(a), ::cos(a)
135 );
136
137 forAll(points, i)
138 {
139 pointsWedge[i] = (rotateZ & points[i]);
140 pointsWedge[i + nPointsij] =
142 (
143 vector(1.0, 1.0, -1.0),
144 pointsWedge[i]
145 );
146 }
147
148 // Write (overwrite) points
149 {
150 IOobject iopoints
151 (
152 "points",
153 runTime.constant(), // instance
154 polyMesh::meshSubDir, // local
155 runTime
156 );
157
158 if (!exists(iopoints.path()))
159 {
160 mkDir(iopoints.path());
161 }
162
163 OFstream os(iopoints.objectPath(), runTime.writeStreamOption());
164
165 Info<< "Writing points to: " << nl
166 << " " << os.name() << endl;
167
168 iopoints.writeHeader(os, vectorIOField::typeName);
169 os << pointsWedge;
170 IOobject::writeEndDivider(os);
171 }
172
173 Info<< "End\n" << endl;
174
175 return 0;
176}
177
178
179// ************************************************************************* //
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
Input from string buffer, using a ISstream. Always UNCOMPRESSED.
Definition: StringStream.H:112
Output to file stream, using an OSstream.
Definition: OFstream.H:57
Extract command arguments and options from the supplied argc and argv parameters.
Definition: argList.H:124
T get(const label index) const
Get a value from the argument at index.
Definition: argListI.H:278
bool check(bool checkArgs=argList::argsMandatory(), bool checkOpts=true) const
Definition: argList.C:1913
Creates a single block of cells from point coordinates, numbers of cells in each direction and an exp...
Definition: block.H:61
void exit(const int errNo=1)
Exit : can be called for any error to exit program.
Definition: error.C:331
A class for handling file names.
Definition: fileName.H:76
A line primitive.
Definition: line.H:68
Tensor of scalars, i.e. Tensor<scalar>.
A token holds an item read from Istream.
Definition: token.H:69
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
A class for handling words, derived from Foam::string.
Definition: word.H:68
volScalarField & p
engineTime & runTime
OBJstream os(runTime.globalPath()/outputName)
const pointField & points
Namespace for OpenFOAM.
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)
constexpr scalar degToRad(const scalar deg) noexcept
Conversion from degrees to radians.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
uint8_t direction
Definition: direction.H:56
dimensioned< Type > cmptMultiply(const dimensioned< Type > &, const dimensioned< Type > &)
error FatalError
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
Foam::argList args(argc, argv)
volScalarField & e
Definition: createFields.H:11
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
Unit conversion functions.
mkDir(pdfPath)