kinematicSingleLayer.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2017 OpenFOAM Foundation
9  Copyright (C) 2019-2020 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "kinematicSingleLayer.H"
30 #include "fvm.H"
31 #include "fvcDiv.H"
32 #include "fvcLaplacian.H"
33 #include "fvcSnGrad.H"
34 #include "fvcReconstruct.H"
35 #include "fvcVolumeIntegrate.H"
36 #include "fvcFlux.H"
38 #include "mappedWallPolyPatch.H"
39 #include "mapDistribute.H"
40 #include "filmThermoModel.H"
41 
42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
43 
44 namespace Foam
45 {
46 namespace regionModels
47 {
48 namespace surfaceFilmModels
49 {
50 
51 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
52 
53 defineTypeNameAndDebug(kinematicSingleLayer, 0);
54 
55 addToRunTimeSelectionTable(surfaceFilmRegionModel, kinematicSingleLayer, mesh);
56 
57 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
58 
60 {
62  {
63  const dictionary& solution = this->solution().subDict("PISO");
64  solution.readEntry("momentumPredictor", momentumPredictor_);
65  solution.readIfPresent("nOuterCorr", nOuterCorr_);
66  solution.readEntry("nCorr", nCorr_);
67  solution.readEntry("nNonOrthCorr", nNonOrthCorr_);
68 
69  return true;
70  }
71 
72  return false;
73 }
74 
75 
77 {
78  rho_ == filmThermo_->rho();
79  mu_ == filmThermo_->mu();
80  sigma_ == filmThermo_->sigma();
81 }
82 
83 
85 {
87 
88  rhoSpPrimary_ == dimensionedScalar(rhoSp_.dimensions(), Zero);
89  USpPrimary_ == dimensionedVector(USp_.dimensions(), Zero);
90  pSpPrimary_ == dimensionedScalar(pSp_.dimensions(), Zero);
91 }
92 
93 
95 {
97 
98  // Update fields from primary region via direct mapped
99  // (coupled) boundary conditions
104 }
105 
106 
108 {
110 
111  volScalarField::Boundary& rhoSpPrimaryBf =
113 
114  volVectorField::Boundary& USpPrimaryBf =
116 
117  volScalarField::Boundary& pSpPrimaryBf =
119 
120  // Convert accumulated source terms into per unit area per unit time
121  const scalar deltaT = time_.deltaTValue();
123  {
124  scalarField rpriMagSfdeltaT
125  (
126  (1.0/deltaT)
127  /primaryMesh().magSf().boundaryField()[patchi]
128  );
129 
130  rhoSpPrimaryBf[patchi] *= rpriMagSfdeltaT;
131  USpPrimaryBf[patchi] *= rpriMagSfdeltaT;
132  pSpPrimaryBf[patchi] *= rpriMagSfdeltaT;
133  }
134 
135  // Retrieve the source fields from the primary region via direct mapped
136  // (coupled) boundary conditions
137  // - fields require transfer of values for both patch AND to push the
138  // values into the first layer of internal cells
142 
143  // update addedMassTotal counter
144  if (time().writeTime())
145  {
146  if (debug)
147  {
148  rhoSp_.write();
149  USp_.write();
150  pSp_.write();
151  }
152 
153  scalar addedMassTotal = 0.0;
154  outputProperties().readIfPresent("addedMassTotal", addedMassTotal);
155  addedMassTotal += returnReduce(addedMassTotal_, sumOp<scalar>());
156  outputProperties().add("addedMassTotal", addedMassTotal, true);
157  addedMassTotal_ = 0.0;
158  }
159 }
160 
161 
163 {
164  return tmp<volScalarField>
165  (
166  new volScalarField
167  (
168  IOobject
169  (
170  typeName + ":pu",
171  time_.timeName(),
172  regionMesh(),
175  ),
176  pPrimary_ // pressure (mapped from primary region)
177  - pSp_ // accumulated particle impingement
178  - fvc::laplacian(sigma_, delta_) // surface tension
179  )
180  );
181 }
182 
183 
185 {
186  return tmp<volScalarField>
187  (
188  new volScalarField
189  (
190  IOobject
191  (
192  typeName + ":pp",
193  time_.timeName(),
194  regionMesh(),
197  ),
198  -rho_*gNormClipped() // hydrostatic effect only
199  )
200  );
201 }
202 
203 
205 {
207 }
208 
209 
211 {
213 
214  // Update injection model - mass returned is mass available for injection
216 
217  // Update transfer model - mass returned is mass available for transfer
219 
220  // Update mass source field
222 
223  turbulence_->correct();
224 }
225 
226 
228 {
229  const volScalarField deltaRho0(deltaRho_);
230 
231  solveContinuity();
232 
233  if (debug)
234  {
238  + dimensionedScalar("SMALL", dimMass*dimVolume, ROOTVSMALL);
239 
240  const scalar sumLocalContErr =
241  (
243  ).value();
244 
245  const scalar globalContErr =
246  (
248  ).value();
249 
251 
253  << "Surface film: " << type() << nl
254  << " time step continuity errors: sum local = "
255  << sumLocalContErr << ", global = " << globalContErr
256  << ", cumulative = " << cumulativeContErr_ << endl;
257  }
258 }
259 
260 
262 {
264 
265  solve
266  (
268  + fvc::div(phi_)
269  ==
270  - rhoSp_
271  );
272 }
273 
274 
276 {
277  // Push boundary film velocity values into internal field
278  for (label i=0; i<intCoupledPatchIDs_.size(); i++)
279  {
280  label patchi = intCoupledPatchIDs_[i];
281  const polyPatch& pp = regionMesh().boundaryMesh()[patchi];
282  UIndirectList<vector>(Uw_, pp.faceCells()) =
283  U_.boundaryField()[patchi];
284  }
285  Uw_ -= nHat()*(Uw_ & nHat());
287 
288  Us_ = turbulence_->Us();
289 }
290 
291 
293 (
294  const volScalarField& pu,
295  const volScalarField& pp
296 )
297 {
299 
300  // Momentum
302  (
303  fvm::ddt(deltaRho_, U_)
304  + fvm::div(phi_, U_)
305  ==
306  - USp_
307  // - fvm::SuSp(rhoSp_, U_)
308  - rhoSp_*U_
309  + forces_.correct(U_)
310  + turbulence_->Su(U_)
311  );
312 
313  fvVectorMatrix& UEqn = tUEqn.ref();
314 
315  UEqn.relax();
316 
317  if (momentumPredictor_)
318  {
319  solve
320  (
321  UEqn
322  ==
324  (
325  - fvc::interpolate(delta_)
326  * (
327  regionMesh().magSf()
328  * (
329  fvc::snGrad(pu, "snGrad(p)")
330  + fvc::snGrad(pp, "snGrad(p)")*fvc::interpolate(delta_)
331  + fvc::snGrad(delta_)*fvc::interpolate(pp)
332  )
333  - fvc::flux(rho_*gTan())
334  )
335  )
336  );
337 
338  // Remove any patch-normal components of velocity
339  U_ -= nHat()*(nHat() & U_);
340  U_.correctBoundaryConditions();
341  }
342 
343  return tUEqn;
344 }
345 
346 
348 (
349  const volScalarField& pu,
350  const volScalarField& pp,
352 )
353 {
355 
356  volScalarField rUA(1.0/UEqn.A());
357  U_ = rUA*UEqn.H();
358 
359  surfaceScalarField deltarUAf(fvc::interpolate(delta_*rUA));
361 
362  surfaceScalarField phiAdd
363  (
364  "phiAdd",
365  regionMesh().magSf()
366  * (
367  fvc::snGrad(pu, "snGrad(p)")
368  + fvc::snGrad(pp, "snGrad(p)")*fvc::interpolate(delta_)
369  )
370  - fvc::flux(rho_*gTan())
371  );
372  constrainFilmField(phiAdd, 0.0);
373 
375  (
376  "phid",
377  fvc::flux(U_*rho_) - deltarUAf*phiAdd*rhof
378  );
379  constrainFilmField(phid, 0.0);
380 
381  surfaceScalarField ddrhorUAppf
382  (
383  "deltaCoeff",
384  fvc::interpolate(delta_)*deltarUAf*rhof*fvc::interpolate(pp)
385  );
386 
387  regionMesh().setFluxRequired(delta_.name());
388 
389  for (int nonOrth=0; nonOrth<=nNonOrthCorr_; nonOrth++)
390  {
391  // Film thickness equation
392  fvScalarMatrix deltaEqn
393  (
394  fvm::ddt(rho_, delta_)
395  + fvm::div(phid, delta_)
396  - fvm::laplacian(ddrhorUAppf, delta_)
397  ==
398  - rhoSp_
399  );
400 
401  deltaEqn.solve();
402 
403  if (nonOrth == nNonOrthCorr_)
404  {
405  phiAdd +=
406  fvc::interpolate(pp)
407  * fvc::snGrad(delta_)
408  * regionMesh().magSf();
409 
410  phi_ == deltaEqn.flux();
411  }
412  }
413 
414  // Bound film thickness by a minimum of zero
415  delta_.max(0.0);
416 
417  // Update U field
418  U_ -= fvc::reconstruct(deltarUAf*phiAdd);
419 
420  // Remove any patch-normal components of velocity
421  U_ -= nHat()*(nHat() & U_);
422 
423  U_.correctBoundaryConditions();
424 
425  // Update film wall and surface velocities
426  updateSurfaceVelocities();
427 
428  // Continuity check
429  continuityCheck();
430 }
431 
432 
433 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
434 
435 kinematicSingleLayer::kinematicSingleLayer
436 (
437  const word& modelType,
438  const fvMesh& mesh,
439  const dimensionedVector& g,
440  const word& regionType,
441  const bool readFields
442 )
443 :
444  surfaceFilmRegionModel(modelType, mesh, g, regionType),
445 
446  momentumPredictor_(solution().subDict("PISO").lookup("momentumPredictor")),
447  nOuterCorr_(solution().subDict("PISO").getOrDefault("nOuterCorr", 1)),
448  nCorr_(solution().subDict("PISO").get<label>("nCorr")),
449  nNonOrthCorr_
450  (
451  solution().subDict("PISO").get<label>("nNonOrthCorr")
452  ),
453 
454  cumulativeContErr_(0.0),
455 
456  deltaSmall_("deltaSmall", dimLength, SMALL),
457  deltaCoLimit_(solution().getOrDefault<scalar>("deltaCoLimit", 1e-4)),
458 
459  rho_
460  (
461  IOobject
462  (
463  "rhof",
464  time().timeName(),
465  regionMesh(),
468  ),
469  regionMesh(),
471  zeroGradientFvPatchScalarField::typeName
472  ),
473  mu_
474  (
475  IOobject
476  (
477  "muf",
478  time().timeName(),
479  regionMesh(),
482  ),
483  regionMesh(),
485  zeroGradientFvPatchScalarField::typeName
486  ),
487  sigma_
488  (
489  IOobject
490  (
491  "sigmaf",
492  time().timeName(),
493  regionMesh(),
496  ),
497  regionMesh(),
499  zeroGradientFvPatchScalarField::typeName
500  ),
501 
502  delta_
503  (
504  IOobject
505  (
506  "deltaf",
507  time().timeName(),
508  regionMesh(),
511  ),
512  regionMesh()
513  ),
514  alpha_
515  (
516  IOobject
517  (
518  "alpha",
519  time().timeName(),
520  regionMesh(),
523  ),
524  regionMesh(),
526  zeroGradientFvPatchScalarField::typeName
527  ),
528  U_
529  (
530  IOobject
531  (
532  "Uf",
533  time().timeName(),
534  regionMesh(),
537  ),
538  regionMesh()
539  ),
540  Us_
541  (
542  IOobject
543  (
544  "Usf",
545  time().timeName(),
546  regionMesh(),
549  ),
550  U_,
551  zeroGradientFvPatchScalarField::typeName
552  ),
553  Uw_
554  (
555  IOobject
556  (
557  "Uwf",
558  time().timeName(),
559  regionMesh(),
562  ),
563  U_,
564  zeroGradientFvPatchScalarField::typeName
565  ),
566  deltaRho_
567  (
568  IOobject
569  (
570  delta_.name() + "*" + rho_.name(),
571  time().timeName(),
572  regionMesh(),
575  ),
576  regionMesh(),
577  dimensionedScalar(delta_.dimensions()*rho_.dimensions(), Zero),
578  zeroGradientFvPatchScalarField::typeName
579  ),
580 
581  phi_
582  (
583  IOobject
584  (
585  "phi",
586  time().timeName(),
587  regionMesh(),
590  ),
591  regionMesh(),
593  ),
594 
595  primaryMassTrans_
596  (
597  IOobject
598  (
599  "primaryMassTrans",
600  time().timeName(),
601  regionMesh(),
604  ),
605  regionMesh(),
607  zeroGradientFvPatchScalarField::typeName
608  ),
609  cloudMassTrans_
610  (
611  IOobject
612  (
613  "cloudMassTrans",
614  time().timeName(),
615  regionMesh(),
618  ),
619  regionMesh(),
621  zeroGradientFvPatchScalarField::typeName
622  ),
623  cloudDiameterTrans_
624  (
625  IOobject
626  (
627  "cloudDiameterTrans",
628  time().timeName(),
629  regionMesh(),
632  ),
633  regionMesh(),
634  dimensionedScalar("minus1", dimLength, -1.0),
635  zeroGradientFvPatchScalarField::typeName
636  ),
637 
638  USp_
639  (
640  IOobject
641  (
642  "USpf",
643  time().timeName(),
644  regionMesh(),
647  ),
648  regionMesh(),
650  this->mappedPushedFieldPatchTypes<vector>()
651  ),
652  pSp_
653  (
654  IOobject
655  (
656  "pSpf",
657  time_.timeName(),
658  regionMesh(),
661  ),
662  regionMesh(),
664  this->mappedPushedFieldPatchTypes<scalar>()
665  ),
666  rhoSp_
667  (
668  IOobject
669  (
670  "rhoSpf",
671  time_.timeName(),
672  regionMesh(),
675  ),
676  regionMesh(),
678  this->mappedPushedFieldPatchTypes<scalar>()
679  ),
680 
681  USpPrimary_
682  (
683  IOobject
684  (
685  USp_.name(), // must have same name as USp_ to enable mapping
686  time().timeName(),
687  primaryMesh(),
690  ),
691  primaryMesh(),
692  dimensionedVector(USp_.dimensions(), Zero)
693  ),
694  pSpPrimary_
695  (
696  IOobject
697  (
698  pSp_.name(), // must have same name as pSp_ to enable mapping
699  time().timeName(),
700  primaryMesh(),
703  ),
704  primaryMesh(),
705  dimensionedScalar(pSp_.dimensions(), Zero)
706  ),
707  rhoSpPrimary_
708  (
709  IOobject
710  (
711  rhoSp_.name(), // must have same name as rhoSp_ to enable mapping
712  time().timeName(),
713  primaryMesh(),
716  ),
717  primaryMesh(),
718  dimensionedScalar(rhoSp_.dimensions(), Zero)
719  ),
720 
721  UPrimary_
722  (
723  IOobject
724  (
725  "U", // must have same name as U to enable mapping
726  time().timeName(),
727  regionMesh(),
730  ),
731  regionMesh(),
733  this->mappedFieldAndInternalPatchTypes<vector>()
734  ),
735  pPrimary_
736  (
737  IOobject
738  (
739  "p", // must have same name as p to enable mapping
740  time().timeName(),
741  regionMesh(),
744  ),
745  regionMesh(),
747  this->mappedFieldAndInternalPatchTypes<scalar>()
748  ),
749  rhoPrimary_
750  (
751  IOobject
752  (
753  "rho", // must have same name as rho to enable mapping
754  time().timeName(),
755  regionMesh(),
758  ),
759  regionMesh(),
761  this->mappedFieldAndInternalPatchTypes<scalar>()
762  ),
763  muPrimary_
764  (
765  IOobject
766  (
767  "thermo:mu", // must have same name as mu to enable mapping
768  time().timeName(),
769  regionMesh(),
772  ),
773  regionMesh(),
775  this->mappedFieldAndInternalPatchTypes<scalar>()
776  ),
777 
778  filmThermo_(filmThermoModel::New(*this, coeffs_)),
779 
780  availableMass_(regionMesh().nCells(), Zero),
781 
782  injection_(*this, coeffs_),
783 
784  transfer_(*this, coeffs_),
785 
786  turbulence_(filmTurbulenceModel::New(*this, coeffs_)),
787 
788  forces_(*this, coeffs_),
789 
790  addedMassTotal_(0.0)
791 {
792  if (readFields)
793  {
794  transferPrimaryRegionThermoFields();
795 
796  correctAlpha();
797 
798  correctThermoFields();
799 
800  deltaRho_ == delta_*rho_;
801 
803  (
804  IOobject
805  (
806  "phi",
807  time().timeName(),
808  regionMesh(),
811  false
812  ),
813  fvc::flux(deltaRho_*U_)
814  );
815 
816  phi_ == phi0;
817  }
818 }
819 
820 
821 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
822 
824 {}
825 
826 
827 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
828 
830 (
831  const label patchi,
832  const label facei,
833  const scalar massSource,
834  const vector& momentumSource,
835  const scalar pressureSource,
836  const scalar energySource
837 )
838 {
840  << "\nSurface film: " << type() << ": adding to film source:" << nl
841  << " mass = " << massSource << nl
842  << " momentum = " << momentumSource << nl
843  << " pressure = " << pressureSource << endl;
844 
845  rhoSpPrimary_.boundaryFieldRef()[patchi][facei] -= massSource;
846  USpPrimary_.boundaryFieldRef()[patchi][facei] -= momentumSource;
847  pSpPrimary_.boundaryFieldRef()[patchi][facei] -= pressureSource;
848 
849  addedMassTotal_ += massSource;
850 }
851 
852 
854 {
856 
858 
860 
862 
864 
866 
867  correctAlpha();
868 
869  // Reset transfer fields
870  availableMass_ = mass();
874 }
875 
876 
878 {
880 
881  // Update sub-models to provide updated source contributions
882  updateSubmodels();
883 
884  // Solve continuity for deltaRho_
885  solveContinuity();
886 
887  // Implicit pressure source coefficient - constant
888  tmp<volScalarField> tpp(this->pp());
889 
890  for (int oCorr=1; oCorr<=nOuterCorr_; oCorr++)
891  {
892  // Explicit pressure source contribution - varies with delta_
893  tmp<volScalarField> tpu(this->pu());
894 
895  // Solve for momentum for U_
896  tmp<fvVectorMatrix> tUEqn = solveMomentum(tpu(), tpp());
897  fvVectorMatrix& UEqn = tUEqn.ref();
898 
899  // Film thickness correction loop
900  for (int corr=1; corr<=nCorr_; corr++)
901  {
902  // Solve thickness for delta_
903  solveThickness(tpu(), tpp(), UEqn);
904  }
905  }
906 
907  // Update deltaRho_ with new delta_
908  deltaRho_ == delta_*rho_;
909 }
910 
911 
913 {
915 
916  // Reset source terms for next time integration
918 }
919 
920 
922 {
923  scalar CoNum = 0.0;
924 
925  if (regionMesh().nInternalFaces() > 0)
926  {
927  const scalarField sumPhi
928  (
929  fvc::surfaceSum(mag(phi_))().primitiveField()
930  / (deltaRho_.primitiveField() + ROOTVSMALL)
931  );
932 
933  forAll(delta_, i)
934  {
935  if ((delta_[i] > deltaCoLimit_) && (alpha_[i] > 0.5))
936  {
937  CoNum = max(CoNum, sumPhi[i]/(delta_[i]*magSf()[i]));
938  }
939  }
940 
941  CoNum *= 0.5*time_.deltaTValue();
942  }
943 
945 
946  Info<< "Film max Courant number: " << CoNum << endl;
947 
948  return CoNum;
949 }
950 
951 
953 {
954  return U_;
955 }
956 
957 
959 {
960  return Us_;
961 }
962 
963 
965 {
966  return Uw_;
967 }
968 
969 
971 {
972  return deltaRho_;
973 }
974 
975 
977 {
978  return phi_;
979 }
980 
981 
983 {
984  return rho_;
985 }
986 
987 
989 {
991  << "T field not available for " << type() << abort(FatalError);
992 
993  return volScalarField::null();
994 }
995 
996 
998 {
1000  << "Ts field not available for " << type() << abort(FatalError);
1001 
1002  return volScalarField::null();
1003 }
1004 
1005 
1007 {
1009  << "Tw field not available for " << type() << abort(FatalError);
1010 
1011  return volScalarField::null();
1012 }
1013 
1014 
1016 {
1018  << "hs field not available for " << type() << abort(FatalError);
1019 
1020  return volScalarField::null();
1021 }
1022 
1023 
1025 {
1027  << "Cp field not available for " << type() << abort(FatalError);
1028 
1029  return volScalarField::null();
1030 }
1031 
1032 
1034 {
1036  << "kappa field not available for " << type() << abort(FatalError);
1037 
1038  return volScalarField::null();
1039 }
1040 
1041 
1043 {
1044  return primaryMassTrans_;
1045 }
1046 
1047 
1049 {
1050  return cloudMassTrans_;
1051 }
1052 
1053 
1055 {
1056  return cloudDiameterTrans_;
1057 }
1058 
1059 
1061 {
1062  Info<< "\nSurface film: " << type() << endl;
1063 
1064  const scalarField& deltaInternal = delta_;
1065  const vectorField& Uinternal = U_;
1066  scalar addedMassTotal = 0.0;
1067  outputProperties().readIfPresent("addedMassTotal", addedMassTotal);
1068  addedMassTotal += returnReduce(addedMassTotal_, sumOp<scalar>());
1069 
1070  Info<< indent << "added mass = " << addedMassTotal << nl
1071  << indent << "current mass = "
1072  << gSum((deltaRho_*magSf())()) << nl
1073  << indent << "min/max(mag(U)) = " << gMin(mag(Uinternal)) << ", "
1074  << gMax(mag(Uinternal)) << nl
1075  << indent << "min/max(delta) = " << gMin(deltaInternal) << ", "
1076  << gMax(deltaInternal) << nl
1077  << indent << "coverage = "
1078  << gSum(alpha_.primitiveField()*magSf())/gSum(magSf()) << nl;
1079 
1080  injection_.info(Info);
1081  transfer_.info(Info);
1082 }
1083 
1084 
1086 {
1088  (
1089  IOobject
1090  (
1091  typeName + ":Srho",
1092  time().timeName(),
1093  primaryMesh(),
1096  false
1097  ),
1098  primaryMesh(),
1100  );
1101 }
1102 
1103 
1106  const label i
1107 ) const
1108 {
1110  (
1111  IOobject
1112  (
1113  typeName + ":Srho(" + Foam::name(i) + ")",
1114  time().timeName(),
1115  primaryMesh(),
1118  false
1119  ),
1120  primaryMesh(),
1122  );
1123 }
1124 
1125 
1127 {
1129  (
1130  IOobject
1131  (
1132  typeName + ":Sh",
1133  time().timeName(),
1134  primaryMesh(),
1137  false
1138  ),
1139  primaryMesh(),
1141  );
1142 }
1143 
1144 
1145 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1146 
1147 } // End namespace surfaceFilmModels
1148 } // End namespace regionModels
1149 } // End namespace Foam
1150 
1151 // ************************************************************************* //
totalMass
dimensionedScalar totalMass
Definition: continuityErrs.H:4
Foam::expressions::patchExpr::debug
int debug
Static debugging option.
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::T
virtual const volScalarField & T() const
Return the film mean temperature [K].
Definition: kinematicSingleLayer.C:988
Foam::IOobject::NO_WRITE
Definition: IOobject.H:195
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::read
virtual bool read()
Read control parameters from dictionary.
Definition: kinematicSingleLayer.C:59
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::correctThermoFields
virtual void correctThermoFields()
Correct the thermo fields.
Definition: kinematicSingleLayer.C:76
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::rhoSpPrimary_
volScalarField rhoSpPrimary_
Mass [kg/m2/s].
Definition: kinematicSingleLayer.H:186
Foam::fvc::reconstruct
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> reconstruct(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcReconstruct.C:56
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::deltaCoLimit_
scalar deltaCoLimit_
Film thickness above which Courant number calculation in valid.
Definition: kinematicSingleLayer.H:107
Foam::dimPressure
const dimensionSet dimPressure
Foam::fvc::snGrad
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcSnGrad.C:47
Foam::fvc::surfaceSum
tmp< GeometricField< Type, fvPatchField, volMesh > > surfaceSum(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcSurfaceIntegrate.C:135
sumLocalContErr
scalar sumLocalContErr
Definition: compressibleContinuityErrs.H:34
Foam::maxOp
Definition: ops.H:223
phid
surfaceScalarField phid("phid", fvc::interpolate(psi) *(fvc::flux(HbyA)+MRF.zeroFilter(rhorAUf *fvc::ddtCorr(rho, U, phi)/fvc::interpolate(rho))))
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
InfoInFunction
#define InfoInFunction
Report an information message using Foam::Info.
Definition: messageStream.H:350
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::~kinematicSingleLayer
virtual ~kinematicSingleLayer()
Destructor.
Definition: kinematicSingleLayer.C:823
filmThermoModel.H
Foam::IOobject::AUTO_WRITE
Definition: IOobject.H:194
mappedWallPolyPatch.H
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::transferPrimaryRegionThermoFields
virtual void transferPrimaryRegionThermoFields()
Transfer thermo fields from the primary region to the film region.
Definition: kinematicSingleLayer.C:94
Foam::fvc::flux
tmp< surfaceScalarField > flux(const volVectorField &vvf)
Return the face-flux field obtained from the given volVectorField.
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::Us_
volVectorField Us_
Velocity - surface [m/s].
Definition: kinematicSingleLayer.H:136
Foam::solution
Selector class for relaxation factors, solver type and solution.
Definition: solution.H:55
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::USpPrimary_
volVectorField USpPrimary_
Momentum [kg/m/s2].
Definition: kinematicSingleLayer.H:180
Foam::returnReduce
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
Definition: PstreamReduceOps.H:94
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::USp_
volVectorField USp_
Momentum [kg/m/s2].
Definition: kinematicSingleLayer.H:167
Foam::TimeState::deltaT
dimensionedScalar deltaT() const
Return time step.
Definition: TimeStateI.H:55
Foam::regionModels::surfaceFilmModels::transferModelList::correct
virtual void correct(scalarField &availableMass, volScalarField &massToTransfer)
Correct kinematic transfers.
Definition: transferModelList.C:101
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::Us
virtual const volVectorField & Us() const
Return the film surface velocity [m/s].
Definition: kinematicSingleLayer.C:958
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::injection_
injectionModelList injection_
Cloud injection.
Definition: kinematicSingleLayer.H:214
Foam::dimLength
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:52
Foam::regionModels::surfaceFilmModels::transferModelList::info
virtual void info(Ostream &os)
Provide some info.
Definition: transferModelList.C:154
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::preEvolveRegion
virtual void preEvolveRegion()
Pre-evolve film hook.
Definition: kinematicSingleLayer.C:853
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::primaryMassTrans_
volScalarField primaryMassTrans_
Film mass available for transfer to the primary region.
Definition: kinematicSingleLayer.H:151
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::dimVelocity
const dimensionSet dimVelocity
Foam::dimEnergy
const dimensionSet dimEnergy
Foam::dimDensity
const dimensionSet dimDensity
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::Ts
virtual const volScalarField & Ts() const
Return the film surface temperature [K].
Definition: kinematicSingleLayer.C:997
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::turbulence_
autoPtr< filmTurbulenceModel > turbulence_
Turbulence model.
Definition: kinematicSingleLayer.H:220
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::filmThermo_
autoPtr< filmThermoModel > filmThermo_
Film thermo model.
Definition: kinematicSingleLayer.H:208
Foam::fvc::domainIntegrate
dimensioned< Type > domainIntegrate(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvcVolumeIntegrate.C:88
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::U
virtual const volVectorField & U() const
Return the film velocity [m/s].
Definition: kinematicSingleLayer.C:952
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::deltaSmall_
const dimensionedScalar deltaSmall_
Small delta.
Definition: kinematicSingleLayer.H:104
Foam::GeometricField::primitiveField
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
Definition: GeometricFieldI.H:53
fvcSnGrad.H
Calculate the snGrad of the given volField.
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::evolveRegion
virtual void evolveRegion()
Evolve the film equations.
Definition: kinematicSingleLayer.C:877
Foam::Time::timeName
static word timeName(const scalar t, const int precision=precision_)
Definition: Time.C:780
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::solveContinuity
virtual void solveContinuity()
Solve continuity equation.
Definition: kinematicSingleLayer.C:261
Foam::regionModels::surfaceFilmModels::injectionModelList::info
virtual void info(Ostream &os)
Provide some info.
Definition: injectionModelList.C:126
fvcDiv.H
Calculate the divergence of the given field.
Foam::regionModels::regionModel::time
const Time & time() const
Return the reference to the time database.
Definition: regionModelI.H:40
Foam::fvc::div
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:49
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::rho_
volScalarField rho_
Density [kg/m3].
Definition: kinematicSingleLayer.H:115
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::sigma_
volScalarField sigma_
Surface tension [m/s2].
Definition: kinematicSingleLayer.H:121
Foam::polyMesh::boundaryMesh
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:444
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::continuityCheck
virtual void continuityCheck()
Continuity check.
Definition: kinematicSingleLayer.C:227
Foam::gSum
Type gSum(const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:594
Foam::regionModels::surfaceFilmModels::addToRunTimeSelectionTable
addToRunTimeSelectionTable(surfaceFilmRegionModel, kinematicSingleLayer, mesh)
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::availableMass_
scalarField availableMass_
Available mass for transfer via sub-models.
Definition: kinematicSingleLayer.H:211
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::pPrimary_
volScalarField pPrimary_
Pressure [Pa].
Definition: kinematicSingleLayer.H:196
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::Srho
virtual tmp< volScalarField::Internal > Srho() const
Return total mass source - Eulerian phase only.
Definition: kinematicSingleLayer.C:1085
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::muPrimary_
volScalarField muPrimary_
Viscosity [Pa.s].
Definition: kinematicSingleLayer.H:202
Foam::sumOp
Definition: ops.H:213
Foam::TimeState::deltaTValue
scalar deltaTValue() const noexcept
Return time step value.
Definition: TimeStateI.H:43
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::dimensionedVector
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
Definition: dimensionedVector.H:50
Foam::regionModels::regionModel::primaryMesh
const fvMesh & primaryMesh() const
Return the reference to the primary mesh database.
Definition: regionModelI.H:34
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::gNormClipped
tmp< volScalarField > gNormClipped() const
Return the gravity normal-to-patch component contribution.
Definition: kinematicSingleLayerI.H:234
Foam::fvc::laplacian
tmp< GeometricField< Type, fvPatchField, volMesh > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcLaplacian.C:47
Foam::regionModels::surfaceFilmModels::surfaceFilmRegionModel::read
virtual bool read()
Read control parameters from dictionary.
Definition: surfaceFilmRegionModel.C:46
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::cloudDiameterTrans_
volScalarField cloudDiameterTrans_
Parcel diameters originating from film to cloud.
Definition: kinematicSingleLayer.H:157
Foam::dimTime
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::deltaRho_
volScalarField deltaRho_
Film thickness*density (helper field) [kg/m2].
Definition: kinematicSingleLayer.H:142
Foam::reduce
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
Definition: PstreamReduceOps.H:51
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::cloudMassTrans
virtual const volScalarField & cloudMassTrans() const
Return the film mass available for transfer to cloud.
Definition: kinematicSingleLayer.C:1048
tUEqn
tmp< fvVectorMatrix > tUEqn(fvm::ddt(rho, U)+fvm::div(phi, U)+MRF.DDt(rho, U)+turbulence->divDevRhoReff(U)==fvOptions(rho, U))
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::addedMassTotal_
scalar addedMassTotal_
Cumulative mass added via sources [kg].
Definition: kinematicSingleLayer.H:229
Foam::dimArea
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:59
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::transferPrimaryRegionSourceFields
virtual void transferPrimaryRegionSourceFields()
Transfer source fields from the primary region to the film region.
Definition: kinematicSingleLayer.C:107
Foam::Field< scalar >
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::delta_
volScalarField delta_
Film thickness [m].
Definition: kinematicSingleLayer.H:127
Foam::regionModels::surfaceFilmModels::surfaceFilmRegionModel
Base class for surface film models.
Definition: surfaceFilmRegionModel.H:55
Foam::Info
messageStream Info
Information stream (stdout output on master, null elsewhere)
DebugInFunction
#define DebugInFunction
Report an information message using Foam::Info.
Definition: messageStream.H:388
Foam::polyPatch
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:68
Foam::fvm::laplacian
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmLaplacian.C:48
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::kappa
virtual const volScalarField & kappa() const
Return the film thermal conductivity [W/m/K].
Definition: kinematicSingleLayer.C:1033
Foam::regionModels::singleLayerRegion::magSf
virtual const volScalarField & magSf() const
Return the face area magnitudes / [m2].
Definition: singleLayerRegion.C:223
Foam::regionModels::regionModel::regionMesh
const fvMesh & regionMesh() const
Return the region mesh database.
Definition: regionModelI.H:64
Foam::dictionary::subDict
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
Definition: dictionary.C:460
Foam::IOobject::READ_IF_PRESENT
Definition: IOobject.H:187
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::pp
virtual tmp< volScalarField > pp()
Implicit pressure source coefficient.
Definition: kinematicSingleLayer.C:184
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::cloudDiameterTrans
virtual const volScalarField & cloudDiameterTrans() const
Return the parcel diameters originating from film to cloud.
Definition: kinematicSingleLayer.C:1054
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::info
virtual void info()
Provide some feedback.
Definition: kinematicSingleLayer.C:1060
Foam::regionModels::regionModel::preEvolveRegion
virtual void preEvolveRegion()
Pre-evolve region.
Definition: regionModel.C:500
Foam::constant::electromagnetic::phi0
const dimensionedScalar phi0
Magnetic flux quantum: default SI units: [Wb].
Foam::dictionary::readEntry
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
Definition: dictionaryTemplates.C:302
Foam::solve
SolverPerformance< Type > solve(faMatrix< Type > &, Istream &)
Solve returning the solution statistics given convergence tolerance.
Foam::dimensionedScalar
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Definition: dimensionedScalarFwd.H:42
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::Uw_
volVectorField Uw_
Velocity - wall [m/s].
Definition: kinematicSingleLayer.H:139
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::primaryMassTrans
virtual tmp< volScalarField > primaryMassTrans() const
Return mass transfer source - Eulerian phase only.
Definition: kinematicSingleLayer.C:1042
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::resetPrimaryRegionSourceTerms
virtual void resetPrimaryRegionSourceTerms()
Reset source term fields.
Definition: kinematicSingleLayer.C:84
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::solveMomentum
virtual tmp< fvVectorMatrix > solveMomentum(const volScalarField &pu, const volScalarField &pp)
Solve for film velocity.
Definition: kinematicSingleLayer.C:293
Foam::max
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
Foam::regionModels::regionModel::outputProperties
const IOdictionary & outputProperties() const
Return const access to the output properties dictionary.
Definition: regionModelI.H:106
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::phi
virtual const surfaceScalarField & phi() const
Return the film flux [kg.m/s].
Definition: kinematicSingleLayer.C:976
Foam::radiation::lookup
Lookup type of boundary radiation properties.
Definition: lookup.H:63
Foam::regionModels::regionModel::solution
const dictionary & solution() const
Return the solution dictionary.
Definition: regionModelI.H:99
timeName
word timeName
Definition: getTimeIndex.H:3
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::correctAlpha
virtual void correctAlpha()
Correct film coverage field.
Definition: kinematicSingleLayer.C:204
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::momentumPredictor_
Switch momentumPredictor_
Momentum predictor.
Definition: kinematicSingleLayer.H:89
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::hs
virtual const volScalarField & hs() const
Return the film surface enthalpy [J/kg].
Definition: kinematicSingleLayer.C:1015
fvm.H
Foam::FatalError
error FatalError
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
kinematicSingleLayer.H
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam::dimensioned< scalar >
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::U_
volVectorField U_
Velocity - mean [m/s].
Definition: kinematicSingleLayer.H:133
Foam::fvMesh
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:85
g
const uniformDimensionedVectorField & g
Definition: createFluidFields.H:26
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::dimMass
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:51
Foam::abort
errorManip< error > abort(error &err)
Definition: errorManip.H:144
Foam::regionModels::surfaceFilmModels::filmThermoModel::New
static autoPtr< filmThermoModel > New(surfaceFilmRegionModel &film, const dictionary &dict)
Return a reference to the selected phase change model.
Definition: filmThermoModelNew.C:43
Foam::fvc::interpolate
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.
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::nCorr_
label nCorr_
Number of PISO-like correctors.
Definition: kinematicSingleLayer.H:95
Foam::indent
Ostream & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:339
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::mass
tmp< volScalarField > mass() const
Return the current film mass.
Definition: kinematicSingleLayerI.H:200
CoNum
scalar CoNum
Definition: compressibleCourantNo.H:31
globalContErr
scalar globalContErr
Definition: compressibleContinuityErrs.H:37
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::pu
virtual tmp< volScalarField > pu()
Explicit pressure source contribution.
Definition: kinematicSingleLayer.C:162
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::UPrimary_
volVectorField UPrimary_
Velocity [m/s].
Definition: kinematicSingleLayer.H:193
Foam::readFields
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const wordHashSet &selectedFields, LIFOStack< regIOobject * > &storedObjects)
Read the selected GeometricFields of the templated type.
Definition: ReadFieldsTemplates.C:312
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::pSp_
volScalarField pSp_
Pressure [Pa].
Definition: kinematicSingleLayer.H:170
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::cumulativeContErr_
scalar cumulativeContErr_
Cumulative continuity error.
Definition: kinematicSingleLayer.H:101
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::phi_
surfaceScalarField phi_
Mass flux (includes film thickness) [kg.m/s].
Definition: kinematicSingleLayer.H:145
Foam::GeometricField::correctBoundaryConditions
void correctBoundaryConditions()
Correct boundary field.
Definition: GeometricField.C:940
Foam::regionModels::surfaceFilmModels::defineTypeNameAndDebug
defineTypeNameAndDebug(kinematicSingleLayer, 0)
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::addSources
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.
Definition: kinematicSingleLayer.C:830
Foam::fvm::ddt
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvmDdt.C:48
Foam::sqr
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Definition: dimensionedSymmTensor.C:51
Foam::regionModels::regionModel::time_
const Time & time_
Reference to the time database.
Definition: regionModel.H:90
fvcLaplacian.H
Calculate the laplacian of the given field.
Foam::nl
constexpr char nl
Definition: Ostream.H:404
mapDistribute.H
Foam::regionModels::surfaceFilmModels::injectionModelList::correct
virtual void correct(scalarField &availableMass, volScalarField &massToInject, volScalarField &diameterToInject)
Correct.
Definition: injectionModelList.C:98
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::alpha_
volScalarField alpha_
Film coverage indicator, 1 = covered, 0 = uncovered [].
Definition: kinematicSingleLayer.H:130
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::pSpPrimary_
volScalarField pSpPrimary_
Pressure [Pa].
Definition: kinematicSingleLayer.H:183
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::mu_
volScalarField mu_
Dynamic viscosity [Pa.s].
Definition: kinematicSingleLayer.H:118
Foam::GeometricField::boundaryFieldRef
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
Definition: GeometricField.C:783
Foam::Vector< scalar >
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::nNonOrthCorr_
label nNonOrthCorr_
Number of non-orthogonal correctors.
Definition: kinematicSingleLayer.H:98
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::Uw
virtual const volVectorField & Uw() const
Return the film wall velocity [m/s].
Definition: kinematicSingleLayer.C:964
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
fvcFlux.H
Calculate the face-flux of the given field.
Foam::type
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:590
Foam::regionModels::singleLayerRegion::nHat
virtual const volVectorField & nHat() const
Return the patch normal vectors.
Definition: singleLayerRegion.C:210
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::rhoPrimary_
volScalarField rhoPrimary_
Density [kg/m3].
Definition: kinematicSingleLayer.H:199
Foam::constant::electromagnetic::e
const dimensionedScalar e
Elementary charge.
Definition: createFields.H:11
fvcVolumeIntegrate.H
Volume integrate volField creating a volField.
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::updateSurfaceVelocities
virtual void updateSurfaceVelocities()
Update film surface velocities.
Definition: kinematicSingleLayer.C:275
Foam::tmp::New
static tmp< T > New(Args &&... args)
Construct tmp of T with forwarding arguments.
Foam::fvMatrix
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvPatchField.H:68
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::Cp
virtual const volScalarField & Cp() const
Return the film specific heat capacity [J/kg/K].
Definition: kinematicSingleLayer.C:1024
Foam::fvMatrix::solve
SolverPerformance< Type > solve(const dictionary &)
Solve returning the solution statistics.
Definition: fvMatrixSolve.C:319
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::transfer_
transferModelList transfer_
Transfer with the continuous phase.
Definition: kinematicSingleLayer.H:217
Foam::regionModels::regionModel::intCoupledPatchIDs_
labelList intCoupledPatchIDs_
List of patch IDs internally coupled with the primary region.
Definition: regionModel.H:114
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::Sh
virtual tmp< volScalarField::Internal > Sh() const
Return enthalpy source - Eulerian phase only.
Definition: kinematicSingleLayer.C:1126
Foam::name
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
fvcReconstruct.H
Reconstruct volField from a face flux field.
Foam::UIndirectList
A List with indirect addressing.
Definition: faMatrix.H:60
UEqn
fvVectorMatrix & UEqn
Definition: UEqn.H:13
Foam::dictionary::add
entry * add(entry *entryPtr, bool mergeEntry=false)
Add a new entry.
Definition: dictionary.C:640
Foam::GeometricField< scalar, fvPatchField, volMesh >::null
static const GeometricField< scalar, fvPatchField, volMesh > & null()
Return a null geometric field.
Definition: GeometricFieldI.H:32
Foam::gMin
Type gMin(const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:593
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::Tw
virtual const volScalarField & Tw() const
Return the film wall temperature [K].
Definition: kinematicSingleLayer.C:1006
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::postEvolveRegion
virtual void postEvolveRegion()
Post-evolve film hook.
Definition: kinematicSingleLayer.C:912
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::rho
virtual const volScalarField & rho() const
Return the film density [kg/m3].
Definition: kinematicSingleLayer.C:982
Foam::dimVolume
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:60
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::deltaRho
virtual const volScalarField & deltaRho() const
Return the film thickness*density (helper field) [kg/m3].
Definition: kinematicSingleLayer.C:970
rhof
surfaceScalarField rhof(fvc::interpolate(rho, "div(phi,rho)"))
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::updateSubmodels
virtual void updateSubmodels()
Update the film sub-models.
Definition: kinematicSingleLayer.C:210
Foam::GeometricField< scalar, fvPatchField, volMesh >
Foam::IOobject::NO_READ
Definition: IOobject.H:188
Foam::fvm::div
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmDiv.C:48
Foam::regionModels::surfaceFilmModels::filmTurbulenceModel::New
static autoPtr< filmTurbulenceModel > New(surfaceFilmRegionModel &film, const dictionary &dict)
Return a reference to the selected injection model.
Definition: filmTurbulenceModelNew.C:43
Foam::fvMatrix::flux
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > flux() const
Return the face-flux field from the matrix.
Definition: fvMatrix.C:1534
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::rhoSp_
volScalarField rhoSp_
Mass [kg/m2/s].
Definition: kinematicSingleLayer.H:173
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::solveThickness
virtual void solveThickness(const volScalarField &pu, const volScalarField &pp, fvVectorMatrix &UEqn)
Solve coupled velocity-thickness equations.
Definition: kinematicSingleLayer.C:348
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::cloudMassTrans_
volScalarField cloudMassTrans_
Film mass available for transfer to cloud.
Definition: kinematicSingleLayer.H:154
Foam::gMax
Type gMax(const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:592
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::nOuterCorr_
label nOuterCorr_
Number of outer correctors.
Definition: kinematicSingleLayer.H:92
Foam::dictionary::readIfPresent
bool readIfPresent(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX) const
Definition: dictionaryTemplates.C:405
Foam::GeometricField::boundaryField
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Definition: GeometricFieldI.H:62
Foam::dimless
const dimensionSet dimless
Dimensionless.
Definition: dimensionSets.C:189
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::CourantNumber
virtual scalar CourantNumber() const
Courant number evaluation.
Definition: kinematicSingleLayer.C:921
Foam::IOobject::MUST_READ
Definition: IOobject.H:185
Foam::pos
dimensionedScalar pos(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:177