Skip to content
Snippets Groups Projects
Commit 709fc12c authored by andy's avatar andy
Browse files

ENH: Robustness updates to pointMVC interpolation

parent 1e77abbe
Branches
Tags
No related merge requests found
Showing with 45 additions and 44 deletions
...@@ -204,8 +204,8 @@ $(interpolation)/interpolationCellPoint/makeInterpolationCellPoint.C ...@@ -204,8 +204,8 @@ $(interpolation)/interpolationCellPoint/makeInterpolationCellPoint.C
$(interpolation)/interpolationCellPointFace/makeInterpolationCellPointFace.C $(interpolation)/interpolationCellPointFace/makeInterpolationCellPointFace.C
$(interpolation)/interpolationCellPointWallModified/cellPointWeightWallModified/cellPointWeightWallModified.C $(interpolation)/interpolationCellPointWallModified/cellPointWeightWallModified/cellPointWeightWallModified.C
$(interpolation)/interpolationCellPointWallModified/makeInterpolationCellPointWallModified.C $(interpolation)/interpolationCellPointWallModified/makeInterpolationCellPointWallModified.C
$(interpolation)/interpolationPoint/pointMVCWeight.C $(interpolation)/interpolationPointMVC/pointMVCWeight.C
$(interpolation)/interpolationPoint/makeInterpolationPoint.C $(interpolation)/interpolationPointMVC/makeInterpolationPointMVC.C
volPointInterpolation = interpolation/volPointInterpolation volPointInterpolation = interpolation/volPointInterpolation
$(volPointInterpolation)/volPointInterpolation.C $(volPointInterpolation)/volPointInterpolation.C
......
EXE_INC = \ EXE_INC = \
-DFULLDEBUG -g -O0 \
-I$(LIB_SRC)/triSurface/lnInclude \ -I$(LIB_SRC)/triSurface/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude -I$(LIB_SRC)/meshTools/lnInclude
......
...@@ -2,7 +2,7 @@ ...@@ -2,7 +2,7 @@
========= | ========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | \\ / O peration |
\\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd. \\ / A nd | Copyright (C) 2010-2011 OpenCFD Ltd.
\\/ M anipulation | \\/ M anipulation |
------------------------------------------------------------------------------- -------------------------------------------------------------------------------
License License
...@@ -23,13 +23,13 @@ License ...@@ -23,13 +23,13 @@ License
\*---------------------------------------------------------------------------*/ \*---------------------------------------------------------------------------*/
#include "interpolationPoint.H" #include "interpolationPointMVC.H"
#include "volPointInterpolation.H" #include "volPointInterpolation.H"
// * * * * * * * * * * * * * * * * Constructor * * * * * * * * * * * * * * * // // * * * * * * * * * * * * * * * * Constructor * * * * * * * * * * * * * * * //
template<class Type> template<class Type>
Foam::interpolationPoint<Type>::interpolationPoint Foam::interpolationPointMVC<Type>::interpolationPointMVC
( (
const GeometricField<Type, fvPatchField, volMesh>& psi const GeometricField<Type, fvPatchField, volMesh>& psi
) )
......
...@@ -2,7 +2,7 @@ ...@@ -2,7 +2,7 @@
========= | ========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | \\ / O peration |
\\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd. \\ / A nd | Copyright (C) 2010-2011 OpenCFD Ltd.
\\/ M anipulation | \\/ M anipulation |
------------------------------------------------------------------------------- -------------------------------------------------------------------------------
License License
...@@ -22,7 +22,7 @@ License ...@@ -22,7 +22,7 @@ License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>. along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Class Class
Foam::interpolationPoint Foam::interpolationPointMVC
Description Description
Given cell centre values interpolates to vertices and uses these to Given cell centre values interpolates to vertices and uses these to
...@@ -30,8 +30,8 @@ Description ...@@ -30,8 +30,8 @@ Description
\*---------------------------------------------------------------------------*/ \*---------------------------------------------------------------------------*/
#ifndef interpolationPoint_H #ifndef interpolationPointMVC_H
#define interpolationPoint_H #define interpolationPointMVC_H
#include "interpolation.H" #include "interpolation.H"
#include "pointMVCWeight.H" #include "pointMVCWeight.H"
...@@ -42,11 +42,11 @@ namespace Foam ...@@ -42,11 +42,11 @@ namespace Foam
{ {
/*---------------------------------------------------------------------------*\ /*---------------------------------------------------------------------------*\
Class interpolationPoint Declaration Class interpolationPointMVC Declaration
\*---------------------------------------------------------------------------*/ \*---------------------------------------------------------------------------*/
template<class Type> template<class Type>
class interpolationPoint class interpolationPointMVC
: :
public interpolation<Type> public interpolation<Type>
{ {
...@@ -57,16 +57,17 @@ protected: ...@@ -57,16 +57,17 @@ protected:
//- Interpolated volfield //- Interpolated volfield
const GeometricField<Type, pointPatchField, pointMesh> psip_; const GeometricField<Type, pointPatchField, pointMesh> psip_;
public: public:
//- Runtime type information //- Runtime type information
TypeName("point"); TypeName("pointMVC");
// Constructors // Constructors
//- Construct from components //- Construct from components
interpolationPoint interpolationPointMVC
( (
const GeometricField<Type, fvPatchField, volMesh>& psi const GeometricField<Type, fvPatchField, volMesh>& psi
); );
...@@ -93,12 +94,12 @@ public: ...@@ -93,12 +94,12 @@ public:
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#include "interpolationPointI.H" #include "interpolationPointMVCI.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository #ifdef NoRepository
# include "interpolationPoint.C" # include "interpolationPointMVC.C"
#endif #endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
......
...@@ -2,7 +2,7 @@ ...@@ -2,7 +2,7 @@
========= | ========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | \\ / O peration |
\\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd. \\ / A nd | Copyright (C) 2010-2011 OpenCFD Ltd.
\\/ M anipulation | \\/ M anipulation |
------------------------------------------------------------------------------- -------------------------------------------------------------------------------
License License
...@@ -26,7 +26,7 @@ License ...@@ -26,7 +26,7 @@ License
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class Type> template<class Type>
inline Type Foam::interpolationPoint<Type>::interpolate inline Type Foam::interpolationPointMVC<Type>::interpolate
( (
const pointMVCWeight& cpw const pointMVCWeight& cpw
) const ) const
...@@ -36,7 +36,7 @@ inline Type Foam::interpolationPoint<Type>::interpolate ...@@ -36,7 +36,7 @@ inline Type Foam::interpolationPoint<Type>::interpolate
template<class Type> template<class Type>
inline Type Foam::interpolationPoint<Type>::interpolate inline Type Foam::interpolationPointMVC<Type>::interpolate
( (
const vector& position, const vector& position,
const label cellI, const label cellI,
......
...@@ -2,7 +2,7 @@ ...@@ -2,7 +2,7 @@
========= | ========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | \\ / O peration |
\\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd. \\ / A nd | Copyright (C) 2010-2011 OpenCFD Ltd.
\\/ M anipulation | \\/ M anipulation |
------------------------------------------------------------------------------- -------------------------------------------------------------------------------
License License
...@@ -23,13 +23,13 @@ License ...@@ -23,13 +23,13 @@ License
\*---------------------------------------------------------------------------*/ \*---------------------------------------------------------------------------*/
#include "interpolationPoint.H" #include "interpolationPointMVC.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam namespace Foam
{ {
makeInterpolation(interpolationPoint); makeInterpolation(interpolationPointMVC);
} }
// ************************************************************************* // // ************************************************************************* //
...@@ -2,7 +2,7 @@ ...@@ -2,7 +2,7 @@
========= | ========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | \\ / O peration |
\\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd. \\ / A nd | Copyright (C) 2010-2011 OpenCFD Ltd.
\\/ M anipulation | \\/ M anipulation |
------------------------------------------------------------------------------- -------------------------------------------------------------------------------
License License
...@@ -56,7 +56,7 @@ void Foam::pointMVCWeight::calcWeights ...@@ -56,7 +56,7 @@ void Foam::pointMVCWeight::calcWeights
forAll(f, j) forAll(f, j)
{ {
label jPlus1 = f.fcIndex(j); label jPlus1 = f.fcIndex(j);
scalar l = mag(u[j]-u[jPlus1]); scalar l = mag(u[j] - u[jPlus1]);
theta[j] = 2.0*Foam::asin(l/2.0); theta[j] = 2.0*Foam::asin(l/2.0);
} }
...@@ -68,7 +68,7 @@ void Foam::pointMVCWeight::calcWeights ...@@ -68,7 +68,7 @@ void Foam::pointMVCWeight::calcWeights
weights[pid] = weights[pid] =
1.0 1.0
/ dist[pid] / dist[pid]
* (Foam::tan(theta[jMin1]/2.0)+Foam::tan(theta[j]/2.0)); * (Foam::tan(theta[jMin1]/2.0) + Foam::tan(theta[j]/2.0));
sumWeight += weights[pid]; sumWeight += weights[pid];
} }
...@@ -132,11 +132,10 @@ void Foam::pointMVCWeight::calcWeights ...@@ -132,11 +132,10 @@ void Foam::pointMVCWeight::calcWeights
//Pout<< " uj:" << u[j] << " ujPlus1:" << u[jPlus1] //Pout<< " uj:" << u[j] << " ujPlus1:" << u[jPlus1]
// << " temp:" << temp << endl; // << " temp:" << temp << endl;
scalar l = mag(u[j]-u[jPlus1]); scalar l = min(mag(u[j] - u[jPlus1]), 2.0);
scalar angle = 2.0*Foam::asin(l/2.0); scalar angle = 2.0*Foam::asin(l/2.0);
//Pout<< " j:" << j << " l:" << l //Pout<< " j:" << j << " l:" << l << " angle:" << angle << endl;
// << " angle:" << angle << endl;
v += 0.5*angle*temp; v += 0.5*angle*temp;
} }
...@@ -145,14 +144,14 @@ void Foam::pointMVCWeight::calcWeights ...@@ -145,14 +144,14 @@ void Foam::pointMVCWeight::calcWeights
v /= vNorm; v /= vNorm;
// Make sure v points towards the polygon // Make sure v points towards the polygon
//if (((v&u[0]) < 0) != (mesh.faceOwner()[faceI] != cellIndex_)) //if (((v&u[0]) < 0) != (mesh.faceOwner()[faceI] != cellIndex_))
//{ //{
// FatalErrorIn("pointMVCWeight::calcWeights(..)") // FatalErrorIn("pointMVCWeight::calcWeights(..)")
// << "v:" << v << " u[0]:" << u[0] // << "v:" << v << " u[0]:" << u[0]
// << exit(FatalError); // << exit(FatalError);
//} //}
if ((v&u[0]) < 0) if ((v & u[0]) < 0)
{ {
v = -v; v = -v;
} }
...@@ -165,22 +164,22 @@ void Foam::pointMVCWeight::calcWeights ...@@ -165,22 +164,22 @@ void Foam::pointMVCWeight::calcWeights
label jPlus1 = f.fcIndex(j); label jPlus1 = f.fcIndex(j);
//Pout<< " uj:" << u[j] << " ujPlus1:" << u[jPlus1] << endl; //Pout<< " uj:" << u[j] << " ujPlus1:" << u[jPlus1] << endl;
vector n0 = u[j] ^ v; vector n0 = u[j]^v;
n0 /= mag(n0); n0 /= mag(n0);
vector n1 = u[jPlus1] ^ v; vector n1 = u[jPlus1]^v;
n1 /= mag(n1); n1 /= mag(n1);
scalar l = mag(n0-n1); scalar l = min(mag(n0 - n1), 2.0);
//Pout<< " l:" << l << endl; //Pout<< " l:" << l << endl;
alpha(j) = 2.0*Foam::asin(l/2.0); alpha(j) = 2.0*Foam::asin(l/2.0);
vector temp = n0 ^ n1; vector temp = n0^n1;
if ((temp&v) < 0.0) if ((temp&v) < 0.0)
{ {
alpha[j] = -alpha[j]; alpha[j] = -alpha[j];
} }
l = mag(u[j]-v); l = min(mag(u[j] - v), 2.0);
//Pout<< " l:" << l << endl; //Pout<< " l:" << l << endl;
theta(j) = 2.0*Foam::asin(l/2.0); theta(j) = 2.0*Foam::asin(l/2.0);
} }
...@@ -211,7 +210,7 @@ void Foam::pointMVCWeight::calcWeights ...@@ -211,7 +210,7 @@ void Foam::pointMVCWeight::calcWeights
sum += sum +=
1.0 1.0
/ Foam::tan(theta[j]) / Foam::tan(theta[j])
* (Foam::tan(alpha[j]/2.0)+Foam::tan(alpha[jMin1]/2.0)); * (Foam::tan(alpha[j]/2.0) + Foam::tan(alpha[jMin1]/2.0));
} }
// The special case when x lies on the polygon, handle it using 2D mvc. // The special case when x lies on the polygon, handle it using 2D mvc.
...@@ -234,7 +233,7 @@ void Foam::pointMVCWeight::calcWeights ...@@ -234,7 +233,7 @@ void Foam::pointMVCWeight::calcWeights
/ sum / sum
/ dist[pid] / dist[pid]
/ Foam::sin(theta[j]) / Foam::sin(theta[j])
* (Foam::tan(alpha[j]/2.0)+Foam::tan(alpha[jMin1]/2.0)); * (Foam::tan(alpha[j]/2.0) + Foam::tan(alpha[jMin1]/2.0));
} }
} }
......
...@@ -2,7 +2,7 @@ ...@@ -2,7 +2,7 @@
========= | ========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | \\ / O peration |
\\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd. \\ / A nd | Copyright (C) 2010-2011 OpenCFD Ltd.
\\/ M anipulation | \\/ M anipulation |
------------------------------------------------------------------------------- -------------------------------------------------------------------------------
License License
......
...@@ -2,7 +2,7 @@ ...@@ -2,7 +2,7 @@
========= | ========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | \\ / O peration |
\\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd. \\ / A nd | Copyright (C) 2010-2011 OpenCFD Ltd.
\\/ M anipulation | \\/ M anipulation |
------------------------------------------------------------------------------- -------------------------------------------------------------------------------
License License
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment