zhyang-dev/OpenFOAM-12
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Merge branch 'master' of github.com-OpenFOAM:OpenFOAM/OpenFOAM-dev

Browse Source
This commit is contained in:
Henry Weller
2022-11-03 14:58:20 +00:00
parent cec0359871 1a7da6391e
commit 392b1fd911
7 changed files with 227 additions and 204 deletions
Download Patch File Download Diff File Expand all files Collapse all files

205
src/OpenFOAM/meshes/primitiveShapes/plane/plane.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2021 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -28,29 +28,37 @@ License
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
void Foam::plane::calcPntAndVec(const scalarList& C)
void Foam::plane::calcPntAndVec
(
const scalar a,
const scalar b,
const scalar c,
const scalar d
)
{
normal_ = vector(C[0], C[1], C[2]);
normal_ = vector(a, b, c);
const scalar magNormal = mag(normal_);
if (magNormal == 0)
// Normalise the normal if possible. Set to invalid if not.
if (magNormal > 0)
{
FatalErrorInFunction
<< "Plane normal has zero length"
<< abort(FatalError);
normal_ /= magNormal;
}
else
{
normal_ = vector::zero;
}
normal_ /= magNormal;
if (magNormal < mag(C[3])*vSmall)
// Construct the point if possible. Set to far away if not.
if (magNormal > mag(d)*vSmall)
{
FatalErrorInFunction
<< "Plane is too far from the origin"
<< abort(FatalError);
point_ = - d/magNormal*normal_;
}
else
{
point_ = point::max;
}
point_ = - C[3]/magNormal*normal_;
}
@ -61,34 +69,9 @@ void Foam::plane::calcPntAndVec
const point& point3
)
{
normal_ = normalised((point1 - point2) ^ (point2 - point3));
point_ = (point1 + point2 + point3)/3;
vector line12 = point1 - point2;
vector line23 = point2 - point3;
if
(
mag(line12) < vSmall
|| mag(line23) < vSmall
|| mag(point3-point1) < vSmall
)
{
FatalErrorInFunction
<< "Bad points:" << point1 << ' ' << point2 << ' ' << point3
<< abort(FatalError);
}
normal_ = line12 ^ line23;
scalar magUnitVector(mag(normal_));
if (magUnitVector < vSmall)
{
FatalErrorInFunction
<< "Plane normal defined with zero length" << nl
<< "Bad points:" << point1 << ' ' << point2 << ' ' << point3
<< abort(FatalError);
}
normal_ /= magUnitVector;
}
@ -96,58 +79,38 @@ void Foam::plane::calcPntAndVec
Foam::plane::plane(const vector& normalVector)
:
normal_(normalVector),
normal_(normalised(normalVector)),
point_(Zero)
{
scalar magUnitVector(mag(normal_));
if (magUnitVector > vSmall)
{
normal_ /= magUnitVector;
}
else
{
FatalErrorInFunction
<< "plane normal has zero length. basePoint:" << point_
<< abort(FatalError);
}
}
{}
Foam::plane::plane(const point& basePoint, const vector& normalVector)
:
normal_(normalVector),
normal_(normalised(normalVector)),
point_(basePoint)
{}
Foam::plane::plane
(
const scalar a,
const scalar b,
const scalar c,
const scalar d
)
{
scalar magUnitVector(mag(normal_));
if (magUnitVector > vSmall)
{
normal_ /= magUnitVector;
}
else
{
FatalErrorInFunction
<< "plane normal has zero length. basePoint:" << point_
<< abort(FatalError);
}
}
Foam::plane::plane(const scalarList& C)
{
calcPntAndVec(C);
calcPntAndVec(a, b, c, d);
}
Foam::plane::plane
(
const point& a,
const point& b,
const point& c
const point& point1,
const point& point2,
const point& point3
)
{
calcPntAndVec(a, b, c);
calcPntAndVec(point1, point2, point3);
}
@ -158,38 +121,35 @@ Foam::plane::plane(const dictionary& dict)
{
const word planeType(dict.lookup("planeType"));
const dictionary& subDict = dict.optionalSubDict(planeType + "Dict");
if (planeType == "planeEquation")
{
const dictionary& subDict = dict.optionalSubDict("planeEquationDict");
scalarList C(4);
C[0] = subDict.lookup<scalar>("a");
C[1] = subDict.lookup<scalar>("b");
C[2] = subDict.lookup<scalar>("c");
C[3] = subDict.lookup<scalar>("d");
calcPntAndVec(C);
calcPntAndVec
(
subDict.lookup<scalar>("a"),
subDict.lookup<scalar>("b"),
subDict.lookup<scalar>("c"),
subDict.lookup<scalar>("d")
);
}
else if (planeType == "embeddedPoints")
{
const dictionary& subDict = dict.optionalSubDict("embeddedPointsDict");
point point1(subDict.lookup("point1"));
point point2(subDict.lookup("point2"));
point point3(subDict.lookup("point3"));
calcPntAndVec(point1, point2, point3);
calcPntAndVec
(
subDict.lookup<point>("point1"),
subDict.lookup<point>("point2"),
subDict.lookup<point>("point3")
);
}
else if (planeType == "pointAndNormal")
{
const dictionary& subDict = dict.optionalSubDict("pointAndNormalDict");
point_ =
subDict.lookupBackwardsCompatible<point>
(
{"point", "basePoint"}
);
normal_ =
normalised
(
@ -203,47 +163,36 @@ Foam::plane::plane(const dictionary& dict)
{
FatalIOErrorInFunction(dict)
<< "Invalid plane type: " << planeType << nl
<< "Valid options include: planeEquation, embeddedPoints and "
<< "pointAndNormal"
<< "Valid options include: "
<< "planeEquation, embeddedPoints and pointAndNormal"
<< abort(FatalIOError);
}
if (normal_ == vector::zero)
{
FatalIOErrorInFunction(subDict)
<< "Plane normal has zero length"
<< exit(FatalIOError);
}
if (point_ == point::max)
{
FatalIOErrorInFunction(subDict)
<< "Plane is too far from the origin"
<< exit(FatalIOError);
}
}
Foam::plane::plane(Istream& is)
:
normal_(is),
normal_(normalised(vector(is))),
point_(is)
{
scalar magUnitVector(mag(normal_));
if (magUnitVector > vSmall)
{
normal_ /= magUnitVector;
}
else
{
FatalErrorInFunction
<< "plane normal has zero length. basePoint:" << point_
<< abort(FatalError);
}
}
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
const Foam::vector& Foam::plane::normal() const
{
return normal_;
}
const Foam::point& Foam::plane::refPoint() const
{
return point_;
}
Foam::FixedList<Foam::scalar, 4> Foam::plane::planeCoeffs() const
{
FixedList<scalar, 4> C(4);

53
src/OpenFOAM/meshes/primitiveShapes/plane/plane.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2020 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -73,6 +73,7 @@ public:
class ray
{
point refPoint_;
vector dir_;
public:
@ -109,10 +110,15 @@ private:
// Private Member Functions
//- Calculates basePoint and normal vector given plane coefficients
void calcPntAndVec(const scalarList& C);
void calcPntAndVec
(
const scalar C1,
const scalar C2,
const scalar C3,
const scalar C4
);
//- Calculates basePoint and normal vector given three points
//- Normal vector determined using right hand rule
void calcPntAndVec
(
const point& point1,
@ -139,27 +145,44 @@ public:
const point& point3
);
//- Construct from coefficients for the
// plane equation: ax + by + cz + d = 0
explicit plane(const scalarList& C);
//- Construct from coefficients for the plane equation,
// ax + by + cz + d = 0
explicit plane
(
const scalar a,
const scalar b,
const scalar c,
const scalar d
);
//- Construct from dictionary
explicit plane(const dictionary& planeDict);
//- Construct from Istream. Assumes the base + normal notation.
//- Construct from Istream
explicit plane(Istream& is);
// Member Functions
//- Return plane normal
const vector& normal() const;
//- Return whether or not this plane is valid
inline bool valid() const
{
return normal_ != vector::zero && point_ != point::max;
}
//- Return or return plane base point
const point& refPoint() const;
//- Return the plane normal
inline const vector& normal() const
{
return normal_;
}
//- Return coefficients for the
// plane equation: ax + by + cz + d = 0
//- Return the plane base point
inline const point& refPoint() const
{
return point_;
}
//- Return coefficients for the plane equation, ax + by + cz + d = 0
FixedList<scalar, 4> planeCoeffs() const;
//- Return a point on the plane
@ -176,7 +199,7 @@ public:
scalar normalIntersect(const point& pnt0, const vector& dir) const;
//- Return cut coefficient for plane and ray
scalar normalIntersect(const ray& r) const
inline scalar normalIntersect(const ray& r) const
{
return normalIntersect(r.refPoint(), r.dir());
}
@ -184,7 +207,7 @@ public:
//- Return the cutting point between the plane and
// a line passing through the supplied points
template<class Point, class PointRef>
scalar lineIntersect(const line<Point, PointRef>& l) const
inline scalar lineIntersect(const line<Point, PointRef>& l) const
{
return normalIntersect(l.start(), l.vec());
}

63
src/finiteVolume/finiteVolume/fvc/fvcDdt.C
View File

@ -102,6 +102,18 @@ ddt
}
template<class Type>
tmp<GeometricField<Type, fvPatchField, volMesh>>
ddt
(
const one&,
const GeometricField<Type, fvPatchField, volMesh>& vf
)
{
return ddt(vf);
}
template<class Type>
tmp<GeometricField<Type, fvPatchField, volMesh>>
ddt
@ -125,25 +137,11 @@ ddt
}
template<class Type>
tmp<GeometricField<Type, fvsPatchField, surfaceMesh>>
ddt
(
const GeometricField<Type, fvsPatchField, surfaceMesh>& sf
)
{
return fv::ddtScheme<Type>::New
(
sf.mesh(),
sf.mesh().schemes().ddt("ddt(" + sf.name() + ')')
).ref().fvcDdt(sf);
}
template<class Type>
tmp<GeometricField<Type, fvPatchField, volMesh>>
ddt
(
const one&,
const one&,
const GeometricField<Type, fvPatchField, volMesh>& vf
)
@ -156,11 +154,40 @@ template<class Type>
tmp<GeometricField<Type, fvPatchField, volMesh>>
ddt
(
const GeometricField<Type, fvPatchField, volMesh>& vf,
const one&
const one&,
const volScalarField& rho,
const GeometricField<Type, fvPatchField, volMesh>& vf
)
{
return ddt(vf);
return ddt(rho, vf);
}
template<class Type>
tmp<GeometricField<Type, fvPatchField, volMesh>>
ddt
(
const volScalarField& alpha,
const one&,
const GeometricField<Type, fvPatchField, volMesh>& vf
)
{
return ddt(alpha, vf);
}
template<class Type>
tmp<GeometricField<Type, fvsPatchField, surfaceMesh>>
ddt
(
const GeometricField<Type, fvsPatchField, surfaceMesh>& sf
)
{
return fv::ddtScheme<Type>::New
(
sf.mesh(),
sf.mesh().schemes().ddt("ddt(" + sf.name() + ')')
).ref().fvcDdt(sf);
}

56
src/finiteVolume/finiteVolume/fvc/fvcDdt.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2020 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -80,6 +80,13 @@ namespace fvc
const GeometricField<Type, fvPatchField, volMesh>&
);
template<class Type>
tmp<GeometricField<Type, fvPatchField, volMesh>> ddt
(
const one&,
const GeometricField<Type, fvPatchField, volMesh>&
);
template<class Type>
tmp<GeometricField<Type, fvPatchField, volMesh>> ddt
(
@ -88,35 +95,36 @@ namespace fvc
const GeometricField<Type, fvPatchField, volMesh>&
);
template<class Type>
tmp<GeometricField<Type, fvPatchField, volMesh>> ddt
(
const one&,
const one&,
const GeometricField<Type, fvPatchField, volMesh>&
);
template<class Type>
tmp<GeometricField<Type, fvPatchField, volMesh>> ddt
(
const one&,
const volScalarField&,
const GeometricField<Type, fvPatchField, volMesh>&
);
template<class Type>
tmp<GeometricField<Type, fvPatchField, volMesh>> ddt
(
const volScalarField&,
const one&,
const GeometricField<Type, fvPatchField, volMesh>&
);
template<class Type>
tmp<GeometricField<Type, fvsPatchField, surfaceMesh>> ddt
(
const GeometricField<Type, fvsPatchField, surfaceMesh>&
);
template<class Type>
tmp<GeometricField<Type, fvPatchField, volMesh>> ddt
(
const one&,
const GeometricField<Type, fvPatchField, volMesh>&
);
template<class Type>
tmp<GeometricField<Type, fvPatchField, volMesh>> ddt
(
const GeometricField<Type, fvPatchField, volMesh>&,
const one&
);
inline geometricZeroField ddt
(
const one&,
const one&
)
{
return geometricZeroField();
}
template<class Type>
tmp
<

24
src/finiteVolume/finiteVolume/fvm/fvmDdt.C
View File

@ -55,18 +55,6 @@ ddt
}
template<class Type>
tmp<fvMatrix<Type>>
ddt
(
const one&,
const GeometricField<Type, fvPatchField, volMesh>& vf
)
{
return ddt(vf);
}
template<class Type>
tmp<fvMatrix<Type>>
ddt
@ -99,6 +87,18 @@ ddt
}
template<class Type>
tmp<fvMatrix<Type>>
ddt
(
const one&,
const GeometricField<Type, fvPatchField, volMesh>& vf
)
{
return ddt(vf);
}
template<class Type>
tmp<fvMatrix<Type>>
ddt

14
src/finiteVolume/finiteVolume/fvm/fvmDdt.H
View File

@ -57,13 +57,6 @@ namespace fvm
const GeometricField<Type, fvPatchField, volMesh>&
);
template<class Type>
tmp<fvMatrix<Type>> ddt
(
const one&,
const GeometricField<Type, fvPatchField, volMesh>&
);
template<class Type>
tmp<fvMatrix<Type>> ddt
(
@ -78,6 +71,13 @@ namespace fvm
const GeometricField<Type, fvPatchField, volMesh>&
);
template<class Type>
tmp<fvMatrix<Type>> ddt
(
const one&,
const GeometricField<Type, fvPatchField, volMesh>&
);
template<class Type>
tmp<fvMatrix<Type>> ddt
(

16
src/meshTools/tetOverlapVolume/tetOverlapVolume.C
View File

@ -61,6 +61,10 @@ Foam::scalar Foam::tetOverlapVolume::tetTetOverlapVol
// face 0
const plane pl0(tetB.b(), tetB.d(), tetB.c());
if (!pl0.valid())
{
return 0;
}
const FixedList<point, 4> t({tetA.a(), tetA.b(), tetA.c(), tetA.d()});
cutTetList1.clear();
tetCut(t, pl0, cut::appendOp<tetListType>(cutTetList1), cut::noOp());
@ -71,6 +75,10 @@ Foam::scalar Foam::tetOverlapVolume::tetTetOverlapVol
// face 1
const plane pl1(tetB.a(), tetB.c(), tetB.d());
if (!pl1.valid())
{
return 0;
}
cutTetList2.clear();
for (label i = 0; i < cutTetList1.size(); i++)
{
@ -84,6 +92,10 @@ Foam::scalar Foam::tetOverlapVolume::tetTetOverlapVol
// face 2
const plane pl2(tetB.a(), tetB.d(), tetB.b());
if (!pl2.valid())
{
return 0;
}
cutTetList1.clear();
for (label i = 0; i < cutTetList2.size(); i++)
{
@ -97,6 +109,10 @@ Foam::scalar Foam::tetOverlapVolume::tetTetOverlapVol
// face 3
const plane pl3(tetB.a(), tetB.b(), tetB.c());
if (!pl3.valid())
{
return 0;
}
scalar v = 0;
for (label i = 0; i < cutTetList1.size(); i++)
{