src/lagrangian/molecularDynamics/molecule/referredCell/referredCell.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 \*----------------------------------------------------------------------------*/
  26
  27 #include "referredCell.H"
  28 #include "moleculeCloud.H"
  29
  30 namespace Foam
  31 {
  32
  33 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  34
  35 void referredCell::setConstructionData
  36 (
  37     const polyMesh& mesh,
  38     const label sourceCell
  39 )
  40 {
  41     // * * * * * * * * * * * Points * * * * * * * * * * *
  42
  43     const labelList& points = mesh.cellPoints()[sourceCell];
  44
  45     vectorList sourceCellVertices(points.size());
  46
  47     forAll(sourceCellVertices, sCV)
  48     {
  49         sourceCellVertices[sCV] = mesh.points()[points[sCV]];
  50     }
  51
  52     vertexPositions_ = referPositions(sourceCellVertices);
  53
  54     // * * * * * * * * * * * Edges * * * * * * * * * * *
  55
  56     const labelList& edges = mesh.cellEdges()[sourceCell];
  57
  58     edgeList sourceCellEdges(edges.size());
  59
  60     forAll(sourceCellEdges, sCE)
  61     {
  62         sourceCellEdges[sCE] = mesh.edges()[edges[sCE]];
  63     }
  64
  65     locallyMapEdgeList(points, sourceCellEdges);
  66
  67     // * * * * * * * * * * * Faces * * * * * * * * * * *
  68
  69     labelList faces(mesh.cells()[sourceCell]);
  70
  71     vectorList sourceCellFaceCentres(faces.size());
  72
  73     vectorList sourceCellFaceAreas(faces.size());
  74
  75     labelListList sourceCellFaces(faces.size());
  76
  77     forAll(faces, f)
  78     {
  79         sourceCellFaces[f] = mesh.faces()[faces[f]];
  80
  81         sourceCellFaceCentres[f] = mesh.faceCentres()[faces[f]];
  82
  83         sourceCellFaceAreas[f] = mesh.faceAreas()[faces[f]];
  84     }
  85
  86     locallyMapFaceList(points, sourceCellFaces);
  87
  88     faceCentres_ = referPositions(sourceCellFaceCentres);
  89
  90     faceAreas_ = rotateVectors(sourceCellFaceAreas);
  91 }
  92
  93
  94 void referredCell::locallyMapEdgeList
  95 (
  96     const labelList& points,
  97     const edgeList& sourceCellEdges
  98 )
  99 {
 100     edges_.setSize(sourceCellEdges.size());
 101
 102     forAll(sourceCellEdges, sCE)
 103     {
 104         const edge& e(sourceCellEdges[sCE]);
 105
 106         edges_[sCE].start() = findIndex(points, e.start());
 107
 108         edges_[sCE].end() = findIndex(points, e.end());
 109
 110         if
 111         (
 112             edges_[sCE].start() == -1
 113          || edges_[sCE].end() == -1
 114         )
 115         {
 116             FatalErrorIn("Foam::referredCell::locallyMapEdgeList")
 117                 << "edgeList and points labelList for "
 118                 << "referred cell do not match: "
 119                 << nl << "points: " << points
 120                 << nl << "egdes: " << sourceCellEdges
 121                 << abort(FatalError);
 122         }
 123     }
 124 }
 125
 126
 127 void referredCell::locallyMapFaceList
 128 (
 129     const labelList& points,
 130     const labelListList& sourceCellFaces
 131 )
 132 {
 133     faces_.setSize(sourceCellFaces.size());
 134
 135     forAll(sourceCellFaces, sCF)
 136     {
 137         const labelList& sourceCellFace(sourceCellFaces[sCF]);
 138
 139         labelList& localFace(faces_[sCF]);
 140
 141         localFace.setSize(sourceCellFace.size());
 142
 143         forAll(sourceCellFace, p)
 144         {
 145             localFace[p] = findIndex(points, sourceCellFace[p]);
 146
 147             if (localFace[p] == -1)
 148             {
 149                 FatalErrorIn("Foam::referredCell::locallyMapEdgeList")
 150                     << "edgeList and points labelList for "
 151                     << "referred cell do not match: "
 152                     << nl << "points: " << points
 153                     << nl << "faces: " << sourceCellFaces
 154                     << abort(FatalError);
 155             }
 156         }
 157     }
 158 }
 159
 160
 161 vector referredCell::referPosition(const vector& positionToRefer)
 162 {
 163     return offset_ + (rotation_ & positionToRefer);
 164 }
 165
 166
 167 vectorList referredCell::referPositions(const vectorList& positionsToRefer)
 168 {
 169     return offset_ + (rotation_ & positionsToRefer);
 170 }
 171
 172
 173 vector referredCell::rotateVector(const vector& vectorToRotate)
 174 {
 175     return rotation_ & vectorToRotate;
 176 }
 177
 178
 179 vectorList referredCell::rotateVectors(const vectorList& vectorsToRotate)
 180 {
 181     return rotation_ & vectorsToRotate;
 182 }
 183
 184
 185 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 186
 187 referredCell::referredCell()
 188 :
 189     DynamicList<referredMolecule>(),
 190     sourceProc_(-1),
 191     sourceCell_(-1),
 192     vertexPositions_(),
 193     offset_(vector::zero),
 194     rotation_(I)
 195 {}
 196
 197
 198 referredCell::referredCell
 199 (
 200     const polyMesh& mesh,
 201     const label sourceProc,
 202     const label sourceCell,
 203     const vector& offset,
 204     const tensor& rotation
 205 )
 206 :
 207     DynamicList<referredMolecule>(),
 208     sourceProc_(sourceProc),
 209     sourceCell_(sourceCell),
 210     offset_(offset),
 211     rotation_(rotation)
 212 {
 213     setConstructionData(mesh, sourceCell);
 214 }
 215
 216
 217 referredCell::referredCell
 218 (
 219     const label sourceProc,
 220     const label sourceCell,
 221     const vectorList& vertexPositions,
 222     const edgeList& localEdges,
 223     const labelListList& localFaces,
 224     const vectorList& faceCentres,
 225     const vectorList& faceAreas,
 226     const vector& offset,
 227     const tensor& rotation
 228 )
 229 :
 230     DynamicList<referredMolecule>(),
 231     sourceProc_(sourceProc),
 232     sourceCell_(sourceCell),
 233     edges_(localEdges),
 234     faces_(localFaces),
 235     offset_(offset),
 236     rotation_(rotation)
 237 {
 238     // Supplied vertexPositions, faceCentres and faceAreas are of the
 239     // "original" cell, and need to be transformed to the referred
 240     // locations on construction
 241
 242     vertexPositions_ = referPositions(vertexPositions);
 243
 244     faceCentres_ = referPositions(faceCentres);
 245
 246     faceAreas_ = rotateVectors(faceAreas);
 247 }
 248
 249
 250 referredCell::referredCell
 251 (
 252     const polyMesh& mesh,
 253     const label sourceProc,
 254     const label sourceCell,
 255     const vector& cS,
 256     const vector& cD,
 257     const vector& nS,
 258     const vector& nD
 259 )
 260 :
 261     DynamicList<referredMolecule>(),
 262     sourceProc_(sourceProc),
 263     sourceCell_(sourceCell)
 264 {
 265     // It is assumed that the vectors originating from the faces being referred
 266     // here are correct periodic faces - i.e. they have the same area etc.
 267
 268     vector nA = -nS/mag(nS);
 269     vector nB = nD/mag(nD);
 270
 271     rotation_ = rotationTensor(nA, nB);
 272
 273     offset_ = cD - (rotation_ & cS);
 274
 275     // Allow sourceCell = -1 to create a dummy referredCell
 276     // to obtain the transformation
 277
 278     if(sourceCell >= 0)
 279     {
 280         setConstructionData(mesh, sourceCell);
 281     }
 282 }
 283
 284
 285 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 286
 287 referredCell::~referredCell()
 288 {}
 289
 290
 291 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 292
 293 referredCell referredCell::reRefer
 294 (
 295     const vector& cS,
 296     const vector& cD,
 297     const vector& nS,
 298     const vector& nD
 299 )
 300 {
 301     vector nA = -nS/mag(nS);
 302     vector nB = nD/mag(nD);
 303
 304     tensor newRotation = rotationTensor(nA, nB);
 305
 306     vector newOffset = cD - (newRotation & cS);
 307
 308     tensor reReferredRotation = newRotation & rotation_;
 309
 310     vector reReferredOffset = newOffset + (newRotation & offset_);
 311
 312     return referredCell
 313     (
 314         sourceProc_,
 315         sourceCell_,
 316         rotation_.T() & (vertexPositions_ - offset_),
 317         edges_,
 318         faces_,
 319         rotation_.T() & (faceCentres_ - offset_),
 320         rotation_.T() & (faceAreas_),
 321         reReferredOffset,
 322         reReferredRotation
 323     );
 324 }
 325
 326
 327 vector referredCell::referPosition(const vector& positionToRefer) const
 328 {
 329     return offset_ + (rotation_ & positionToRefer);
 330 }
 331
 332
 333 vectorList referredCell::referPositions
 334 (
 335     const vectorList& positionsToRefer
 336 ) const
 337 {
 338     return offset_ + (rotation_ & positionsToRefer);
 339 }
 340
 341
 342 vector referredCell::rotateVector(const vector& vectorToRotate) const
 343 {
 344     return rotation_ & vectorToRotate;
 345 }
 346
 347
 348 vectorList referredCell::rotateVectors(const vectorList& vectorsToRotate) const
 349 {
 350     return rotation_ & vectorsToRotate;
 351 }
 352
 353
 354 void referredCell::referInMols(const List<referredMolecule>& incomingMols)
 355 {
 356     clear();
 357
 358     forAll(incomingMols, iM)
 359     {
 360         append
 361         (
 362             referredMolecule
 363             (
 364                 incomingMols[iM].id(),
 365                 referPosition
 366                 (
 367                     incomingMols[iM].position()
 368                 )
 369             )
 370         );
 371     }
 372 }
 373
 374
 375 bool referredCell::duplicate(const referredCell& refCellDupl) const
 376 {
 377     return
 378     (
 379         sourceProc_ == refCellDupl.sourceProc()
 380         && sourceCell_ == refCellDupl.sourceCell()
 381         && mag(offset_ - refCellDupl.offset()) < moleculeCloud::transTol
 382     );
 383 }
 384
 385
 386 bool referredCell::duplicate(const label procNo,const label nCells) const
 387 {
 388     return
 389     (
 390         sourceProc_ == procNo
 391         && sourceCell_ < nCells
 392         && mag(offset_) < moleculeCloud::transTol
 393     );
 394 }
 395
 396
 397 // * * * * * * * * * * * * * * * Friend Functions  * * * * * * * * * * * * * //
 398
 399 Istream& operator>>(Istream& is, referredCell& rC)
 400 {
 401
 402     is >> rC.sourceProc_
 403         >> rC.sourceCell_
 404         >> rC.vertexPositions_
 405         >> rC.edges_
 406         >> rC.faces_
 407         >> rC.faceCentres_
 408         >> rC.faceAreas_
 409         >> rC.offset_
 410         >> rC.rotation_;
 411
 412     is.check("Istream& operator<<(Istream& f, const referredCell& rC");
 413
 414     return is;
 415 }
 416
 417
 418 Ostream& operator<<(Ostream& os, const referredCell& rC)
 419 {
 420
 421     os << rC.sourceProc()
 422         << token::SPACE << rC.sourceCell()
 423         << token::SPACE << rC.vertexPositions()
 424         << token::SPACE << rC.edges()
 425         << token::SPACE << rC.faces()
 426         << token::SPACE << rC.faceCentres()
 427         << token::SPACE << rC.faceAreas()
 428         << token::SPACE << rC.offset()
 429         << token::SPACE << rC.rotation();
 430
 431     os.check("Ostream& operator<<(Ostream& f, const referredCell& rC");
 432
 433     return os;
 434 }
 435
 436
 437 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 438
 439 }  // End namespace Foam
 440
 441 // ************************************************************************* //