src/postProcessing/functionObjects/check/mixingPlaneCheck/mixingPlaneCheckFunctionObject.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright held by original author
   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 Author
  26     Martin Beaudoin, Hydro-Quebec, 2009.  All rights reserved
  27
  28 \*---------------------------------------------------------------------------*/
  29
  30 #include "mixingPlaneCheckFunctionObject.H"
  31 #include "addToRunTimeSelectionTable.H"
  32 #include "mixingPlaneFvsPatchFields.H"
  33 #include "surfaceFields.H"
  34 #include "mixingPlanePolyPatch.H"
  35
  36 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
  37
  38 namespace Foam
  39 {
  40     defineTypeNameAndDebug(mixingPlaneCheckFunctionObject, 0);
  41
  42     addToRunTimeSelectionTable
  43     (
  44         functionObject,
  45         mixingPlaneCheckFunctionObject,
  46         dictionary
  47     );
  48 }
  49
  50
  51 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
  52
  53 Foam::mixingPlaneCheckFunctionObject::mixingPlaneCheckFunctionObject
  54 (
  55     const word& name,
  56     const Time& t,
  57     const dictionary& dict
  58 )
  59 :
  60     functionObject(name),
  61     time_(t),
  62     regionName_(polyMesh::defaultRegion),
  63     phiName_(dict.lookup("phi"))
  64 {
  65     if (dict.found("region"))
  66     {
  67         dict.lookup("region") >> regionName_;
  68     }
  69
  70     Info << "Creating mixingPlane check functionObject" << endl;
  71 }
  72
  73
  74 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
  75
  76 bool Foam::mixingPlaneCheckFunctionObject::start()
  77 {
  78     return true;
  79 }
  80
  81
  82 bool Foam::mixingPlaneCheckFunctionObject::execute()
  83 {
  84     const polyMesh& mesh =
  85         time_.lookupObject<polyMesh>(regionName_);
  86
  87     boolList visited(mesh.boundaryMesh().size(), false);
  88
  89     forAll (mesh.boundaryMesh(), patchI)
  90     {
  91         if (
  92             isA<mixingPlanePolyPatch>(mesh.boundaryMesh()[patchI]) &&
  93             mesh.boundaryMesh()[patchI].size()
  94         )
  95         {
  96             if (!visited[patchI])
  97             {
  98                 visited[patchI] = true;
  99
 100                 const mixingPlanePolyPatch& mixingMaster =
 101                     refCast<const mixingPlanePolyPatch>
 102                     (
 103                         mesh.boundaryMesh()[patchI]
 104                     );
 105
 106                 const label shadowPatchI = mixingMaster.shadowIndex();
 107
 108                 visited[shadowPatchI] = true;
 109
 110                 const mixingPlanePolyPatch& mixingShadow =
 111                     mixingMaster.shadow();
 112
 113                 // Get access to the mixing plane patch
 114                 const standAlonePatch& mixingPlanePatch =
 115                     mixingMaster.patchToPatch().mixingPlanePatch();
 116
 117                 // Calculate areas of the mixing patch
 118                 scalarField mixingPlanePatchAreas(mixingPlanePatch.size());
 119                 const vectorField& mixingPlanePatchPoints =
 120                     mixingPlanePatch.points();
 121
 122                 forAll (mixingPlanePatch, faceI)
 123                 {
 124                     mixingPlanePatchAreas[faceI] =
 125                         mixingPlanePatch[faceI].mag(mixingPlanePatchPoints);
 126                 }
 127
 128                 const scalarField masterAreas = mag(mixingMaster.faceAreas());
 129                 const scalarField shadowAreas = mag(mixingShadow.faceAreas());
 130
 131                 // Until the mixingPlane is fully parallelized, we stick with
 132                 // the serial version of sum. The interface is residing on a
 133                 // single processor when running in parallel
 134                 scalar sumMasterAreas = sum(masterAreas);
 135                 scalar sumShadowAreas = sum(shadowAreas);
 136 //                 scalar sumMixingAreas = sum(mixingPlanePatchAreas);
 137
 138 #if 0   // Remove this for now
 139                 Info<< "Mixing plane functionObject: area check " << nl
 140                     << "     "  << "Master " << mixingMaster.name()
 141                     << " = " << sumMasterAreas << nl
 142                     << "     "  << "Shadow  " << mixingMaster.shadow().name()
 143                     << " = " << sumShadowAreas << nl
 144                     << "     "  << "Mixing = " << sumMixingAreas << nl
 145                     << "     "  << "mixingPlanePatchAreas = " << mixingPlanePatchAreas
 146                     //<< "     "  << "mixingPlanePatchPoints = " << mixingPlanePatchPoints
 147                     << endl;
 148 #endif
 149
 150                 // Calculate master to strip sum
 151                 scalarField masterToStripsAreas(mixingPlanePatch.size(), 0);
 152
 153                 const labelListList& mppAddr =
 154                     mixingMaster.patchToPatch().masterProfileToPatchAddr();
 155
 156                 const scalarListList& mppWeights =
 157                     mixingMaster.patchToPatch().masterProfileToPatchWeights();
 158
 159                 forAll (masterAreas, masterI)
 160                 {
 161                     const labelList& curMppAddr = mppAddr[masterI];
 162                     const scalarList& curMppWeights = mppWeights[masterI];
 163
 164                     forAll (curMppAddr, i)
 165                     {
 166                         masterToStripsAreas[curMppAddr[i]] +=
 167                             curMppWeights[i]*masterAreas[masterI];
 168                     }
 169                 }
 170
 171 #if 0
 172                 // Does not work when in cylindrical coordinates
 173                 Info<< "Master scaling = "
 174                     << mixingPlanePatchAreas/masterToStripsAreas << endl;
 175 #endif
 176
 177                 // Calculate shadow to strip sum
 178                 scalarField shadowToStripsAreas(mixingPlanePatch.size(), 0);
 179
 180                 const labelListList& sppAddr =
 181                     mixingMaster.patchToPatch().slaveProfileToPatchAddr();
 182
 183                 const scalarListList& sppWeights =
 184                     mixingMaster.patchToPatch().slaveProfileToPatchWeights();
 185
 186                 forAll (shadowAreas, shadowI)
 187                 {
 188                     const labelList& curSppAddr = sppAddr[shadowI];
 189                     const scalarList& curSppWeights = sppWeights[shadowI];
 190
 191                     forAll (curSppAddr, i)
 192                     {
 193                         shadowToStripsAreas[curSppAddr[i]] +=
 194                             curSppWeights[i]*shadowAreas[shadowI];
 195                     }
 196                 }
 197
 198 #if 0
 199                 // Does not work when in cylindrical coordinates
 200                 Info<< "Shadow scaling = "
 201                     << mixingPlanePatchAreas/shadowToStripsAreas << endl;
 202 #endif
 203
 204                // Old way of computing phi balance
 205
 206                 if( mesh.foundObject<surfaceScalarField>(phiName_) )
 207                 {
 208                     const surfaceScalarField& phi =
 209                         mesh.lookupObject<surfaceScalarField>(phiName_);
 210
 211                     // Calculate local and shadow flux
 212                     scalar masterPatchScaleFactor_ = 1.0;
 213                     scalar shadowPatchScaleFactor_ = sumMasterAreas/sumShadowAreas;
 214                     scalar localFlux    = masterPatchScaleFactor_ * sum(phi.boundaryField()[patchI]);
 215                     scalar localFluxMag = mag(localFlux);
 216
 217                     scalar shadowFlux    = shadowPatchScaleFactor_ * sum(phi.boundaryField()[shadowPatchI]);
 218 //                     scalar shadowFluxMag = mag(shadowFlux);
 219
 220                     Info<< "Mixing plane pair "
 221                         << "(" << mixingMaster.name() << ", "
 222                         << mixingMaster.shadow().name() << ") : "
 223                         << localFlux << " " << shadowFlux
 224                         << " Diff = " << localFlux + shadowFlux << " or "
 225                         << mag(localFlux + shadowFlux)/(localFluxMag + SMALL)*100
 226                         << " %" << endl;
 227                 }
 228             }
 229         }
 230     }
 231
 232
 233 //         {
 234 //             word patchName = phi.boundaryField()[patchI].patch().name();
 235
 236 //             if (patchName == masterPatchName_ && !visited[patchI])
 237 //             {
 238 //                 visited[patchI] = true;
 239
 240 //                 // Calculate local and shadow flux
 241 //                 scalar localFlux    = masterPatchScaleFactor_ * sum(phi.boundaryField()[patchI]);
 242 //                 //scalar localFluxMag = masterPatchScaleFactor_ * sumMag(phi.boundaryField()[patchI]);
 243 //                 scalar localFluxMag = mag(localFlux);
 244
 245 //                 const mixingPlanePolyPatch& mixingPlanePatch =
 246 //                     refCast<const mixingPlanePolyPatch>
 247 //                     (
 248 //                         phi.boundaryField()[patchI].patch().patch()
 249 //                     );
 250
 251 //                 const label shadowPatchI = mixingPlanePatch.shadowIndex();
 252
 253 //                 visited[shadowPatchI] = true;
 254
 255 //                 scalar shadowFlux    = shadowPatchScaleFactor_ * sum(phi.boundaryField()[shadowPatchI]);
 256 //                 //scalar shadowFluxMag = shadowPatchScaleFactor_ * sumMag(phi.boundaryField()[shadowPatchI]);
 257 //                 scalar shadowFluxMag = mag(shadowFlux);
 258
 259 //                 Info<< "mixingPlane pair " << name_ << " (" << mixingPlanePatch.name() << ", " << mixingPlanePatch.shadow().name() << ") : "
 260 //                     << " flux: " << localFlux << " " << shadowFlux
 261 //                     << " : mag: " <<  localFluxMag << " " << shadowFluxMag
 262 //                     << " Diff = " << localFlux + shadowFlux << " or "
 263 //                     << mag(localFlux + shadowFlux)/(localFluxMag + SMALL)*100
 264 //                     << " %" << endl;
 265 //             }
 266 //         }
 267 //     }
 268
 269     return true;
 270 }
 271
 272
 273 bool Foam::mixingPlaneCheckFunctionObject::read(const dictionary& dict)
 274 {
 275     return false;
 276 }
 277
 278
 279 // ************************************************************************* //