src/lagrangian/intermediate/submodels/CloudFunctionObjects/FacePostProcessing/FacePostProcessing.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
   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
  13     the Free Software Foundation, either version 3 of the License, or
  14     (at your 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, see <http://www.gnu.org/licenses/>.
  23
  24 \*---------------------------------------------------------------------------*/
  25
  26 #include "FacePostProcessing.H"
  27 #include "Pstream.H"
  28 #include "ListListOps.H"
  29 #include "surfaceWriter.H"
  30 #include "globalIndex.H"
  31
  32 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  33
  34 template<class CloudType>
  35 void Foam::FacePostProcessing<CloudType>::makeLogFile
  36 (
  37     const word& zoneName,
  38     const label zoneI,
  39     const label nFaces,
  40     const scalar totArea
  41 )
  42 {
  43     // Create the output file if not already created
  44     if (log_)
  45     {
  46         if (debug)
  47         {
  48             Info<< "Creating output file." << endl;
  49         }
  50
  51         if (Pstream::master())
  52         {
  53             const fileName logDir = outputDir_/this->owner().time().timeName();
  54
  55             // Create directory if does not exist
  56             mkDir(logDir);
  57
  58             // Open new file at start up
  59             outputFilePtr_.set
  60             (
  61                 zoneI,
  62                 new OFstream(logDir/(type() + '_' + zoneName + ".dat"))
  63             );
  64
  65             outputFilePtr_[zoneI]
  66                 << "# Source    : " << type() << nl
  67                 << "# Face zone : " << zoneName << nl
  68                 << "# Faces     : " << nFaces << nl
  69                 << "# Area      : " << totArea << nl
  70                 << "# Time" << tab << "mass" << tab << "massFlux" << endl;
  71         }
  72     }
  73 }
  74
  75
  76 template<class CloudType>
  77 void Foam::FacePostProcessing<CloudType>::write()
  78 {
  79     const fvMesh& mesh = this->owner().mesh();
  80     const Time& time = mesh.time();
  81     const faceZoneMesh& fzm = mesh.faceZones();
  82     const scalar dt = time.deltaTValue();
  83
  84     totalTime_ += dt;
  85
  86     const scalar alpha = (totalTime_ - dt)/totalTime_;
  87     const scalar beta = dt/totalTime_;
  88
  89     forAll(faceZoneIDs_, zoneI)
  90     {
  91         massTotal_[zoneI] += mass_[zoneI];
  92         massFlux_[zoneI] = alpha*massFlux_[zoneI] + beta*mass_[zoneI]/dt;
  93     }
  94
  95     const label procI = Pstream::myProcNo();
  96
  97     Info<< "particleFaceFlux output:" << nl;
  98
  99     List<scalarField> zoneMassTotal(mass_.size());
 100     List<scalarField> zoneMassFlux(massFlux_.size());
 101     forAll(faceZoneIDs_, zoneI)
 102     {
 103         const word& zoneName = fzm[faceZoneIDs_[zoneI]].name();
 104
 105         scalarListList allProcMass(Pstream::nProcs());
 106         allProcMass[procI] = massTotal_[zoneI];
 107         Pstream::gatherList(allProcMass);
 108         zoneMassTotal[zoneI] =
 109             ListListOps::combine<scalarList>
 110             (
 111                 allProcMass, accessOp<scalarList>()
 112             );
 113         const scalar sumMassTotal = sum(zoneMassTotal[zoneI]);
 114
 115         scalarListList allProcMassFlux(Pstream::nProcs());
 116         allProcMassFlux[procI] = massFlux_[zoneI];
 117         Pstream::gatherList(allProcMassFlux);
 118         zoneMassFlux[zoneI] =
 119             ListListOps::combine<scalarList>
 120             (
 121                 allProcMassFlux, accessOp<scalarList>()
 122             );
 123         const scalar sumMassFlux = sum(zoneMassFlux[zoneI]);
 124
 125         Info<< "    " << zoneName
 126             << ": total mass = " << sumMassTotal
 127             << "; average mass flux = " << sumMassFlux
 128             << nl;
 129
 130         if (outputFilePtr_.set(zoneI))
 131         {
 132             OFstream& os = outputFilePtr_[zoneI];
 133             os  << time.timeName() << token::TAB << sumMassTotal << token::TAB
 134                 <<  sumMassFlux<< endl;
 135         }
 136     }
 137
 138     Info<< endl;
 139
 140
 141     if (surfaceFormat_ != "none")
 142     {
 143         forAll(faceZoneIDs_, zoneI)
 144         {
 145             const faceZone& fZone = fzm[faceZoneIDs_[zoneI]];
 146
 147             labelList pointToGlobal;
 148             labelList uniqueMeshPointLabels;
 149             autoPtr<globalIndex> globalPointsPtr =
 150                 mesh.globalData().mergePoints
 151                 (
 152                     fZone().meshPoints(),
 153                     fZone().meshPointMap(),
 154                     pointToGlobal,
 155                     uniqueMeshPointLabels
 156                 );
 157
 158             pointField uniquePoints(mesh.points(), uniqueMeshPointLabels);
 159             List<pointField> allProcPoints(Pstream::nProcs());
 160             allProcPoints[procI] = uniquePoints;
 161             Pstream::gatherList(allProcPoints);
 162
 163             faceList faces(fZone().localFaces());
 164             forAll(faces, i)
 165             {
 166                 inplaceRenumber(pointToGlobal, faces[i]);
 167             }
 168             List<faceList> allProcFaces(Pstream::nProcs());
 169             allProcFaces[procI] = faces;
 170             Pstream::gatherList(allProcFaces);
 171
 172             if (Pstream::master())
 173             {
 174                 pointField allPoints
 175                 (
 176                     ListListOps::combine<pointField>
 177                     (
 178                         allProcPoints, accessOp<pointField>()
 179                     )
 180                 );
 181
 182                 faceList allFaces
 183                 (
 184                     ListListOps::combine<faceList>
 185                     (
 186                         allProcFaces, accessOp<faceList>()
 187                     )
 188                 );
 189
 190                 autoPtr<surfaceWriter> writer
 191                 (
 192                     surfaceWriter::New(surfaceFormat_)
 193                 );
 194
 195                 writer->write
 196                 (
 197                     outputDir_/time.timeName(),
 198                     fZone.name(),
 199                     allPoints,
 200                     allFaces,
 201                     "massTotal",
 202                     zoneMassTotal[zoneI],
 203                     false
 204                 );
 205
 206                 writer->write
 207                 (
 208                     outputDir_/time.timeName(),
 209                     fZone.name(),
 210                     allPoints,
 211                     allFaces,
 212                     "massFlux",
 213                     zoneMassFlux[zoneI],
 214                     false
 215                 );
 216             }
 217         }
 218     }
 219
 220
 221     if (resetOnWrite_)
 222     {
 223         forAll(faceZoneIDs_, zoneI)
 224         {
 225             massFlux_[zoneI] = 0.0;
 226         }
 227         totalTime_ = 0.0;
 228     }
 229
 230     forAll(mass_, zoneI)
 231     {
 232         mass_[zoneI] = 0.0;
 233     }
 234
 235     // writeProperties();
 236 }
 237
 238
 239 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 240
 241 template<class CloudType>
 242 Foam::FacePostProcessing<CloudType>::FacePostProcessing
 243 (
 244     const dictionary& dict,
 245     CloudType& owner
 246 )
 247 :
 248     CloudFunctionObject<CloudType>(dict, owner, typeName),
 249     faceZoneIDs_(),
 250     surfaceFormat_(this->coeffDict().lookup("surfaceFormat")),
 251     resetOnWrite_(this->coeffDict().lookup("resetOnWrite")),
 252     totalTime_(0.0),
 253     mass_(),
 254     massTotal_(),
 255     massFlux_(),
 256     log_(this->coeffDict().lookup("log")),
 257     outputFilePtr_(),
 258     outputDir_(owner.mesh().time().path())
 259 {
 260     wordList faceZoneNames(this->coeffDict().lookup("faceZones"));
 261     mass_.setSize(faceZoneNames.size());
 262     massTotal_.setSize(faceZoneNames.size());
 263     massFlux_.setSize(faceZoneNames.size());
 264
 265     outputFilePtr_.setSize(faceZoneNames.size());
 266
 267     if (Pstream::parRun())
 268     {
 269         // Put in undecomposed case (Note: gives problems for
 270         // distributed data running)
 271         outputDir_ =
 272             outputDir_/".."/"postProcessing"/cloud::prefix/owner.name();
 273     }
 274     else
 275     {
 276         outputDir_ =
 277             outputDir_/"postProcessing"/cloud::prefix/owner.name();
 278     }
 279
 280     DynamicList<label> zoneIDs;
 281     const faceZoneMesh& fzm = owner.mesh().faceZones();
 282     const surfaceScalarField& magSf = owner.mesh().magSf();
 283     const polyBoundaryMesh& pbm = owner.mesh().boundaryMesh();
 284     forAll(faceZoneNames, i)
 285     {
 286         const word& zoneName = faceZoneNames[i];
 287         label zoneI = fzm.findZoneID(zoneName);
 288         if (zoneI != -1)
 289         {
 290             zoneIDs.append(zoneI);
 291             const faceZone& fz = fzm[zoneI];
 292             mass_[i].setSize(fz.size(), 0.0);
 293             massTotal_[i].setSize(fz.size(), 0.0);
 294             massFlux_[i].setSize(fz.size(), 0.0);
 295
 296             label nFaces = returnReduce(fz.size(), sumOp<label>());
 297             Info<< "        " << zoneName << " faces: " << nFaces << nl;
 298
 299             scalar totArea = 0.0;
 300             forAll(fz, j)
 301             {
 302                 label faceI = fz[j];
 303                 if (faceI < owner.mesh().nInternalFaces())
 304                 {
 305                     totArea += magSf[fz[j]];
 306                 }
 307                 else
 308                 {
 309                     label bFaceI = faceI - owner.mesh().nInternalFaces();
 310                     label patchI = pbm.patchID()[bFaceI];
 311                     const polyPatch& pp = pbm[patchI];
 312
 313                     if
 314                     (
 315                         !pp.coupled()
 316                      || refCast<const coupledPolyPatch>(pp).owner()
 317                     )
 318                     {
 319                         label localFaceI = pp.whichFace(faceI);
 320                         totArea += magSf.boundaryField()[patchI][localFaceI];
 321                     }
 322                 }
 323             }
 324             totArea = returnReduce(totArea, sumOp<scalar>());
 325
 326             makeLogFile(zoneName, i, nFaces, totArea);
 327         }
 328     }
 329
 330     faceZoneIDs_.transfer(zoneIDs);
 331
 332     // readProperties(); AND initialise mass... fields
 333 }
 334
 335
 336 template<class CloudType>
 337 Foam::FacePostProcessing<CloudType>::FacePostProcessing
 338 (
 339     const FacePostProcessing<CloudType>& pff
 340 )
 341 :
 342     CloudFunctionObject<CloudType>(pff),
 343     faceZoneIDs_(pff.faceZoneIDs_),
 344     surfaceFormat_(pff.surfaceFormat_),
 345     resetOnWrite_(pff.resetOnWrite_),
 346     totalTime_(pff.totalTime_),
 347     mass_(pff.mass_),
 348     massTotal_(pff.massTotal_),
 349     massFlux_(pff.massFlux_),
 350     log_(pff.log_),
 351     outputFilePtr_(),
 352     outputDir_(pff.outputDir_)
 353 {}
 354
 355
 356 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 357
 358 template<class CloudType>
 359 Foam::FacePostProcessing<CloudType>::~FacePostProcessing()
 360 {}
 361
 362
 363 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 364
 365 template<class CloudType>
 366 void Foam::FacePostProcessing<CloudType>::postFace
 367 (
 368     const parcelType& p,
 369     const label faceI
 370 )
 371 {
 372     if
 373     (
 374         this->owner().solution().output()
 375      || this->owner().solution().transient()
 376     )
 377     {
 378         const faceZoneMesh& fzm = this->owner().mesh().faceZones();
 379
 380         forAll(faceZoneIDs_, i)
 381         {
 382             const faceZone& fz = fzm[faceZoneIDs_[i]];
 383
 384             label faceId = -1;
 385             forAll(fz, j)
 386             {
 387                 if (fz[j] == faceI)
 388                 {
 389                     faceId = j;
 390                     break;
 391                 }
 392             }
 393
 394             if (faceId != -1)
 395             {
 396                 mass_[i][faceId] += p.mass()*p.nParticle();
 397             }
 398         }
 399     }
 400 }
 401
 402
 403 // ************************************************************************* //