Commit 9a155dd0 authored by Henry Weller's avatar Henry Weller
Browse files

blockMesh: Added edge projection

New functionality contributed by Mattijs Janssens:
  - new edge projection: projectCurve for use with new geometry
    'searchableCurve'
  - new tutorial 'pipe'
  - naming of vertices and blocks (see pipe tutorial). Including back
    substitution for error messages.
parent d8e9decd
......@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
......@@ -92,6 +92,14 @@ public:
const bool doCollapse = false
);
//- Construct from components
inline cellShape
(
const word& model,
const labelList&,
const bool doCollapse = false
);
//- Construct from Istream
inline cellShape(Istream& is);
......
......@@ -25,6 +25,7 @@ License
#include "Istream.H"
#include "cell.H"
#include "cellModeller.H"
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
......@@ -51,6 +52,23 @@ inline Foam::cellShape::cellShape
}
inline Foam::cellShape::cellShape
(
const word& model,
const labelList& l,
const bool doCollapse
)
:
labelList(l),
m(cellModeller::lookup(model))
{
if (doCollapse)
{
collapse();
}
}
inline Foam::cellShape::cellShape(Istream& is)
{
is >> *this;
......
blockVertices/blockVertex/blockVertex.C
blockVertices/pointVertex/pointVertex.C
blockVertices/projectVertex/projectVertex.C
blockVertices/namedVertex/namedVertex.C
blockEdges/blockEdge/blockEdge.C
blockEdges/lineDivide/lineDivide.C
......@@ -13,6 +14,7 @@ blockEdges/BSplineEdge/BSplineEdge.C
blockEdges/splineEdge/CatmullRomSpline.C
blockEdges/splineEdge/splineEdge.C
blockEdges/projectEdge/projectEdge.C
blockEdges/projectCurveEdge/projectCurveEdge.C
blockFaces/blockFace/blockFace.C
blockFaces/projectFace/projectFace.C
......@@ -23,8 +25,9 @@ gradingDescriptor/gradingDescriptors.C
blockDescriptor/blockDescriptor.C
blockDescriptor/blockDescriptorEdges.C
block/block.C
block/blockCreate.C
blocks/block/block.C
blocks/block/blockCreate.C
blocks/namedBlock/namedBlock.C
blockMesh/blockMesh.C
blockMesh/blockMeshCreate.C
......@@ -33,4 +36,6 @@ blockMesh/blockMeshCheck.C
blockMesh/blockMeshMerge.C
blockMesh/blockMeshMergeFast.C
searchableCurve/searchableCurve.C
LIB = $(FOAM_LIBBIN)/libblockMesh
EXE_INC = \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/fileFormats/lnInclude \
-I$(LIB_SRC)/edgeMesh/lnInclude \
-I$(LIB_SRC)/surfMesh/lnInclude
LIB_LIBS = \
-lmeshTools \
-lfileFormats \
-ledgeMesh \
-lsurfMesh
......@@ -128,6 +128,90 @@ void Foam::blockDescriptor::findCurvedFaces()
}
void Foam::blockDescriptor::read
(
Istream& is,
label& val,
const dictionary& dict
)
{
token t(is);
if (t.isLabel())
{
val = t.labelToken();
}
else if (t.isWord())
{
const word& varName = t.wordToken();
const entry* ePtr = dict.lookupScopedEntryPtr
(
varName,
true,
true
);
if (ePtr)
{
// Read as label
val = Foam::readLabel(ePtr->stream());
}
else
{
FatalIOErrorInFunction(is)
<< "Undefined variable "
<< varName << ". Valid variables are " << dict
<< exit(FatalIOError);
}
}
else
{
FatalIOErrorInFunction(is)
<< "Illegal token " << t.info()
<< " when trying to read label"
<< exit(FatalIOError);
}
is.fatalCheck
(
"operator>>(Istream&, List<T>&) : reading entry"
);
}
Foam::label Foam::blockDescriptor::read
(
Istream& is,
const dictionary& dict
)
{
label val;
read(is, val, dict);
return val;
}
void Foam::blockDescriptor::write
(
Ostream& os,
const label val,
const dictionary& dict
)
{
forAllConstIter(dictionary, dict, iter)
{
if (iter().isStream())
{
label keyVal(Foam::readLabel(iter().stream()));
if (keyVal == val)
{
os << iter().keyword();
return;
}
}
}
os << val;
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::blockDescriptor::blockDescriptor
......@@ -164,6 +248,8 @@ Foam::blockDescriptor::blockDescriptor
Foam::blockDescriptor::blockDescriptor
(
const dictionary& dict,
const label index,
const pointField& vertices,
const blockEdgeList& edges,
const blockFaceList& faces,
......@@ -173,13 +259,24 @@ Foam::blockDescriptor::blockDescriptor
vertices_(vertices),
edges_(edges),
faces_(faces),
blockShape_(is),
density_(),
expand_(12, gradingDescriptors()),
zoneName_(),
curvedFaces_(-1),
nCurvedFaces_(0)
{
// Read cell model and list of vertices (potentially with variables)
word model(is);
blockShape_ = cellShape
(
model,
read<label>
(
is,
dict.subOrEmptyDict("namedVertices")
)
);
// Examine next token
token t(is);
......
......@@ -152,6 +152,8 @@ public:
//- Construct from Istream
blockDescriptor
(
const dictionary& dict,
const label index,
const pointField& vertices,
const blockEdgeList&,
const blockFaceList&,
......@@ -161,6 +163,9 @@ public:
// Member Functions
//- Reference to point field defining the block mesh
inline const pointField& vertices() const;
//- Return reference to the list of curved faces
inline const blockFaceList& faces() const;
......@@ -233,6 +238,23 @@ public:
// to lie on the faces of the block
void correctFacePoints(FixedList<pointField, 6>&) const;
//- In-place read with dictionary lookup
static void read(Istream&, label&, const dictionary&);
//- In-place read with dictionary lookup
template<class T>
static void read(Istream&, List<T>&, const dictionary&);
//- Return-read with dictionary lookup
static label read(Istream&, const dictionary&);
//- Return-read with dictionary lookup
template<class T>
static List<T> read(Istream& is, const dictionary&);
//- Write with dictionary lookup
static void write(Ostream&, const label, const dictionary&);
// IOstream Operators
......@@ -248,6 +270,10 @@ public:
#include "blockDescriptorI.H"
#ifdef NoRepository
#include "blockDescriptorTemplates.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
......
......@@ -25,6 +25,12 @@ License
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
inline const Foam::pointField& Foam::blockDescriptor::vertices() const
{
return vertices_;
}
inline const Foam::blockFaceList& Foam::blockDescriptor::faces() const
{
return faces_;
......
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class T>
void Foam::blockDescriptor::read
(
Istream& is,
List<T>& L,
const dictionary& dict
)
{
token firstToken(is);
if (firstToken.isLabel())
{
label s = firstToken.labelToken();
// Set list length to that read
L.setSize(s);
// Read list contents depending on data format
// Read beginning of contents
char delimiter = is.readBeginList("List");
if (s)
{
if (delimiter == token::BEGIN_LIST)
{
for (label i=0; i<s; i++)
{
read(is, L[i], dict);
}
}
}
// Read end of contents
is.readEndList("List");
}
else if (firstToken.isPunctuation())
{
if (firstToken.pToken() != token::BEGIN_LIST)
{
FatalIOErrorInFunction(is)
<< "incorrect first token, expected '(', found "
<< firstToken.info()
<< exit(FatalIOError);
}
SLList<T> sll;
while (true)
{
token t(is);
if (t.isPunctuation() && t.pToken() == token::END_LIST)
{
break;
}
is.putBack(t);
T elem;
read(is, elem, dict);
sll.append(elem);
}
// Convert the singly-linked list to this list
L = sll;
}
else
{
FatalIOErrorInFunction(is)
<< "incorrect first token, expected <int> or '(', found "
<< firstToken.info()
<< exit(FatalIOError);
}
}
template<class T>
Foam::List<T> Foam::blockDescriptor::read
(
Istream& is,
const dictionary& dict
)
{
List<T> L;
read(is, L, dict);
return L;
}
// ************************************************************************* //
......@@ -62,12 +62,14 @@ Foam::blockEdges::BSplineEdge::BSplineEdge
Foam::blockEdges::BSplineEdge::BSplineEdge
(
const dictionary& dict,
const label index,
const searchableSurfaces& geometry,
const pointField& points,
Istream& is
)
:
blockEdge(points, is),
blockEdge(dict, index, points, is),
BSpline(appendEndPoints(points, start_, end_, pointField(is)))
{
token t(is);
......
......@@ -83,6 +83,8 @@ public:
//- Construct from Istream, setting pointsList
BSplineEdge
(
const dictionary& dict,
const label index,
const searchableSurfaces& geometry,
const pointField&,
Istream&
......
......@@ -123,12 +123,14 @@ Foam::blockEdges::arcEdge::arcEdge
Foam::blockEdges::arcEdge::arcEdge
(
const dictionary& dict,
const label index,
const searchableSurfaces& geometry,
const pointField& points,
Istream& is
)
:
blockEdge(points, is),
blockEdge(dict, index, points, is),
p1_(points_[start_]),
p2_(is),
p3_(points_[end_]),
......
......@@ -91,6 +91,8 @@ public:
//- Construct from Istream setting pointsList
arcEdge
(
const dictionary& dict,
const label index,
const searchableSurfaces& geometry,
const pointField& points,
Istream&
......
......@@ -24,6 +24,7 @@ License
\*---------------------------------------------------------------------------*/
#include "blockEdge.H"
#include "blockDescriptor.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
......@@ -51,13 +52,15 @@ Foam::blockEdge::blockEdge
Foam::blockEdge::blockEdge
(
const dictionary& dict,
const label index,
const pointField& points,
Istream& is
)
:
points_(points),
start_(readLabel(is)),
end_(readLabel(is))
start_(blockDescriptor::read(is, dict.subOrEmptyDict("namedVertices"))),
end_(blockDescriptor::read(is, dict.subOrEmptyDict("namedVertices")))
{}
......@@ -70,6 +73,8 @@ Foam::autoPtr<Foam::blockEdge> Foam::blockEdge::clone() const
Foam::autoPtr<Foam::blockEdge> Foam::blockEdge::New
(
const dictionary& dict,
const label index,
const searchableSurfaces& geometry,
const pointField& points,
Istream& is
......@@ -95,7 +100,7 @@ Foam::autoPtr<Foam::blockEdge> Foam::blockEdge::New
<< abort(FatalError);
}
return autoPtr<blockEdge>(cstrIter()(geometry, points, is));
return autoPtr<blockEdge>(cstrIter()(dict, index, geometry, points, is));
}
......@@ -139,6 +144,25 @@ Foam::blockEdge::position(const scalarList& lambdas) const
}
void Foam::blockEdge::write(Ostream& os, const dictionary& d) const
{
const dictionary* varDictPtr = d.subDictPtr("namedVertices");
if (varDictPtr)
{
const dictionary& varDict = *varDictPtr;
blockDescriptor::write(os, start_, varDict);
os << tab;
blockDescriptor::write(os, end_, varDict);
os << endl;
}
else
{
os << start_ << tab << end_ << endl;
}
}
// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
Foam::Ostream& Foam::operator<<(Ostream& os, const blockEdge& p)
......
......@@ -92,11 +92,13 @@ public:
blockEdge,
Istream,
(
const dictionary& dict,
const label index,
const searchableSurfaces& geometry,
const pointField& points,
Istream& is
),
(geometry, points, is)
(dict, index, geometry, points, is)
);
......@@ -113,6 +115,8 @@ public:
//- Construct from Istream setting pointsList
blockEdge
(
const dictionary& dict,
const label index,
const pointField&,
Istream&
);
......@@ -123,6 +127,8 @@ public:
//- New function which constructs and returns pointer to a blockEdge
static autoPtr<blockEdge> New
(
const dictionary& dict,
const label index,
const searchableSurfaces& geometry,
const pointField&,
Istream&
......@@ -132,20 +138,29 @@ public:
// PtrLists of blockEdge
class iNew
{
const dictionary& dict_;
const searchableSurfaces& geometry_;
const pointField& points_;
mutable label index_;
public:
iNew(const searchableSurfaces& geometry, const pointField& points)
iNew
(
const dictionary& dict,
const searchableSurfaces& geometry,
const pointField& points
)
:
dict_(dict),
geometry_(geometry),
points_(points)
points_(points),
index_(0)
{}
autoPtr<blockEdge> operator()(Istream& is) const
{
return blockEdge::New(geometry_, points_, is);
return blockEdge::New(dict_, index_++, geometry_, points_, is);
}
};
......@@ -195,6 +210,9 @@ public:
//- Return the length of the curve
virtual scalar length() const = 0;
//- Write edge with variable backsubstitution
void write(Ostream&, const dictionary&) const;
// Ostream operator
......
......@@ -53,12 +53,14 @@ Foam::blockEdges::lineEdge::lineEdge
Foam::blockEdges::lineEdge::lineEdge
(
const dictionary& dict,
const label index,
const searchableSurfaces& geometry,