kinematicSingleLayer.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) 2011-2018 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 "kinematicSingleLayer.H"
27 #include "fvm.H"
28 #include "fvcDiv.H"
29 #include "fvcLaplacian.H"
30 #include "fvcSnGrad.H"
31 #include "fvcReconstruct.H"
32 #include "fvcVolumeIntegrate.H"
33 #include "fvcFlux.H"
35 #include "mappedWallPolyPatch.H"
36 #include "mapDistribute.H"
37 #include "filmThermoModel.H"
38 
39 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
40 
41 namespace Foam
42 {
43 namespace regionModels
44 {
45 namespace surfaceFilmModels
46 {
47 
48 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
49 
50 defineTypeNameAndDebug(kinematicSingleLayer, 0);
51 
52 addToRunTimeSelectionTable(surfaceFilmRegionModel, kinematicSingleLayer, mesh);
53 
54 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
55 
57 {
59  {
60  const dictionary& solution = this->solution().subDict("PISO");
61  solution.lookup("momentumPredictor") >> momentumPredictor_;
62  solution.readIfPresent("nOuterCorr", nOuterCorr_);
63  solution.lookup("nCorr") >> nCorr_;
64  solution.lookup("nNonOrthCorr") >> nNonOrthCorr_;
65 
66  return true;
67  }
68  else
69  {
70  return false;
71  }
72 }
73 
74 
76 {
77  rho_ == filmThermo_->rho();
78  mu_ == filmThermo_->mu();
79  sigma_ == filmThermo_->sigma();
80 }
81 
82 
84 {
85  if (debug)
86  {
88  }
89 
92  pSpPrimary_ == dimensionedScalar("zero", pSp_.dimensions(), 0.0);
93 }
94 
95 
97 {
98  if (debug)
99  {
100  InfoInFunction << endl;
101  }
102 
103  // Update fields from primary region via direct mapped
104  // (coupled) boundary conditions
109 }
110 
111 
113 {
114  if (debug)
115  {
116  InfoInFunction << endl;
117  }
118 
119  volScalarField::Boundary& rhoSpPrimaryBf =
121 
122  volVectorField::Boundary& USpPrimaryBf =
124 
125  volScalarField::Boundary& pSpPrimaryBf =
127 
128  // Convert accummulated source terms into per unit area per unit time
129  const scalar deltaT = time_.deltaTValue();
131  {
132  scalarField rpriMagSfdeltaT
133  (
134  (1.0/deltaT)
135  /primaryMesh().magSf().boundaryField()[patchi]
136  );
137 
138  rhoSpPrimaryBf[patchi] *= rpriMagSfdeltaT;
139  USpPrimaryBf[patchi] *= rpriMagSfdeltaT;
140  pSpPrimaryBf[patchi] *= rpriMagSfdeltaT;
141  }
142 
143  // Retrieve the source fields from the primary region via direct mapped
144  // (coupled) boundary conditions
145  // - fields require transfer of values for both patch AND to push the
146  // values into the first layer of internal cells
150 
151  // update addedMassTotal counter
152  if (time().writeTime())
153  {
154  scalar addedMassTotal = 0.0;
155  outputProperties().readIfPresent("addedMassTotal", addedMassTotal);
156  addedMassTotal += returnReduce(addedMassTotal_, sumOp<scalar>());
157  outputProperties().add("addedMassTotal", addedMassTotal, true);
158  addedMassTotal_ = 0.0;
159  }
160 }
161 
162 
164 {
165  return tmp<volScalarField>
166  (
167  new volScalarField
168  (
169  IOobject
170  (
171  typeName + ":pu",
172  time_.timeName(),
173  regionMesh(),
176  ),
177  pPrimary_ // pressure (mapped from primary region)
178  - pSp_ // accumulated particle impingement
179  - fvc::laplacian(sigma_, delta_) // surface tension
180  )
181  );
182 }
183 
184 
186 {
187  return tmp<volScalarField>
188  (
189  new volScalarField
190  (
191  IOobject
192  (
193  typeName + ":pp",
194  time_.timeName(),
195  regionMesh(),
198  ),
199  -rho_*gNormClipped() // hydrostatic effect only
200  )
201  );
202 }
203 
204 
206 {
208 }
209 
210 
212 {
213  if (debug)
214  {
215  InfoInFunction << endl;
216  }
217 
218  // Update injection model - mass returned is mass available for injection
220 
221  // Update transfer model - mass returned is mass available for transfer
223 
224  // Update mass source field
226 
227  turbulence_->correct();
228 }
229 
230 
232 {
233  const volScalarField deltaRho0(deltaRho_);
234 
235  solveContinuity();
236 
237  if (debug)
238  {
242  + dimensionedScalar("small", dimMass*dimVolume, rootVSmall);
243 
244  const scalar sumLocalContErr =
245  (
246  fvc::domainIntegrate(mag(mass - magSf()*deltaRho0))/totalMass
247  ).value();
248 
249  const scalar globalContErr =
250  (
251  fvc::domainIntegrate(mass - magSf()*deltaRho0)/totalMass
252  ).value();
253 
255 
257  << "Surface film: " << type() << nl
258  << " time step continuity errors: sum local = "
259  << sumLocalContErr << ", global = " << globalContErr
260  << ", cumulative = " << cumulativeContErr_ << endl;
261  }
262 }
263 
264 
266 {
267  if (debug)
268  {
269  InfoInFunction << endl;
270  }
271 
272  solve
273  (
275  + fvc::div(phi_)
276  ==
277  - rhoSp_
278  );
279 }
280 
281 
283 {
284  // Push boundary film velocity values into internal field
285  for (label i=0; i<intCoupledPatchIDs_.size(); i++)
286  {
288  const polyPatch& pp = regionMesh().boundaryMesh()[patchi];
291  }
292  Uw_ -= nHat()*(Uw_ & nHat());
294 
295  Us_ = turbulence_->Us();
296 }
297 
298 
300 (
301  const volScalarField& pu,
302  const volScalarField& pp
303 )
304 {
305  if (debug)
306  {
307  InfoInFunction << endl;
308  }
309 
310  // Momentum
312  (
314  + fvm::div(phi_, U_)
315  ==
316  - USp_
317  // - fvm::SuSp(rhoSp_, U_)
318  - rhoSp_*U_
319  + forces_.correct(U_)
320  + turbulence_->Su(U_)
321  );
322 
323  fvVectorMatrix& UEqn = tUEqn.ref();
324 
325  UEqn.relax();
326 
327  if (momentumPredictor_)
328  {
329  solve
330  (
331  UEqn
332  ==
334  (
336  * (
337  regionMesh().magSf()
338  * (
339  fvc::snGrad(pu, "snGrad(p)")
340  + fvc::snGrad(pp, "snGrad(p)")*fvc::interpolate(delta_)
342  )
343  - fvc::flux(rho_*gTan())
344  )
345  )
346  );
347 
348  // Remove any patch-normal components of velocity
349  U_ -= nHat()*(nHat() & U_);
351  }
352 
353  return tUEqn;
354 }
355 
356 
358 (
359  const volScalarField& pu,
360  const volScalarField& pp,
361  const fvVectorMatrix& UEqn
362 )
363 {
364  if (debug)
365  {
366  InfoInFunction << endl;
367  }
368 
369  volScalarField rUA(1.0/UEqn.A());
370  U_ = rUA*UEqn.H();
371 
374 
375  surfaceScalarField phiAdd
376  (
377  "phiAdd",
378  regionMesh().magSf()
379  * (
380  fvc::snGrad(pu, "snGrad(p)")
381  + fvc::snGrad(pp, "snGrad(p)")*fvc::interpolate(delta_)
382  )
383  - fvc::flux(rho_*gTan())
384  );
385  constrainFilmField(phiAdd, 0.0);
386 
388  (
389  "phid",
390  fvc::flux(U_*rho_) - deltarUAf*phiAdd*rhof
391  );
392  constrainFilmField(phid, 0.0);
393 
394  surfaceScalarField ddrhorUAppf
395  (
396  "deltaCoeff",
397  fvc::interpolate(delta_)*deltarUAf*rhof*fvc::interpolate(pp)
398  );
399 
401 
402  for (int nonOrth=0; nonOrth<=nNonOrthCorr_; nonOrth++)
403  {
404  // Film thickness equation
405  fvScalarMatrix deltaEqn
406  (
407  fvm::ddt(rho_, delta_)
408  + fvm::div(phid, delta_)
409  - fvm::laplacian(ddrhorUAppf, delta_)
410  ==
411  - rhoSp_
412  );
413 
414  deltaEqn.solve();
415 
416  if (nonOrth == nNonOrthCorr_)
417  {
418  phiAdd +=
419  fvc::interpolate(pp)
421  * regionMesh().magSf();
422 
423  phi_ == deltaEqn.flux();
424  }
425  }
426 
427  // Bound film thickness by a minimum of zero
428  delta_.max(0.0);
429 
430  // Update U field
431  U_ -= fvc::reconstruct(deltarUAf*phiAdd);
432 
433  // Remove any patch-normal components of velocity
434  U_ -= nHat()*(nHat() & U_);
435 
437 
438  // Continuity check
439  continuityCheck();
440 }
441 
442 
443 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
444 
445 kinematicSingleLayer::kinematicSingleLayer
446 (
447  const word& modelType,
448  const fvMesh& mesh,
449  const dimensionedVector& g,
450  const word& regionType,
451  const bool readFields
452 )
453 :
454  surfaceFilmRegionModel(modelType, mesh, g, regionType),
455 
456  momentumPredictor_(solution().subDict("PISO").lookup("momentumPredictor")),
457  nOuterCorr_(solution().subDict("PISO").lookupOrDefault("nOuterCorr", 1)),
458  nCorr_(readLabel(solution().subDict("PISO").lookup("nCorr"))),
460  (
461  readLabel(solution().subDict("PISO").lookup("nNonOrthCorr"))
462  ),
463 
464  cumulativeContErr_(0.0),
465 
466  deltaSmall_("deltaSmall", dimLength, small),
467  deltaCoLimit_(solution().lookupOrDefault("deltaCoLimit", 1e-4)),
468 
469  rho_
470  (
471  IOobject
472  (
473  "rhof",
474  time().timeName(),
475  regionMesh(),
478  ),
479  regionMesh(),
480  dimensionedScalar("zero", dimDensity, 0.0),
481  zeroGradientFvPatchScalarField::typeName
482  ),
483  mu_
484  (
485  IOobject
486  (
487  "muf",
488  time().timeName(),
489  regionMesh(),
492  ),
493  regionMesh(),
494  dimensionedScalar("zero", dimPressure*dimTime, 0.0),
495  zeroGradientFvPatchScalarField::typeName
496  ),
497  sigma_
498  (
499  IOobject
500  (
501  "sigmaf",
502  time().timeName(),
503  regionMesh(),
506  ),
507  regionMesh(),
508  dimensionedScalar("zero", dimMass/sqr(dimTime), 0.0),
509  zeroGradientFvPatchScalarField::typeName
510  ),
511 
512  delta_
513  (
514  IOobject
515  (
516  "deltaf",
517  time().timeName(),
518  regionMesh(),
521  ),
522  regionMesh()
523  ),
524  alpha_
525  (
526  IOobject
527  (
528  "alpha",
529  time().timeName(),
530  regionMesh(),
533  ),
534  regionMesh(),
535  dimensionedScalar("zero", dimless, 0.0),
536  zeroGradientFvPatchScalarField::typeName
537  ),
538  U_
539  (
540  IOobject
541  (
542  "Uf",
543  time().timeName(),
544  regionMesh(),
547  ),
548  regionMesh()
549  ),
550  Us_
551  (
552  IOobject
553  (
554  "Usf",
555  time().timeName(),
556  regionMesh(),
559  ),
560  U_,
561  zeroGradientFvPatchScalarField::typeName
562  ),
563  Uw_
564  (
565  IOobject
566  (
567  "Uwf",
568  time().timeName(),
569  regionMesh(),
572  ),
573  U_,
574  zeroGradientFvPatchScalarField::typeName
575  ),
576  deltaRho_
577  (
578  IOobject
579  (
580  delta_.name() + "*" + rho_.name(),
581  time().timeName(),
582  regionMesh(),
585  ),
586  regionMesh(),
588  zeroGradientFvPatchScalarField::typeName
589  ),
590 
591  phi_
592  (
593  IOobject
594  (
595  "phi",
596  time().timeName(),
597  regionMesh(),
600  ),
601  regionMesh(),
602  dimensionedScalar("0", dimLength*dimMass/dimTime, 0.0)
603  ),
604 
606  (
607  IOobject
608  (
609  "primaryMassTrans",
610  time().timeName(),
611  regionMesh(),
614  ),
615  regionMesh(),
616  dimensionedScalar("zero", dimMass, 0.0),
617  zeroGradientFvPatchScalarField::typeName
618  ),
620  (
621  IOobject
622  (
623  "cloudMassTrans",
624  time().timeName(),
625  regionMesh(),
628  ),
629  regionMesh(),
630  dimensionedScalar("zero", dimMass, 0.0),
631  zeroGradientFvPatchScalarField::typeName
632  ),
634  (
635  IOobject
636  (
637  "cloudDiameterTrans",
638  time().timeName(),
639  regionMesh(),
642  ),
643  regionMesh(),
644  dimensionedScalar("zero", dimLength, -1.0),
645  zeroGradientFvPatchScalarField::typeName
646  ),
647 
648  USp_
649  (
650  IOobject
651  (
652  "USpf",
653  time().timeName(),
654  regionMesh(),
657  ),
658  regionMesh(),
660  (
661  "zero", dimMass*dimVelocity/dimArea/dimTime, Zero
662  ),
663  this->mappedPushedFieldPatchTypes<vector>()
664  ),
665  pSp_
666  (
667  IOobject
668  (
669  "pSpf",
670  time_.timeName(),
671  regionMesh(),
674  ),
675  regionMesh(),
676  dimensionedScalar("zero", dimPressure, 0.0),
677  this->mappedPushedFieldPatchTypes<scalar>()
678  ),
679  rhoSp_
680  (
681  IOobject
682  (
683  "rhoSpf",
684  time_.timeName(),
685  regionMesh(),
688  ),
689  regionMesh(),
690  dimensionedScalar("zero", dimMass/dimTime/dimArea, 0.0),
691  this->mappedPushedFieldPatchTypes<scalar>()
692  ),
693 
695  (
696  IOobject
697  (
698  USp_.name(), // must have same name as USp_ to enable mapping
699  time().timeName(),
700  primaryMesh(),
703  ),
704  primaryMesh(),
706  ),
708  (
709  IOobject
710  (
711  pSp_.name(), // must have same name as pSp_ to enable mapping
712  time().timeName(),
713  primaryMesh(),
716  ),
717  primaryMesh(),
718  dimensionedScalar("zero", pSp_.dimensions(), 0.0)
719  ),
721  (
722  IOobject
723  (
724  rhoSp_.name(), // must have same name as rhoSp_ to enable mapping
725  time().timeName(),
726  primaryMesh(),
729  ),
730  primaryMesh(),
731  dimensionedScalar("zero", rhoSp_.dimensions(), 0.0)
732  ),
733 
734  UPrimary_
735  (
736  IOobject
737  (
738  "U", // must have same name as U to enable mapping
739  time().timeName(),
740  regionMesh(),
743  ),
744  regionMesh(),
746  this->mappedFieldAndInternalPatchTypes<vector>()
747  ),
748  pPrimary_
749  (
750  IOobject
751  (
752  "p", // must have same name as p to enable mapping
753  time().timeName(),
754  regionMesh(),
757  ),
758  regionMesh(),
759  dimensionedScalar("zero", dimPressure, 0.0),
760  this->mappedFieldAndInternalPatchTypes<scalar>()
761  ),
763  (
764  IOobject
765  (
766  "rho", // must have same name as rho to enable mapping
767  time().timeName(),
768  regionMesh(),
771  ),
772  regionMesh(),
773  dimensionedScalar("zero", dimDensity, 0.0),
774  this->mappedFieldAndInternalPatchTypes<scalar>()
775  ),
776  muPrimary_
777  (
778  IOobject
779  (
780  "thermo:mu", // must have same name as mu to enable mapping
781  time().timeName(),
782  regionMesh(),
785  ),
786  regionMesh(),
787  dimensionedScalar("zero", dimPressure*dimTime, 0.0),
788  this->mappedFieldAndInternalPatchTypes<scalar>()
789  ),
790 
792 
793  availableMass_(regionMesh().nCells(), 0.0),
794 
795  injection_(*this, coeffs_),
796 
797  transfer_(*this, coeffs_),
798 
800 
801  forces_(*this, coeffs_),
802 
803  addedMassTotal_(0.0)
804 {
805  if (readFields)
806  {
808 
809  correctAlpha();
810 
812 
813  deltaRho_ == delta_*rho_;
814 
816  (
817  IOobject
818  (
819  "phi",
820  time().timeName(),
821  regionMesh(),
824  false
825  ),
827  );
828 
829  phi_ == phi0;
830  }
831 }
832 
833 
834 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
835 
837 {}
838 
839 
840 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
841 
843 (
844  const label patchi,
845  const label facei,
846  const scalar massSource,
847  const vector& momentumSource,
848  const scalar pressureSource,
849  const scalar energySource
850 )
851 {
852  if (debug)
853  {
855  << "\nSurface film: " << type() << ": adding to film source:" << nl
856  << " mass = " << massSource << nl
857  << " momentum = " << momentumSource << nl
858  << " pressure = " << pressureSource << endl;
859  }
860 
861  rhoSpPrimary_.boundaryFieldRef()[patchi][facei] -= massSource;
862  USpPrimary_.boundaryFieldRef()[patchi][facei] -= momentumSource;
863  pSpPrimary_.boundaryFieldRef()[patchi][facei] -= pressureSource;
864 
865  addedMassTotal_ += massSource;
866 }
867 
868 
870 {
871  if (debug)
872  {
873  InfoInFunction << endl;
874  }
875 
877 
879 
881 
883 
884  // Reset transfer fields
885  availableMass_ = mass();
886  cloudMassTrans_ == dimensionedScalar("zero", dimMass, 0.0);
889 }
890 
891 
893 {
894  if (debug)
895  {
896  InfoInFunction << endl;
897  }
898 
899  // Update film coverage indicator
900  correctAlpha();
901 
902  // Update film wall and surface velocities
904 
905  // Update sub-models to provide updated source contributions
906  updateSubmodels();
907 
908  // Solve continuity for deltaRho_
909  solveContinuity();
910 
911  // Implicit pressure source coefficient - constant
912  tmp<volScalarField> tpp(this->pp());
913 
914  for (int oCorr=1; oCorr<=nOuterCorr_; oCorr++)
915  {
916  // Explicit pressure source contribution - varies with delta_
917  tmp<volScalarField> tpu(this->pu());
918 
919  // Solve for momentum for U_
920  tmp<fvVectorMatrix> UEqn = solveMomentum(tpu(), tpp());
921 
922  // Film thickness correction loop
923  for (int corr=1; corr<=nCorr_; corr++)
924  {
925  // Solve thickness for delta_
926  solveThickness(tpu(), tpp(), UEqn());
927  }
928  }
929 
930  // Update deltaRho_ with new delta_
931  deltaRho_ == delta_*rho_;
932 
933  // Reset source terms for next time integration
935 }
936 
937 
939 {
940  scalar CoNum = 0.0;
941 
942  if (regionMesh().nInternalFaces() > 0)
943  {
944  const scalarField sumPhi
945  (
946  fvc::surfaceSum(mag(phi_))().primitiveField()
947  / (deltaRho_.primitiveField() + rootVSmall)
948  );
949 
950  forAll(delta_, i)
951  {
952  if (delta_[i] > deltaCoLimit_)
953  {
954  CoNum = max(CoNum, sumPhi[i]/(delta_[i]*magSf()[i]));
955  }
956  }
957 
958  CoNum *= 0.5*time_.deltaTValue();
959  }
960 
961  reduce(CoNum, maxOp<scalar>());
962 
963  Info<< "Film max Courant number: " << CoNum << endl;
964 
965  return CoNum;
966 }
967 
968 
970 {
971  return U_;
972 }
973 
974 
976 {
977  return Us_;
978 }
979 
980 
982 {
983  return Uw_;
984 }
985 
986 
988 {
989  return deltaRho_;
990 }
991 
992 
994 {
995  return phi_;
996 }
997 
998 
1000 {
1001  return rho_;
1002 }
1003 
1004 
1006 {
1008  << "T field not available for " << type() << abort(FatalError);
1009 
1010  return volScalarField::null();
1011 }
1012 
1013 
1015 {
1017  << "Ts field not available for " << type() << abort(FatalError);
1018 
1019  return volScalarField::null();
1020 }
1021 
1022 
1024 {
1026  << "Tw field not available for " << type() << abort(FatalError);
1027 
1028  return volScalarField::null();
1029 }
1030 
1031 
1033 {
1035  << "hs field not available for " << type() << abort(FatalError);
1036 
1037  return volScalarField::null();
1038 }
1039 
1040 
1042 {
1044  << "Cp field not available for " << type() << abort(FatalError);
1045 
1046  return volScalarField::null();
1047 }
1048 
1049 
1051 {
1053  << "kappa field not available for " << type() << abort(FatalError);
1054 
1055  return volScalarField::null();
1056 }
1057 
1058 
1060 {
1061  return primaryMassTrans_;
1062 }
1063 
1064 
1066 {
1067  return cloudMassTrans_;
1068 }
1069 
1070 
1072 {
1073  return cloudDiameterTrans_;
1074 }
1075 
1076 
1078 {
1079  Info<< "\nSurface film: " << type() << endl;
1080 
1081  const scalarField& deltaInternal = delta_;
1082  const vectorField& Uinternal = U_;
1083  scalar addedMassTotal = 0.0;
1084  outputProperties().readIfPresent("addedMassTotal", addedMassTotal);
1085  addedMassTotal += returnReduce(addedMassTotal_, sumOp<scalar>());
1086 
1087  Info<< indent << "added mass = " << addedMassTotal << nl
1088  << indent << "current mass = "
1089  << gSum((deltaRho_*magSf())()) << nl
1090  << indent << "min/max(mag(U)) = " << gMin(mag(Uinternal)) << ", "
1091  << gMax(mag(Uinternal)) << nl
1092  << indent << "min/max(delta) = " << gMin(deltaInternal) << ", "
1093  << gMax(deltaInternal) << nl
1094  << indent << "coverage = "
1095  << gSum(alpha_.primitiveField()*magSf())/gSum(magSf()) << nl;
1096 
1097  injection_.info(Info);
1098  transfer_.info(Info);
1099 }
1100 
1101 
1103 {
1105  (
1107  (
1108  IOobject
1109  (
1110  typeName + ":Srho",
1111  time().timeName(),
1112  primaryMesh(),
1115  false
1116  ),
1117  primaryMesh(),
1118  dimensionedScalar("zero", dimMass/dimVolume/dimTime, 0.0)
1119  )
1120  );
1121 }
1122 
1123 
1126  const label i
1127 ) const
1128 {
1130  (
1132  (
1133  IOobject
1134  (
1135  typeName + ":Srho(" + Foam::name(i) + ")",
1136  time().timeName(),
1137  primaryMesh(),
1140  false
1141  ),
1142  primaryMesh(),
1143  dimensionedScalar("zero", dimMass/dimVolume/dimTime, 0.0)
1144  )
1145  );
1146 }
1147 
1148 
1150 {
1152  (
1154  (
1155  IOobject
1156  (
1157  typeName + ":Sh",
1158  time().timeName(),
1159  primaryMesh(),
1162  false
1163  ),
1164  primaryMesh(),
1165  dimensionedScalar("zero", dimEnergy/dimVolume/dimTime, 0.0)
1166  )
1167  );
1168 }
1169 
1170 
1171 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1172 
1173 } // End namespace surfaceFilmModels
1174 } // End namespace regionModels
1175 } // End namespace Foam
1176 
1177 // ************************************************************************* //
virtual bool read()
Read control parameters from dictionary.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:424
autoPtr< filmThermoModel > filmThermo_
Film thermo model.
virtual bool read()
Read control parameters from dictionary.
virtual const volVectorField & U() const
Return the film velocity [m/s].
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
autoPtr< filmTurbulenceModel > turbulence_
Turbulence model.
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
const word & name() const
Return name.
Definition: IOobject.H:297
virtual void correct(scalarField &availableMass, volScalarField &massToInject, volScalarField &diameterToInject)
Correct.
Ostream & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:226
virtual const volScalarField & deltaRho() const
Return the film thickness*density (helper field) [kg/m3].
virtual const volScalarField & Ts() const
Return the film surface temperature [K].
volScalarField primaryMassTrans_
Film mass available for transfer to the primary region.
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
Type gMin(const FieldField< Field, Type > &f)
const Boundary & boundaryField() const
Return const-reference to the boundary field.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
Reconstruct volField from a face flux field.
virtual void transferPrimaryRegionThermoFields()
Transfer thermo fields from the primary region to the film region.
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
tmp< GeometricField< Type, fvPatchField, volMesh > > H() const
Return the H operation source.
Definition: fvMatrix.C:766
virtual const volVectorField & Us() const
Return the film surface velocity [m/s].
scalar globalContErr
const dimensionedScalar phi0
Magnetic flux quantum: default SI units: [Wb].
void setFluxRequired(const word &name) const
Definition: fvSchemes.C:497
static autoPtr< filmThermoModel > New(surfaceFilmRegionModel &film, const dictionary &dict)
Return a reference to the selected phase change model.
tmp< fvVectorMatrix > tUEqn(fvm::ddt(rho, U)+fvm::div(phi, U)+MRF.DDt(rho, U)+turbulence->divDevRhoReff(U)==fvOptions(rho, U))
dimensioned< Type > domainIntegrate(const GeometricField< Type, fvPatchField, volMesh > &vf)
virtual tmp< volScalarField::Internal > Srho() const
Return total mass source - Eulerian phase only.
Calculate the snGrad of the given volField.
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
scalarField availableMass_
Available mass for transfer via sub-models.
volScalarField cloudDiameterTrans_
Parcel diameters originating from film to cloud.
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:626
scalar addedMassTotal_
Cumulative mass added via sources [kg].
const Time & time() const
Return the reference to the time database.
Definition: regionModelI.H:37
DimensionedField< scalar, volMesh > Internal
Type of the internal field from which this GeometricField is derived.
tmp< GeometricField< Type, fvPatchField, volMesh > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcLaplacian.C:45
virtual void info(Ostream &os)
Provide some info.
virtual tmp< volScalarField > primaryMassTrans() const
Return mass transfer source - Eulerian phase only.
addToRunTimeSelectionTable(surfaceFilmRegionModel, kinematicSingleLayer, mesh)
virtual void solveThickness(const volScalarField &pu, const volScalarField &pp, const fvVectorMatrix &UEqn)
Solve coupled velocity-thickness equations.
virtual void evolveRegion()
Evolve the film equations.
Macros for easy insertion into run-time selection tables.
virtual void solveContinuity()
Solve continuity equation.
surfaceScalarField phid("phid", fvc::interpolate(psi) *(fvc::flux(HbyA)+MRF.zeroFilter(rhorAUf *fvc::ddtCorr(rho, U, phi)/fvc::interpolate(rho))))
bool add(entry *, bool mergeEntry=false)
Add a new entry.
Definition: dictionary.C:814
virtual tmp< volScalarField > pp()
Implicit pressure source coefficient.
tmp< fvVectorMatrix > correct(volVectorField &U)
Return (net) force system.
Definition: forceList.C:78
volScalarField deltaRho_
Film thickness*density (helper field) [kg/m2].
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:692
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
virtual void transferPrimaryRegionSourceFields()
Transfer source fields from the primary region to the film region.
virtual void correct(scalarField &availableMass, volScalarField &massToTransfer)
Correct kinematic transfers.
dimensionedScalar pos(const dimensionedScalar &ds)
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:58
const dimensionSet & dimensions() const
Return dimensions.
const fvMesh & primaryMesh() const
Return the reference to the primary mesh database.
Definition: regionModelI.H:31
Type gSum(const FieldField< Field, Type > &f)
virtual tmp< fvVectorMatrix > solveMomentum(const volScalarField &pu, const volScalarField &pp)
Solve for film velocity.
virtual const volScalarField & cloudMassTrans() const
Return the film mass available for transfer to cloud.
virtual void preEvolveRegion()
Pre-evolve region.
Definition: regionModel.C:519
static autoPtr< filmTurbulenceModel > New(surfaceFilmRegionModel &film, const dictionary &dict)
Return a reference to the selected injection model.
const labelUList & faceCells() const
Return face-cell addressing.
Definition: polyPatch.C:315
tmp< volVectorField > gTan() const
Return the gravity tangential component contributions.
tmp< volScalarField > mass() const
Return the current film mass.
A class for handling words, derived from string.
Definition: word.H:59
SolverPerformance< Type > solve(fvMatrix< Type > &, const dictionary &)
Solve returning the solution statistics given convergence tolerance.
virtual const volScalarField & cloudDiameterTrans() const
Return the parcel diameters originating from film to cloud.
Calculate the face-flux of the given field.
Calculate the laplacian of the given field.
tmp< volScalarField > gNormClipped() const
Return the gravity normal-to-patch component contribution.
const IOdictionary & outputProperties() const
Return const access to the output properties dictionary.
Definition: regionModelI.H:110
const dictionary & solution() const
Return the solution dictionary.
Definition: regionModelI.H:103
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvmDdt.C:46
scalar deltaTValue() const
Return time step value.
Definition: TimeStateI.H:41
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
virtual const volVectorField & nHat() const
Return the patch normal vectors.
virtual const surfaceScalarField & phi() const
Return the film flux [kg.m/s].
const fvMesh & regionMesh() const
Return the region mesh database.
Definition: regionModelI.H:61
Volume integrate volField creating a volField.
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
Definition: fvPatchField.H:72
word timeName
Definition: getTimeIndex.H:3
static const zero Zero
Definition: zero.H:97
const dimensionSet dimPressure
label nNonOrthCorr_
Number of non-orthogonal correctors.
errorManip< error > abort(error &err)
Definition: errorManip.H:131
label readLabel(Istream &is)
Definition: label.H:64
virtual const volScalarField & kappa() const
Return the film thermal conductivity [W/m/K].
virtual tmp< volScalarField > pu()
Explicit pressure source contribution.
Calculate the divergence of the given field.
scalar CoNum
dimensionedScalar totalMass
Definition: continuityErrs.H:4
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
SolverPerformance< Type > solve(const dictionary &)
Solve segregated or coupled returning the solution statistics.
Definition: fvMatrixSolve.C:57
virtual void info(Ostream &os)
Provide some info.
static const char nl
Definition: Ostream.H:265
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmDiv.C:46
Type gMax(const FieldField< Field, Type > &f)
virtual void correctAlpha()
Correct film coverage field.
virtual const volScalarField & hs() const
Return the film surface enthalpy [J/kg].
virtual void resetPrimaryRegionSourceTerms()
Reset source term fields.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
fileName::Type type(const fileName &, const bool followLink=true)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:481
const dimensionSet dimEnergy
transferModelList transfer_
Transfer with the continuous phase.
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
const dimensionSet dimDensity
dictionary coeffs_
Model coefficients dictionary.
Definition: regionModel.H:102
static const GeometricField< scalar, fvPatchField, volMesh > & null()
Return a null geometric field.
virtual void updateSubmodels()
Update the film sub-models.
virtual const volScalarField & magSf() const
Return the face area magnitudes / [m2].
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> reconstruct(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
tmp< volScalarField > A() const
Return the central coefficient.
Definition: fvMatrix.C:737
label patchi
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
virtual const volVectorField & Uw() const
Return the film wall velocity [m/s].
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmLaplacian.C:46
const Time & time_
Reference to the time database.
Definition: regionModel.H:87
void max(const dimensioned< Type > &)
surfaceScalarField phi_
Mass flux (includes film thickness) [kg.m/s].
fvVectorMatrix & UEqn
Definition: UEqn.H:13
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
volScalarField alpha_
Film coverage indicator, 1 = covered, 0 = uncovered [].
A List with indirect addressing.
Definition: fvMatrix.H:106
virtual void updateSurfaceVelocities()
Update film surface velocities.
virtual const volScalarField & Tw() const
Return the film wall temperature [K].
void correctBoundaryConditions()
Correct boundary field.
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
tmp< GeometricField< Type, fvPatchField, volMesh > > surfaceSum(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
virtual void addSources(const label patchi, const label facei, const scalar massSource, const vector &momentumSource, const scalar pressureSource, const scalar energySource=0)
External hook to add sources to the film.
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > flux() const
Return the face-flux field from the matrix.
Definition: fvMatrix.C:876
Selector class for relaxation factors, solver type and solution.
Definition: solution.H:48
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:98
virtual const volScalarField & rho() const
Return the film density [kg/m3].
volScalarField cloudMassTrans_
Film mass available for transfer to cloud.
void constrainFilmField(Type &field, const typename Type::cmptType &value)
Constrain a film region master/slave boundaries of a field to a.
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
labelList intCoupledPatchIDs_
List of patch IDs internally coupled with the primary region.
Definition: regionModel.H:114
A class for managing temporary objects.
Definition: PtrList.H:53
tmp< surfaceScalarField > flux(const volVectorField &vvf)
Return the face-flux field obtained from the given volVectorField.
Definition: fvcFlux.C:32
virtual const volScalarField & Cp() const
Return the film specific heat capacity [J/kg/K].
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
virtual const volScalarField & T() const
Return the film mean temperature [K].
defineTypeNameAndDebug(kinematicSingleLayer, 0)
virtual tmp< volScalarField::Internal > Sh() const
Return enthalpy source - Eulerian phase only.
virtual scalar CourantNumber() const
Courant number evaluation.
surfaceScalarField rhof(fvc::interpolate(rho, "div(phi,rho)"))
dimensionedScalar deltaT() const
Return time step.
Definition: TimeStateI.H:53
scalar deltaCoLimit_
Film thickness above which Courant number calculation in valid.
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcSnGrad.C:45
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:57
scalar sumLocalContErr
Namespace for OpenFOAM.
virtual void correctThermoFields()
Correct the thermo fields.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:576
#define InfoInFunction
Report an information message using Foam::Info.
const dimensionSet dimVelocity