Commit e22f1b35 authored by Henry Weller's avatar Henry Weller
Browse files

refineMesh: Added customizable direction fields to multiDirRefinement

Contribution provided by Bruno Santos
Resolves feature request http://www.openfoam.org/mantisbt/view.php?id=2004
parent 4758c2ac
......@@ -22,8 +22,12 @@ set c0;
// x,y,z axis. Specify in globalCoeffs section below.
// - patchLocal : coordinate system different for every cell. Specify in
// patchLocalCoeffs section below.
// - fieldBased : uses the list of field names from the directions list for
// selecting the directions to cut. Meant to be used with geometricCut, but
// can also be used with useHexTopology.
coordinateSystem global;
//coordinateSystem patchLocal;
//coordinateSystem fieldBased;
// .. and its coefficients. x,y in this case. (normal direction is calculated
// as tan1^tan2)
......@@ -39,7 +43,7 @@ patchLocalCoeffs
tan1 (1 0 0);
}
// List of directions to refine
// List of directions to refine, if global or patchLocal
directions
(
tan1
......@@ -47,6 +51,15 @@ directions
normal
);
// List of directions to refine, if "fieldBased". Keep in mind that these
// fields must be of type "vectorField", not "volVectorField".
//directions
//(
// radialDirectionFieldName
// angularDirectionFieldName
// heightDirectionFieldName
//);
// Whether to use hex topology. This will
// - if patchLocal: all cells on selected patch should be hex
// - split all hexes in 2x2x2 through the middle of edges.
......
......@@ -38,7 +38,7 @@ License
namespace Foam
{
defineTypeNameAndDebug(cellCuts, 0);
defineTypeNameAndDebug(cellCuts, 0);
}
......@@ -120,7 +120,7 @@ void Foam::cellCuts::writeUncutOBJ
const label cellI
) const
{
//- Cell edges
// Cell edges
OFstream cutsStream(dir / "cell_" + name(cellI) + ".obj");
Pout<< "Writing cell for time " << mesh().time().timeName()
......@@ -135,7 +135,7 @@ void Foam::cellCuts::writeUncutOBJ
labelList(1, cellI)
);
//- Loop cutting cell in two
// Loop cutting cell in two
OFstream cutStream(dir / "cellCuts_" + name(cellI) + ".obj");
Pout<< "Writing raw cuts on cell for time " << mesh().time().timeName()
......@@ -180,7 +180,7 @@ void Foam::cellCuts::writeOBJ
const labelList& anchors
) const
{
//- Cell edges
// Cell edges
OFstream cutsStream(dir / "cell_" + name(cellI) + ".obj");
Pout<< "Writing cell for time " << mesh().time().timeName()
......@@ -196,7 +196,7 @@ void Foam::cellCuts::writeOBJ
);
//- Loop cutting cell in two
// Loop cutting cell in two
OFstream loopStream(dir / "cellLoop_" + name(cellI) + ".obj");
Pout<< "Writing loop for time " << mesh().time().timeName()
......@@ -207,7 +207,7 @@ void Foam::cellCuts::writeOBJ
writeOBJ(loopStream, loopPoints, vertI);
//- Anchors for cell
// Anchors for cell
OFstream anchorStream(dir / "anchors_" + name(cellI) + ".obj");
Pout<< "Writing anchors for time " << mesh().time().timeName()
......@@ -2240,7 +2240,7 @@ void Foam::cellCuts::setFromCellCutter
WarningInFunction
<< "Found loop on cell " << cellI
<< " that resulted in an unexpected bad cut."
<< " that resulted in an unexpected bad cut." << nl
<< " Suggestions:" << nl
<< " - Turn on the debug switch for 'cellCuts' to get"
<< " geometry files that identify this cell." << nl
......
......@@ -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
......@@ -57,13 +57,12 @@ const Foam::NamedEnum<Foam::directions::directionType, 3>
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
// For debugging
void Foam::directions::writeOBJ(Ostream& os, const point& pt)
{
os << "v " << pt.x() << ' ' << pt.y() << ' ' << pt.z() << endl;
}
// For debugging
void Foam::directions::writeOBJ
(
Ostream& os,
......@@ -81,7 +80,6 @@ void Foam::directions::writeOBJ
}
// Dump to file.
void Foam::directions::writeOBJ
(
const fileName& fName,
......@@ -137,7 +135,6 @@ void Foam::directions::check2D
}
// Get direction on all cells
Foam::vectorField Foam::directions::propagateDirection
(
const polyMesh& mesh,
......@@ -283,34 +280,35 @@ Foam::directions::directions
List<vectorField>(wordList(dict.lookup("directions")).size())
{
const wordList wantedDirs(dict.lookup("directions"));
const word coordSystem(dict.lookup("coordinateSystem"));
bool wantNormal = false;
bool wantTan1 = false;
bool wantTan2 = false;
label nDirs = 0;
forAll(wantedDirs, i)
if (coordSystem != "fieldBased")
{
directionType wantedDir = directionTypeNames_[wantedDirs[i]];
if (wantedDir == NORMAL)
forAll(wantedDirs, i)
{
wantNormal = true;
}
else if (wantedDir == TAN1)
{
wantTan1 = true;
}
else if (wantedDir == TAN2)
{
wantTan2 = true;
directionType wantedDir = directionTypeNames_[wantedDirs[i]];
if (wantedDir == NORMAL)
{
wantNormal = true;
}
else if (wantedDir == TAN1)
{
wantTan1 = true;
}
else if (wantedDir == TAN2)
{
wantTan2 = true;
}
}
}
label nDirs = 0;
const word coordSystem(dict.lookup("coordinateSystem"));
if (coordSystem == "global")
{
const dictionary& globalDict = dict.subDict("globalCoeffs");
......@@ -424,12 +422,29 @@ Foam::directions::directions
this->operator[](nDirs++) = tan2Dirs;
}
}
else if (coordSystem == "fieldBased")
{
forAll(wantedDirs, i)
{
operator[](nDirs++) =
vectorIOField
(
IOobject
(
mesh.instance()/wantedDirs[i],
mesh,
IOobject::MUST_READ,
IOobject::NO_WRITE
)
);
}
}
else
{
FatalErrorInFunction
<< "Unknown coordinate system "
<< coordSystem << endl
<< "Known types are global and patchLocal"
<< "Known types are global, patchLocal and fieldBased"
<< exit(FatalError);
}
}
......
......@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
......@@ -30,7 +30,8 @@ Description
Used in splitting cells.
Either all cells have similar refinement direction ('global') or
direction is dependent on local cell geometry. Controlled by dictionary.
direction is dependent on local cell geometry, or loads selected fields
by name ('fieldBased'). Controlled by dictionary.
SourceFiles
directions.C
......@@ -65,7 +66,9 @@ class directions
:
public List<vectorField>
{
public:
// Data types
//- Enumeration listing the possible coordinate directions.
......@@ -76,6 +79,7 @@ public:
NORMAL
};
private:
static const NamedEnum<directionType, 3> directionTypeNames_;
......@@ -83,7 +87,6 @@ private:
// Private Member Functions
//- For debugging. Write point coordinate.
static void writeOBJ(Ostream& os, const point& pt);
......@@ -141,7 +144,6 @@ public:
const dictionary& dict,
const twoDPointCorrector* correct2DPtr = NULL
);
};
......
#!/bin/sh
cd ${0%/*} || exit 1 # Run from this directory
# Source tutorial clean functions
. $WM_PROJECT_DIR/bin/tools/CleanFunctions
rm -rf 0 > /dev/null 2>&1
cleanCase
#------------------------------------------------------------------------------
#!/bin/sh
cd ${0%/*} || exit 1 # Run from this directory
# Source tutorial run functions
. $WM_PROJECT_DIR/bin/tools/RunFunctions
wmake calcRadiusField
wclean calcRadiusField
runApplication blockMesh
##### Procedure for special refinement over Z
# We need the 0 folder to exist for these steps
mkdir 0
# Refine over Z, in 6 passes
for index in 1 2 3 4 5 6; do
runApplication calcRadiusField
mv log.calcRadiusField log.calcRadiusField.tier$index
runApplication topoSet -dict system/topoSetDict.tier$index
mv log.topoSet log.topoSet.tier$index
## foamToVTK -cellSet tier$index
runApplication refineMesh -dict system/refineMeshDict.tier$index -overwrite
mv log.refineMesh log.refineMesh.tier$index
rm -r 0/*
done
# Refine over cylindrical coordinates, in 3 passes
for index in 1 2 3; do
runApplication calcRadiusField -calcDirections
mv log.calcRadiusField log.calcRadiusField.range$index
runApplication topoSet -dict system/topoSetDict.range$index
mv log.topoSet log.topoSet.range$index
## foamToVTK -cellSet tier$index
runApplication refineMesh -dict system/refineMeshDict.range$index \
-overwrite
mv log.refineMesh log.refineMesh.range$index
rm -r 0/*
done
#------------------------------------------------------------------------------
calcRadiusField.C
EXE = $(FOAM_USER_APPBIN)/calcRadiusField
EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude
EXE_LIBS = \
-lfiniteVolume
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is derived from 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/>.
Description
Write the volScalarField "radiusFieldXY" that has the distance to the
origin over X,Y.
And also write the direction fields based on the option "-calcDirections".
The resulting fields are:
- radialDirection
- angularDirection
- heightDirection
Derived from:
$FOAM_UTILITIES/postProcessing/miscellaneous/writeCellCentres
\*---------------------------------------------------------------------------*/
#include "argList.H"
#include "timeSelector.H"
#include "Time.H"
#include "fvMesh.H"
#include "vectorIOField.H"
#include "volFields.H"
using namespace Foam;
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
timeSelector::addOptions();
#include "addRegionOption.H"
argList::addBoolOption
(
"calcDirections",
"calculate the direction fields as well"
);
#include "setRootCase.H"
#include "createTime.H"
instantList timeDirs = timeSelector::select0(runTime, args);
const bool calcDirections = args.optionFound("calcDirections");
#include "createNamedMesh.H"
forAll(timeDirs, timeI)
{
runTime.setTime(timeDirs[timeI], timeI);
Info<< "Time = " << runTime.timeName() << endl;
// Check for new mesh
mesh.readUpdate();
Info<< "Writing radius field over X,Y in "
<< runTime.timeName() << endl;
volScalarField radiusFieldXY
(
IOobject
(
"radiusFieldXY",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
sqrt
(
mesh.C().component(0)*mesh.C().component(0)
+ mesh.C().component(1)*mesh.C().component(1)
)
);
radiusFieldXY.write();
if(calcDirections)
{
vectorIOField radialDirection
(
IOobject
(
"radialDirection",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh.C()/magSqr(mesh.C())
);
radialDirection.replace(vector::Z, scalar(0.0));
radialDirection /= sqrt(magSqr(radialDirection));
radialDirection.write();
const tensor transform2Tangencial
(
0, -1, 0,
1, 0, 0,
0, 0, 1
);
vectorIOField angularDirection
(
IOobject
(
"angularDirection",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
transform2Tangencial & mesh.C()
);
angularDirection.replace(vector::Z, scalar(0.0));
angularDirection /= sqrt(magSqr(angularDirection));
angularDirection.write();
vectorIOField heightDirection
(
IOobject
(
"heightDirection",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
radialDirection ^ angularDirection
);
heightDirection.write();
}
}
Info<< "\nEnd\n" << endl;
return 0;
}
// ************************************************************************* //
FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "constant/polyMesh";
object blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
convertToMeters 1;
/*
A ------------------------------- B
/ H /
C\ /
\ /
\ /
\ / G
\ /
\ /
\ /
\/
D
*/
Ax 0.438912; Ay 0.000000;
Bx 18.28800; By 0.000000;
Cx 0.310358; Cy -0.310358;
Dx 12.931569; Dy -12.931569;
Hx 0.405502; Hy -0.167964;
Gx 16.895909; Gy -6.998515;
Z_DB_low -0.88706;
Z_AC_low -1.07;
Z_high 4.39208;
vertices
(
($Cx $Cy $Z_AC_low) //0
($Dx $Dy $Z_DB_low) //1
($Bx $By $Z_DB_low) //2
($Ax $Ay $Z_AC_low) //3
($Cx $Cy $Z_high) //4
($Dx $Dy $Z_high) //5
($Bx $By $Z_high) //6
($Ax $Ay $Z_high) //7
);
blocks
(
hex (0 1 2 3 4 5 6 7) (47 10 4) simpleGrading (41.6669 1 1)
);
edges
(
arc 0 3 ($Hx $Hy $Z_AC_low)
arc 4 7 ($Hx $Hy $Z_high)
arc 1 2 ($Gx $Gy $Z_DB_low)
arc 5 6 ($Gx $Gy $Z_high)
);
patches
(
patch maxX
(
(1 2 6 5)
)
patch minX
(
(0 4 7 3)
)
patch maxY
(
(3 7 6 2)
)
patch minY
(
(1 5 4 0)
)
patch maxZ
(
(4 5 6 7)
)
patch minZ
(
(0 3 2 1)
)
);
mergePatchPairs
(
);
// ************************************************************************* //