VoFPatchTransfer.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2017-2019 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 "VoFPatchTransfer.H"
27 #include "twoPhaseMixtureThermo.H"
28 #include "thermoSingleLayer.H"
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace regionModels
36 {
37 namespace surfaceFilmModels
38 {
39 
40 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
41 
42 defineTypeNameAndDebug(VoFPatchTransfer, 0);
43 addToRunTimeSelectionTable(transferModel, VoFPatchTransfer, dictionary);
44 
45 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
46 
48 (
49  surfaceFilmRegionModel& film,
50  const dictionary& dict
51 )
52 :
53  transferModel(type(), film, dict),
54  deltaFactorToVoF_
55  (
56  coeffDict_.lookupOrDefault<scalar>("deltaFactorToVoF", 1.0)
57  ),
58  deltaFactorToFilm_
59  (
60  coeffDict_.lookupOrDefault<scalar>("deltaFactorToFilm", 0.5)
61  ),
62  alphaToVoF_
63  (
64  coeffDict_.lookupOrDefault<scalar>("alphaToVoF", 0.5)
65  ),
66  alphaToFilm_
67  (
68  coeffDict_.lookupOrDefault<scalar>("alphaToFilm", 0.1)
69  ),
70  transferRateCoeff_
71  (
72  coeffDict_.lookupOrDefault<scalar>("transferRateCoeff", 0.1)
73  )
74 {
75  const polyBoundaryMesh& pbm = film.regionMesh().boundaryMesh();
76  patchIDs_.setSize
77  (
78  pbm.size() - film.regionMesh().globalData().processorPatches().size()
79  );
80 
81  if (coeffDict_.found("patches"))
82  {
83  const wordReList patchNames(coeffDict_.lookup("patches"));
84  const labelHashSet patchSet = pbm.patchSet(patchNames);
85 
86  Info<< " applying to patches:" << nl;
87 
88  label pidi = 0;
89  forAllConstIter(labelHashSet, patchSet, iter)
90  {
91  const label patchi = iter.key();
92  patchIDs_[pidi++] = patchi;
93  Info<< " " << pbm[patchi].name() << endl;
94  }
95  patchIDs_.setSize(pidi);
96  patchTransferredMasses_.setSize(pidi, 0);
97  }
98  else
99  {
100  Info<< " applying to all patches" << endl;
101 
102  forAll(patchIDs_, patchi)
103  {
104  patchIDs_[patchi] = patchi;
105  }
106 
107  patchTransferredMasses_.setSize(patchIDs_.size(), 0);
108  }
109 
110  if (!patchIDs_.size())
111  {
113  << "No patches selected"
114  << exit(FatalError);
115  }
116 }
117 
118 
119 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
120 
121 VoFPatchTransfer::~VoFPatchTransfer()
122 {}
123 
124 
125 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
126 
128 (
129  scalarField& availableMass,
130  scalarField& massToTransfer
131 )
132 {
134 }
135 
136 
138 (
139  scalarField& availableMass,
140  scalarField& massToTransfer,
141  scalarField& energyToTransfer
142 )
143 {
144  // Do not correct if no patches selected
145  if (!patchIDs_.size()) return;
146 
147  const thermoSingleLayer& film = filmType<thermoSingleLayer>();
148 
149  const scalarField& delta = film.delta();
150  const scalarField& rho = film.rho();
151  const scalarField& magSf = film.magSf();
152 
153  const polyBoundaryMesh& pbm = film.regionMesh().boundaryMesh();
154 
155 
156  const twoPhaseMixtureThermo& thermo
157  (
158  film.primaryMesh().lookupObject<twoPhaseMixtureThermo>
159  (
161  )
162  );
163 
164  const volScalarField& heVoF = thermo.thermo1().he();
165  const volScalarField& TVoF = thermo.thermo1().T();
166  const volScalarField CpVoF(thermo.thermo1().Cp());
167  const volScalarField& rhoVoF = thermo.thermo1().rho()();
168  const volScalarField& alphaVoF = thermo.alpha1();
169 
170  forAll(patchIDs_, pidi)
171  {
172  const label patchi = patchIDs_[pidi];
173  label primaryPatchi = -1;
174 
175  forAll(film.intCoupledPatchIDs(), i)
176  {
177  const label filmPatchi = film.intCoupledPatchIDs()[i];
178 
179  if (filmPatchi == patchi)
180  {
181  primaryPatchi = film.primaryPatchIDs()[i];
182  }
183  }
184 
185  if (primaryPatchi != -1)
186  {
187  scalarField deltaCoeffs
188  (
189  film.primaryMesh().boundary()[primaryPatchi].deltaCoeffs()
190  );
191  film.toRegion(patchi, deltaCoeffs);
192 
193  scalarField hp(heVoF.boundaryField()[primaryPatchi]);
194  film.toRegion(patchi, hp);
195 
196  scalarField Tp(TVoF.boundaryField()[primaryPatchi]);
197  film.toRegion(patchi, Tp);
198 
199  scalarField Cpp(CpVoF.boundaryField()[primaryPatchi]);
200  film.toRegion(patchi, Cpp);
201 
202  scalarField rhop(rhoVoF.boundaryField()[primaryPatchi]);
203  film.toRegion(patchi, rhop);
204 
205  scalarField alphap(alphaVoF.boundaryField()[primaryPatchi]);
206  film.toRegion(patchi, alphap);
207 
208  scalarField Vp
209  (
210  film.primaryMesh().boundary()[primaryPatchi]
211  .patchInternalField(film.primaryMesh().V())
212  );
213  film.toRegion(patchi, Vp);
214 
215  const polyPatch& pp = pbm[patchi];
216  const labelList& faceCells = pp.faceCells();
217 
218  // Accumulate the total mass removed from patch
219  scalar dMassPatch = 0;
220 
221  forAll(faceCells, facei)
222  {
223  const label celli = faceCells[facei];
224 
225  scalar dMass = 0;
226 
227  if
228  (
229  delta[celli] > 2*deltaFactorToVoF_/deltaCoeffs[facei]
230  || alphap[facei] > alphaToVoF_
231  )
232  {
233  dMass =
234  transferRateCoeff_*delta[celli]*rho[celli]*magSf[celli];
235 
236  massToTransfer[celli] += dMass;
237  energyToTransfer[celli] += dMass*film.h()[celli];
238  }
239 
240  if
241  (
242  alphap[facei] > 0
243  && delta[celli] < 2*deltaFactorToFilm_/deltaCoeffs[facei]
244  && alphap[facei] < alphaToFilm_
245  )
246  {
247  dMass =
248  -transferRateCoeff_*alphap[facei]*rhop[facei]*Vp[facei];
249 
250  massToTransfer[celli] += dMass;
251  energyToTransfer[celli] += dMass*hp[facei];
252  }
253 
254  availableMass[celli] -= dMass;
255  dMassPatch += dMass;
256  }
257 
258  patchTransferredMasses_[pidi] += dMassPatch;
259  addToTransferredMass(dMassPatch);
260  }
261  }
262 
264 
265  if (writeTime())
266  {
267  scalarField patchTransferredMasses0
268  (
269  getModelProperty<scalarField>
270  (
271  "patchTransferredMasses",
272  scalarField(patchTransferredMasses_.size(), 0)
273  )
274  );
275 
276  scalarField patchTransferredMassTotals(patchTransferredMasses_);
277  Pstream::listCombineGather
278  (
279  patchTransferredMassTotals,
280  plusEqOp<scalar>()
281  );
282  patchTransferredMasses0 += patchTransferredMassTotals;
283 
284  setModelProperty<scalarField>
285  (
286  "patchTransferredMasses",
287  patchTransferredMasses0
288  );
289 
290  patchTransferredMasses_ = 0;
291  }
292 }
293 
294 
295 void VoFPatchTransfer::patchTransferredMassTotals
296 (
297  scalarField& patchMasses
298 ) const
299 {
300  // Do not correct if no patches selected
301  if (!patchIDs_.size()) return;
302 
303  scalarField patchTransferredMasses
304  (
305  getModelProperty<scalarField>
306  (
307  "patchTransferredMasses",
308  scalarField(patchTransferredMasses_.size(), 0)
309  )
310  );
311 
312  scalarField patchTransferredMassTotals(patchTransferredMasses_);
313  Pstream::listCombineGather(patchTransferredMassTotals, plusEqOp<scalar>());
314 
315  forAll(patchIDs_, pidi)
316  {
317  const label patchi = patchIDs_[pidi];
318  patchMasses[patchi] +=
319  patchTransferredMasses[pidi] + patchTransferredMassTotals[pidi];
320  }
321 }
322 
323 
324 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
325 
326 } // End namespace surfaceFilmModels
327 } // End namespace regionModels
328 } // End namespace Foam
329 
330 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
rhoReactionThermo & thermo
Definition: createFields.H:28
addToRunTimeSelectionTable(surfaceFilmRegionModel, kinematicSingleLayer, mesh)
VoFPatchTransfer(surfaceFilmRegionModel &film, const dictionary &dict)
Construct from surface film model.
Macros for easy insertion into run-time selection tables.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
HashSet< label, Hash< label > > labelHashSet
A HashSet with label keys.
Definition: HashSet.H:211
wordList patchNames(nPatches)
List< label > labelList
A List of labels.
Definition: labelList.H:56
const word dictName("particleTrackDict")
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
volScalarField scalarField(fieldObject, mesh)
static const char nl
Definition: Ostream.H:260
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(rho0+psi *p, rhoMin);# 1 "/home/ubuntu/OpenFOAM-8/applications/solvers/multiphase/cavitatingFoam/alphavPsi.H" 1{ alphav=max(min((rho - rholSat)/(rhovSat - rholSat), scalar(1)), scalar(0));alphal=1.0 - alphav;Info<< "max-min alphav: "<< max(alphav).value()<< " "<< min(alphav).value()<< endl;psiModel-> correct()
Definition: pEqn.H:68
label patchi
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
messageStream Info
List< wordRe > wordReList
A List of wordRe (word or regular expression)
Definition: wordReList.H:50
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:366
defineTypeNameAndDebug(kinematicSingleLayer, 0)
Namespace for OpenFOAM.