From 43dd149577a4f4ec321a6d1e298b29306481b44c Mon Sep 17 00:00:00 2001
From: andy <andy>
Date: Mon, 13 Jan 2014 14:54:57 +0000
Subject: [PATCH] ENH: cyclicAMI - enable optional setting of rotation angle
---
.../cyclicAMIPolyPatch/cyclicAMIPolyPatch.C | 169 ++++++++++++++----
.../cyclicAMIPolyPatch/cyclicAMIPolyPatch.H | 12 +-
2 files changed, 143 insertions(+), 38 deletions(-)
diff --git a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C
index 52373c438d0..98ce25d2600 100644
--- a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C
+++ b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2014 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -45,7 +45,7 @@ namespace Foam
// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
-Foam::vector Foam::cyclicAMIPolyPatch::findFaceMaxRadius
+Foam::vector Foam::cyclicAMIPolyPatch::findFaceNormalMaxRadius
(
const pointField& faceCentres
) const
@@ -138,45 +138,115 @@ void Foam::cyclicAMIPolyPatch::calcTransforms
{
case ROTATIONAL:
{
- point n0 = vector::zero;
- if (half0Ctrs.size())
- {
- n0 = findFaceMaxRadius(half0Ctrs);
- }
- point n1 = vector::zero;
- if (half1Ctrs.size())
- {
- n1 = -findFaceMaxRadius(half1Ctrs);
- }
+ tensor revT = tensor::zero;
- reduce(n0, maxMagSqrOp<point>());
- reduce(n1, maxMagSqrOp<point>());
+ if (rotationAngleDefined_)
+ {
+ tensor T(rotationAxis_*rotationAxis_);
- n0 /= mag(n0) + VSMALL;
- n1 /= mag(n1) + VSMALL;
+ tensor S
+ (
+ 0, -rotationAxis_.z(), rotationAxis_.y(),
+ rotationAxis_.z(), 0, -rotationAxis_.x(),
+ -rotationAxis_.y(), rotationAxis_.x(), 0
+ );
- if (debug)
+ tensor RPos
+ (
+ T
+ + cos(rotationAngle_)*(tensor::I + T)
+ + sin(rotationAngle_)*S
+ );
+ tensor RNeg
+ (
+ T
+ + cos(-rotationAngle_)*(tensor::I + T)
+ + sin(-rotationAngle_)*S
+ );
+
+ // check - assume correct angle when difference in face areas
+ // is the smallest
+ vector transformedAreaPos = sum(half0Areas & RPos);
+ vector transformedAreaNeg = sum(half0Areas & RNeg);
+ vector area1 = sum(half1Areas);
+ reduce(transformedAreaPos, sumOp<vector>());
+ reduce(transformedAreaNeg, sumOp<vector>());
+ reduce(area1, sumOp<vector>());
+
+ scalar errorPos = mag(transformedAreaPos - area1);
+ scalar errorNeg = mag(transformedAreaNeg - area1);
+
+ if (errorPos < errorNeg)
+ {
+ revT = RPos;
+ }
+ else
+ {
+ revT = RNeg;
+ rotationAngle_ *= -1;
+ }
+
+ if (debug)
+ {
+ scalar theta = radToDeg(rotationAngle_);
+
+ Pout<< "cyclicAMIPolyPatch::calcTransforms: patch:"
+ << name()
+ << " Specified rotation:"
+ << " swept angle: " << theta << " [deg]"
+ << " reverse transform: " << revT
+ << endl;
+ }
+ }
+ else
{
- Pout<< "cyclicAMIPolyPatch::calcTransforms : patch:" << name()
- << " Specified rotation :"
- << " n0:" << n0 << " n1:" << n1 << endl;
+ point n0 = vector::zero;
+ point n1 = vector::zero;
+ if (half0Ctrs.size())
+ {
+ n0 = findFaceNormalMaxRadius(half0Ctrs);
+ }
+ if (half1Ctrs.size())
+ {
+ n1 = -findFaceNormalMaxRadius(half1Ctrs);
+ }
+
+ reduce(n0, maxMagSqrOp<point>());
+ reduce(n1, maxMagSqrOp<point>());
+
+ n0 /= mag(n0) + VSMALL;
+ n1 /= mag(n1) + VSMALL;
+
+ // Extended tensor from two local coordinate systems calculated
+ // using normal and rotation axis
+ const tensor E0
+ (
+ rotationAxis_,
+ (n0 ^ rotationAxis_),
+ n0
+ );
+ const tensor E1
+ (
+ rotationAxis_,
+ (-n1 ^ rotationAxis_),
+ -n1
+ );
+ revT = E1.T() & E0;
+
+ if (debug)
+ {
+ scalar theta = radToDeg(acos(n0 & n1));
+
+ Pout<< "cyclicAMIPolyPatch::calcTransforms: patch:"
+ << name()
+ << " Specified rotation:"
+ << " n0:" << n0 << " n1:" << n1
+ << " swept angle: " << theta << " [deg]"
+ << " reverse transform: " << revT
+ << endl;
+ }
}
- // Extended tensor from two local coordinate systems calculated
- // using normal and rotation axis
- const tensor E0
- (
- rotationAxis_,
- (n0 ^ rotationAxis_),
- n0
- );
- const tensor E1
- (
- rotationAxis_,
- (-n1 ^ rotationAxis_),
- -n1
- );
- const tensor revT(E1.T() & E0);
const_cast<tensorField&>(forwardT()) = tensorField(1, revT.T());
const_cast<tensorField&>(reverseT()) = tensorField(1, revT);
const_cast<vectorField&>(separation()).setSize(0);
@@ -395,6 +465,8 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
nbrPatchID_(-1),
rotationAxis_(vector::zero),
rotationCentre_(point::zero),
+ rotationAngleDefined_(false),
+ rotationAngle_(0.0),
separationVector_(vector::zero),
AMIPtr_(NULL),
AMIReverse_(false),
@@ -421,6 +493,8 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
nbrPatchID_(-1),
rotationAxis_(vector::zero),
rotationCentre_(point::zero),
+ rotationAngleDefined_(false),
+ rotationAngle_(0.0),
separationVector_(vector::zero),
AMIPtr_(NULL),
AMIReverse_(dict.lookupOrDefault<bool>("flipNormals", false)),
@@ -466,6 +540,17 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
{
dict.lookup("rotationAxis") >> rotationAxis_;
dict.lookup("rotationCentre") >> rotationCentre_;
+ if (dict.readIfPresent("rotationAngle", rotationAngle_))
+ {
+ rotationAngleDefined_ = true;
+ rotationAngle_ = degToRad(rotationAngle_);
+
+ if (debug)
+ {
+ Info<< "rotationAngle: " << rotationAngle_ << " [rad]"
+ << endl;
+ }
+ }
scalar magRot = mag(rotationAxis_);
if (magRot < SMALL)
@@ -516,6 +601,8 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
nbrPatchID_(-1),
rotationAxis_(pp.rotationAxis_),
rotationCentre_(pp.rotationCentre_),
+ rotationAngleDefined_(pp.rotationAngleDefined_),
+ rotationAngle_(pp.rotationAngle_),
separationVector_(pp.separationVector_),
AMIPtr_(NULL),
AMIReverse_(pp.AMIReverse_),
@@ -543,6 +630,8 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
nbrPatchID_(-1),
rotationAxis_(pp.rotationAxis_),
rotationCentre_(pp.rotationCentre_),
+ rotationAngleDefined_(pp.rotationAngleDefined_),
+ rotationAngle_(pp.rotationAngle_),
separationVector_(pp.separationVector_),
AMIPtr_(NULL),
AMIReverse_(pp.AMIReverse_),
@@ -584,6 +673,8 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
nbrPatchID_(-1),
rotationAxis_(pp.rotationAxis_),
rotationCentre_(pp.rotationCentre_),
+ rotationAngleDefined_(pp.rotationAngleDefined_),
+ rotationAngle_(pp.rotationAngle_),
separationVector_(pp.separationVector_),
AMIPtr_(NULL),
AMIReverse_(pp.AMIReverse_),
@@ -874,6 +965,14 @@ void Foam::cyclicAMIPolyPatch::write(Ostream& os) const
<< token::END_STATEMENT << nl;
os.writeKeyword("rotationCentre") << rotationCentre_
<< token::END_STATEMENT << nl;
+
+ if (rotationAngleDefined_)
+ {
+ os.writeKeyword("rotationAngle") << radToDeg(rotationAngle_)
+ << token::END_STATEMENT << nl;
+
+ }
+
break;
}
case TRANSLATIONAL:
diff --git a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.H b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.H
index e3bb1d35d97..fbdffb3deee 100644
--- a/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.H
+++ b/src/meshTools/AMIInterpolation/patches/cyclicAMI/cyclicAMIPolyPatch/cyclicAMIPolyPatch.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2014 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -75,9 +75,15 @@ private:
//- Axis of rotation for rotational cyclics
vector rotationAxis_;
- //- point on axis of rotation for rotational cyclics
+ //- Point on axis of rotation for rotational cyclics
point rotationCentre_;
+ //- Flag to show whether the rotation angle is defined
+ bool rotationAngleDefined_;
+
+ //- Rotation angle
+ scalar rotationAngle_;
+
// For translation
@@ -101,7 +107,7 @@ private:
// Private Member Functions
//- Return normal of face at max distance from rotation axis
- vector findFaceMaxRadius(const pointField& faceCentres) const;
+ vector findFaceNormalMaxRadius(const pointField& faceCentres) const;
void calcTransforms
(
--
GitLab