applications/utilities/mesh/generation/blockMesh/curvedEdges/arcEdge.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 1991-2008 OpenCFD Ltd.
   6      \\/     M anipulation  |
   7 -------------------------------------------------------------------------------
   8 License
   9     This file is part of OpenFOAM.
  10
  11     OpenFOAM is free software; you can redistribute it and/or modify it
  12     under the terms of the GNU General Public License as published by the
  13     Free Software Foundation; either version 2 of the License, or (at your
  14     option) any later version.
  15
  16     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
  17     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
  18     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
  19     for more details.
  20
  21     You should have received a copy of the GNU General Public License
  22     along with OpenFOAM; if not, write to the Free Software Foundation,
  23     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
  24
  25 Description
  26     arcEdge class : defines the arcEdge of a circle in terms of 3 points on its
  27     circumference
  28
  29 \*---------------------------------------------------------------------------*/
  30
  31 #include "arcEdge.H"
  32 #include "mathematicalConstants.H"
  33
  34 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  35
  36 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
  37
  38 namespace Foam
  39 {
  40     defineTypeNameAndDebug(arcEdge, 0);
  41
  42     // Add the curvedEdge constructor functions to the hash tables
  43     curvedEdge::addIstreamConstructorToTable<arcEdge>
  44         addArcEdgeIstreamConstructorToTable_;
  45 }
  46
  47
  48 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  49
  50 Foam::cylindricalCS Foam::arcEdge::calcAngle()
  51 {
  52     vector a = p2_ - p1_;
  53     vector b = p3_ - p1_;
  54
  55     // find centre of arcEdge
  56     scalar asqr = a & a;
  57     scalar bsqr = b & b;
  58     scalar adotb = a & b;
  59
  60     scalar denom = asqr*bsqr - adotb*adotb;
  61
  62     if (mag(denom) < VSMALL)
  63     {
  64         FatalErrorIn("cylindricalCS arcEdge::calcAngle()")
  65             << "Invalid arc definition - are the points co-linear?  Denom ="
  66             << denom
  67             << abort(FatalError);
  68     }
  69
  70     scalar fact = 0.5*(bsqr - adotb)/denom;
  71
  72     vector centre = 0.5*a + fact*((a ^ b) ^ a);
  73
  74     centre += p1_;
  75
  76     // find position vectors w.r.t. the arcEdge centre
  77     vector r1(p1_ - centre);
  78     vector r2(p2_ - centre);
  79     vector r3(p3_ - centre);
  80
  81     // find angles
  82     scalar tmp = (r3&r1)/(mag(r3)*mag(r1));
  83     angle_ = acos(tmp)*180.0/mathematicalConstant::pi;
  84
  85     // check if the vectors define an exterior or an interior arcEdge
  86     if (((r1  ^ r2)&(r1 ^ r3)) < 0.0) angle_ = 360 - angle_;
  87
  88     vector tempAxis(0.0,0.0,0.0);
  89
  90     if (angle_ <= 180.0)
  91     {
  92         tempAxis = r1 ^ r3;
  93
  94         if (mag(tempAxis)/(mag(r1)*mag(r3)) < 0.001) tempAxis = r1 ^ r2;
  95     }
  96     else
  97     {
  98         tempAxis = r3 ^ r1;
  99     }
 100
 101     radius_ = mag(r3);
 102
 103     // set up and return the local coordinate system
 104     return cylindricalCS("tmpCS", centre, tempAxis, r1);
 105 }
 106
 107
 108 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 109
 110 // Construct from components
 111 Foam::arcEdge::arcEdge
 112 (
 113     const pointField& points,
 114     const label start,
 115     const label end,
 116     const vector& P2
 117 )
 118 :
 119     curvedEdge(points, start, end),
 120     p1_(points_[start_]),
 121     p2_(P2),
 122     p3_(points_[end_]),
 123     cs_(calcAngle())
 124 {}
 125
 126
 127 // Construct from Istream
 128 Foam::arcEdge::arcEdge(const pointField& points, Istream& is)
 129 :
 130     curvedEdge(points, is),
 131     p1_(points_[start_]),
 132     p2_(is),
 133     p3_(points_[end_]),
 134     cs_(calcAngle())
 135 {}
 136
 137
 138 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 139
 140 Foam::vector Foam::arcEdge::position(const scalar lambda) const
 141 {
 142     if (lambda < 0 || lambda > 1)
 143     {
 144         FatalErrorIn("arcEdge::position(const scalar lambda) const")
 145             << "Parameter out of range, lambda = " << lambda
 146             << abort(FatalError);
 147     }
 148
 149     if (lambda < SMALL)
 150     {
 151         return p1_;
 152     }
 153     else if (lambda > 1-SMALL)
 154     {
 155         return p3_;
 156     }
 157     else
 158     {
 159         return cs_.globalPosition(vector(radius_, lambda*angle_, 0.0));
 160     }
 161 }
 162
 163
 164 //- Return the length of the curve
 165 Foam::scalar Foam::arcEdge::length() const
 166 {
 167     return angle_*radius_*mathematicalConstant::pi/180.0;
 168 }
 169
 170
 171 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 172
 173 // ************************************************************************* //