Commit 2cecd709 authored by graham's avatar graham
Browse files

New surfaceOffsetLinearDistance cellSizeFunction.

Implemented linearDistance and uniformDistance cellSizeFunctions.

Reading mesh from file as an fvMesh once read, creating a volScalarField to test
the cellSizeFunctions.
parent 254f98f5
......@@ -14,6 +14,7 @@ $(cellSiseFunctions)/cellSizeFunction/cellSizeFunction.C
$(cellSiseFunctions)/uniform/uniform.C
$(cellSiseFunctions)/uniformDistance/uniformDistance.C
$(cellSiseFunctions)/linearDistance/linearDistance.C
$(cellSiseFunctions)/surfaceOffsetLinearDistance/surfaceOffsetLinearDistance.C
initialPointsMethod/initialPointsMethod/initialPointsMethod.C
initialPointsMethod/uniformGrid/uniformGrid.C
......
......@@ -49,17 +49,78 @@ linearDistance::linearDistance
cellSizeFunction(typeName, initialPointsDict, cvMesh, surface),
surfaceCellSize_(readScalar(coeffsDict().lookup("surfaceCellSize"))),
distanceCellSize_(readScalar(coeffsDict().lookup("distanceCellSize"))),
distance_(readScalar(coeffsDict().lookup("distance")))
distance_(readScalar(coeffsDict().lookup("distance"))),
distanceSqr_(sqr(distance_)),
gradient_((distanceCellSize_ - surfaceCellSize_)/distance_)
{}
// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
scalar linearDistance::sizeFunction(scalar d) const
{
return gradient_*d + surfaceCellSize_;
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
bool linearDistance::cellSize(const point& pt, scalar& size) const
{
size = surfaceCellSize_;
size = 0;
List<pointIndexHit> hits;
surface_.findNearest
(
pointField(1, pt),
scalarField(1, distanceSqr_),
hits
);
const pointIndexHit& hitInfo = hits[0];
if (hitInfo.hit())
{
if (sideMode_ == BOTHSIDES)
{
size = sizeFunction(mag(pt - hitInfo.hitPoint()));
return true;
}
pointField ptF(1, pt);
List<searchableSurface::volumeType> vTL;
surface_.getVolumeType(ptF, vTL);
bool functionApplied = false;
if
(
sideMode_ == INSIDE
&& vTL[0] == searchableSurface::INSIDE
)
{
size = sizeFunction(mag(pt - hitInfo.hitPoint()));
functionApplied = true;
}
else if
(
sideMode_ == OUTSIDE
&& vTL[0] == searchableSurface::OUTSIDE
)
{
size = sizeFunction(mag(pt - hitInfo.hitPoint()));
functionApplied = true;
}
return functionApplied;
}
return true;
return false;
}
......
......@@ -64,6 +64,18 @@ private:
//- distance from the surface to control over
scalar distance_;
//- distance squared
scalar distanceSqr_;
//- storing gradient for linear function
scalar gradient_;
// Private Member Functions
//- Calculate the cell size as a function of the given distance
scalar sizeFunction(scalar d) const;
public:
......
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2009-2009 OpenCFD Ltd.
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
\*---------------------------------------------------------------------------*/
#include "surfaceOffsetLinearDistance.H"
#include "addToRunTimeSelectionTable.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
defineTypeNameAndDebug(surfaceOffsetLinearDistance, 0);
addToRunTimeSelectionTable
(
cellSizeFunction,
surfaceOffsetLinearDistance,
dictionary
);
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
surfaceOffsetLinearDistance::surfaceOffsetLinearDistance
(
const dictionary& initialPointsDict,
const conformalVoronoiMesh& cvMesh,
const searchableSurface& surface
)
:
cellSizeFunction(typeName, initialPointsDict, cvMesh, surface),
surfaceCellSize_(readScalar(coeffsDict().lookup("surfaceCellSize"))),
distanceCellSize_(readScalar(coeffsDict().lookup("distanceCellSize"))),
surfaceOffset_(readScalar(coeffsDict().lookup("surfaceOffset"))),
totalDistance_
(
readScalar(coeffsDict().lookup("distance")) + surfaceOffset_
),
totalDistanceSqr_(sqr(totalDistance_)),
gradient_
(
(distanceCellSize_ - surfaceCellSize_)/(totalDistance_ - surfaceOffset_)
),
intercept_(surfaceCellSize_ - gradient_*surfaceOffset_)
{}
// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
scalar surfaceOffsetLinearDistance::sizeFunction(scalar d) const
{
if (d <= surfaceOffset_)
{
return surfaceCellSize_;
}
return gradient_*d + intercept_;
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
bool surfaceOffsetLinearDistance::cellSize(const point& pt, scalar& size) const
{
size = 0;
List<pointIndexHit> hits;
surface_.findNearest
(
pointField(1, pt),
scalarField(1, totalDistanceSqr_),
hits
);
const pointIndexHit& hitInfo = hits[0];
if (hitInfo.hit())
{
if (sideMode_ == BOTHSIDES)
{
size = sizeFunction(mag(pt - hitInfo.hitPoint()));
return true;
}
pointField ptF(1, pt);
List<searchableSurface::volumeType> vTL;
surface_.getVolumeType(ptF, vTL);
bool functionApplied = false;
if
(
sideMode_ == INSIDE
&& vTL[0] == searchableSurface::INSIDE
)
{
size = sizeFunction(mag(pt - hitInfo.hitPoint()));
functionApplied = true;
}
else if
(
sideMode_ == OUTSIDE
&& vTL[0] == searchableSurface::OUTSIDE
)
{
size = sizeFunction(mag(pt - hitInfo.hitPoint()));
functionApplied = true;
}
return functionApplied;
}
return false;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// ************************************************************************* //
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2009-2009 OpenCFD Ltd.
\\/ 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 2 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Class
Foam::surfaceOffsetLinearDistance
Description
SourceFiles
surfaceOffsetLinearDistance.C
\*---------------------------------------------------------------------------*/
#ifndef surfaceOffsetLinearDistance_H
#define surfaceOffsetLinearDistance_H
#include "cellSizeFunction.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class surfaceOffsetLinearDistance Declaration
\*---------------------------------------------------------------------------*/
class surfaceOffsetLinearDistance
:
public cellSizeFunction
{
private:
// Private data
//- cell size at the surface
scalar surfaceCellSize_;
//- cell size at distance_ from the surface
scalar distanceCellSize_;
//- Offset distance from surface for constant size portion
scalar surfaceOffset_;
//- Total distance from the surface to control over (distance +
// surfaceOffset)
scalar totalDistance_;
//- totalDistance squared
scalar totalDistanceSqr_;
//- storing gradient for linear function
scalar gradient_;
//- storing intercept for linear function
scalar intercept_;
// Private Member Functions
//- Calculate the cell size as a function of the given distance
scalar sizeFunction(scalar d) const;
public:
//- Runtime type information
TypeName("surfaceOffsetLinearDistance");
// Constructors
//- Construct from components
surfaceOffsetLinearDistance
(
const dictionary& initialPointsDict,
const conformalVoronoiMesh& cvMesh,
const searchableSurface& surface
);
//- Destructor
virtual ~surfaceOffsetLinearDistance()
{}
// Member Functions
//- Modify scalar argument to the cell size specified by function.
// Return a boolean specifying if the function was used, i.e. false if
// the point was not in range of the surface for a spatially varying
// size.
virtual bool cellSize(const point& pt, scalar& size) const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //
......@@ -55,6 +55,8 @@ uniform::uniform
bool uniform::cellSize(const point& pt, scalar& size) const
{
size = 0;
if (sideMode_ == BOTHSIDES)
{
size = cellSize_;
......@@ -67,14 +69,12 @@ bool uniform::cellSize(const point& pt, scalar& size) const
surface_.getVolumeType(ptF, vTL);
size = 0;
bool functionApplied = false;
if
(
sideMode_ == INSIDE
&& vTL[1] == searchableSurface::INSIDE
&& vTL[0] == searchableSurface::INSIDE
)
{
size = cellSize_;
......@@ -84,7 +84,7 @@ bool uniform::cellSize(const point& pt, scalar& size) const
else if
(
sideMode_ == OUTSIDE
&& vTL[1] == searchableSurface::OUTSIDE
&& vTL[0] == searchableSurface::OUTSIDE
)
{
size = cellSize_;
......
......@@ -48,7 +48,8 @@ uniformDistance::uniformDistance
:
cellSizeFunction(typeName, initialPointsDict, cvMesh, surface),
cellSize_(readScalar(coeffsDict().lookup("cellSize"))),
distance_(readScalar(coeffsDict().lookup("distance")))
distance_(readScalar(coeffsDict().lookup("distance"))),
distanceSqr_(sqr(distance_))
{}
......@@ -56,9 +57,58 @@ uniformDistance::uniformDistance
bool uniformDistance::cellSize(const point& pt, scalar& size) const
{
size = cellSize_;
return true;
size = 0;
List<pointIndexHit> hits;
surface_.findNearest
(
pointField(1, pt),
scalarField(1, distanceSqr_),
hits
);
if (hits[0].hit())
{
if (sideMode_ == BOTHSIDES)
{
size = cellSize_;
return true;
}
pointField ptF(1, pt);
List<searchableSurface::volumeType> vTL;
surface_.getVolumeType(ptF, vTL);
bool functionApplied = false;
if
(
sideMode_ == INSIDE
&& vTL[0] == searchableSurface::INSIDE
)
{
size = cellSize_;
functionApplied = true;
}
else if
(
sideMode_ == OUTSIDE
&& vTL[0] == searchableSurface::OUTSIDE
)
{
size = cellSize_;
functionApplied = true;
}
return functionApplied;
}
return false;
}
......
......@@ -61,6 +61,9 @@ private:
//- Distance
scalar distance_;
//- Distance squared
scalar distanceSqr_;
public:
......
......@@ -32,6 +32,8 @@ License
#include "ulong.H"
#include "surfaceFeatures.H"
#include "zeroGradientPointPatchField.H"
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::conformalVoronoiMesh::conformalVoronoiMesh
......@@ -100,9 +102,6 @@ Foam::conformalVoronoiMesh::conformalVoronoiMesh
{
timeCheck();
point sizeTestPt = vector(0.5, 0.3, 0.6);
Info<< nl << cellSizeControl().cellSize(sizeTestPt) << endl;
createFeaturePoints();
timeCheck();
......@@ -117,6 +116,51 @@ Foam::conformalVoronoiMesh::conformalVoronoiMesh
writeMesh();
timeCheck();
// Test field creation and cell size functions
Info<< "Create fvMesh" << endl;
fvMesh fMesh
(
IOobject
(
Foam::polyMesh::defaultRegion,
runTime_.constant(),
runTime_,
IOobject::MUST_READ
)
);
timeCheck();
Info<< "Create volScalarField" << endl;
volScalarField cellSizeTest
(
IOobject
(
"cellSizeTest",
runTime_.timeName(),
runTime_,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
fMesh,
dimensionedScalar("cellSize", dimLength, 0),
zeroGradientPointPatchField<scalar>::typeName
);
scalarField& cellSize = cellSizeTest.internalField();
const vectorField& C = fMesh.cellCentres();
forAll(cellSize, i)
{
cellSize[i] = cellSizeControl().cellSize(C[i]);
}
cellSizeTest.write();
}
......@@ -128,7 +172,6 @@ Foam::conformalVoronoiMesh::~conformalVoronoiMesh()
// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //