scalarTransport.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) 2012-2024 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 "scalarTransport.H"
27 #include "surfaceFields.H"
28 #include "fvmDdt.H"
29 #include "fvcDdt.H"
30 #include "fvmDiv.H"
31 #include "fvmLaplacian.H"
32 #include "fvmSup.H"
33 #include "fvcFlux.H"
34 #include "fvModels.H"
35 #include "fvConstraints.H"
38 
39 #include "CMULES.H"
40 #include "EulerDdtScheme.H"
41 #include "localEulerDdtScheme.H"
42 #include "CrankNicolsonDdtScheme.H"
43 #include "subCycle.H"
44 #include "interfaceCompression.H"
45 
47 
48 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
49 
50 namespace Foam
51 {
52 namespace functionObjects
53 {
55 
57  (
61  );
62 }
63 }
64 
65 
66 template<>
67 const char* Foam::NamedEnum
68 <
70  3
71 >::names[] =
72 {
73  "none",
74  "constant",
75  "viscosity"
76 };
77 
78 const Foam::NamedEnum
79 <
81  3
82 > Foam::functionObjects::scalarTransport::diffusivityTypeNames_;
83 
84 
85 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
86 
88 Foam::functionObjects::scalarTransport::D() const
89 {
90  const word Dname("D" + fieldName_);
91 
92  if (diffusivity_ == diffusivityType::constant)
93  {
94  return volScalarField::New
95  (
96  Dname,
97  mesh_,
99  );
100  }
101  else
102  {
105 
106  return alphal_*turbulence.nu() + alphat_*turbulence.nut();
107  }
108 }
109 
110 
111 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
112 
114 (
115  const word& name,
116  const Time& runTime,
117  const dictionary& dict
118 )
119 :
120  fvMeshFunctionObject(name, runTime, dict),
121  fieldName_(dict.lookupOrDefault<word>("field", "s")),
122  diffusivity_(diffusivityType::none),
123  D_(0),
124  nCorr_(0),
125  s_
126  (
127  IOobject
128  (
129  fieldName_,
130  time_.name(),
131  mesh_,
132  IOobject::MUST_READ,
133  IOobject::AUTO_WRITE
134  ),
135  mesh_
136  ),
137  MULES_(false),
138  deltaN_
139  (
140  "deltaN",
141  1e-8/pow(average(mesh_.V()), 1.0/3.0)
142  ),
143  sRestart_(false)
144 {
145  read(dict);
146 
147  if (mesh_.solution().solversDict().found(fieldName_))
148  {
149  const dictionary& controls = mesh_.solution().solverDict(fieldName_);
150 
151  if (controls.found("nSubCycles"))
152  {
153  MULES_ = true;
154 
156  (
157  IOobject::groupName("sPhi", s_.group()),
158  runTime.name(),
159  mesh_,
162  );
163 
164  sRestart_ = sPhiHeader.headerOk();
165 
166  if (sRestart_)
167  {
168  Info << "Restarting s" << endl;
169  }
170 
171  const surfaceScalarField& phi =
173 
174  // Scalar volumetric flux
175  tsPhi_ = new surfaceScalarField
176  (
177  sPhiHeader,
178  phi*fvc::interpolate(s_)
179  );
180 
181  if (controls.lookupOrDefault<Switch>("MULESCorr", false))
182  {
183  mesh_.schemes().setFluxRequired(fieldName_);
184  }
185  }
186  }
187 }
188 
189 
190 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
191 
193 {}
194 
195 
196 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
197 
199 {
201 
202  phiName_ = dict.lookupOrDefault<word>("phi", "phi");
203  rhoName_ = dict.lookupOrDefault<word>("rho", "rho");
204  schemesField_ = dict.lookupOrDefault<word>("schemesField", fieldName_);
205 
206  diffusivity_ = diffusivityTypeNames_.read(dict.lookup("diffusivity"));
207 
208  switch(diffusivity_)
209  {
210  case diffusivityType::none:
211  break;
212 
213  case diffusivityType::constant:
214  dict.lookup("D") >> D_;
215  break;
216 
218  dict.lookup("alphal") >> alphal_;
219  dict.lookup("alphat") >> alphat_;
220  break;
221  }
222 
223  dict.readIfPresent("nCorr", nCorr_);
224 
225  return true;
226 }
227 
228 
230 {
231  return wordList{phiName_};
232 }
233 
234 
236 {
237  Info<< type() << " write:" << endl;
238 
239  const surfaceScalarField& phi =
240  mesh_.lookupObject<surfaceScalarField>(phiName_);
241 
242  const word divScheme("div(phi," + schemesField_ + ")");
243 
244  // Set under-relaxation coeff
245  scalar relaxCoeff = 0.0;
246  if (mesh_.solution().relaxEquation(schemesField_))
247  {
248  relaxCoeff = mesh_.solution().equationRelaxationFactor(schemesField_);
249  }
250 
253  (
255  );
256 
257  if (phi.dimensions() == dimVolume/dimTime)
258  {
259  if (MULES_)
260  {
261  subCycleMULES();
262 
264  }
265  else
266  {
267  for (int i=0; i<=nCorr_; i++)
268  {
269  fvScalarMatrix sEqn
270  (
271  fvm::ddt(s_)
272  + fvm::div(phi, s_, divScheme)
273  ==
274  fvModels.source(s_)
275  );
276 
277  if (diffusivity_ != diffusivityType::none)
278  {
279  sEqn -= fvm::laplacian(D(), s_);
280  }
281 
282  sEqn.relax(relaxCoeff);
283 
284  fvConstraints.constrain(sEqn);
285 
286  sEqn.solve(schemesField_);
287 
289  }
290  }
291  }
292  else if (phi.dimensions() == dimMass/dimTime)
293  {
294  const volScalarField& rho =
295  mesh_.lookupObject<volScalarField>(rhoName_);
296 
297  for (int i=0; i<=nCorr_; i++)
298  {
299  fvScalarMatrix sEqn
300  (
301  fvm::ddt(rho, s_)
302  + fvm::div(phi, s_, divScheme)
303  ==
304  fvModels.source(rho, s_)
305  );
306 
307  if (diffusivity_ != diffusivityType::none)
308  {
309  sEqn -= fvm::laplacian(rho*D(), s_);
310  }
311 
312  sEqn.relax(relaxCoeff);
313 
314  fvConstraints.constrain(sEqn);
315 
316  sEqn.solve(schemesField_);
317 
319  }
320  }
321  else
322  {
324  << "Incompatible dimensions for phi: " << phi.dimensions() << nl
325  << "Dimensions should be " << dimMass/dimTime << " or "
327  }
328 
329  Info<< endl;
330 
331  return true;
332 }
333 
334 
335 void Foam::functionObjects::scalarTransport::subCycleMULES()
336 {
337  const dictionary& controls = mesh_.solution().solverDict(fieldName_);
338  const label nSubCycles(controls.lookup<label>("nSubCycles"));
339  const bool LTS = fv::localEulerDdt::enabled(mesh_);
340 
341  if (nSubCycles > 1)
342  {
343  tmp<volScalarField> trSubDeltaT;
344 
345  if (LTS)
346  {
347  trSubDeltaT =
348  fv::localEulerDdt::localRSubDeltaT(mesh_, nSubCycles);
349  }
350 
351  for
352  (
353  subCycle<volScalarField> sSubCycle(s_, nSubCycles);
354  !(++sSubCycle).end();
355  )
356  {
357  solveMULES();
358  }
359  }
360  else
361  {
362  solveMULES();
363  }
364 
365 
366  // Apply the diffusivity term separately to allow implicit solution
367  // and boundedness of the explicit advection
368  if (diffusivity_ != diffusivityType::none)
369  {
370  fvScalarMatrix sEqn
371  (
372  fvm::ddt(s_) - fvc::ddt(s_)
373  - fvm::laplacian(D(), s_)
374  );
375 
376  sEqn.solve(controls.subDict("diffusivity"));
377 
378  Info<< fieldName_ << " volume fraction = "
379  << s_.weightedAverage(mesh_.V()).value()
380  << " Min(" << fieldName_ << ") = " << min(s_).value()
381  << " Max(" << fieldName_ << ") = " << max(s_).value()
382  << endl;
383  }
384 }
385 
386 
387 void Foam::functionObjects::scalarTransport::solveMULES()
388 {
389  const dictionary& controls = mesh_.solution().solverDict(fieldName_);
390  const label nCorr(controls.lookup<label>("nCorr"));
391  const label nSubCycles(controls.lookup<label>("nSubCycles"));
392  const bool MULESCorr(controls.lookupOrDefault<Switch>("MULESCorr", false));
393 
394  // Apply the compression correction from the previous iteration
395  // Improves efficiency for steady-simulations but can only be applied
396  // once the s field is reasonably steady, i.e. fully developed
397  const bool applyPrevCorr
398  (
399  controls.lookupOrDefault<Switch>("applyPrevCorr", false)
400  );
401 
402  const bool LTS = fv::localEulerDdt::enabled(mesh_);
403 
404  const word divScheme("div(phi," + schemesField_ + ")");
405 
406  const surfaceScalarField& phi =
407  mesh_.lookupObject<surfaceScalarField>(phiName_);
408 
409  surfaceScalarField& sPhi_ = tsPhi_.ref();
410 
411  const word sScheme(mesh_.schemes().div(divScheme)[1].wordToken());
412 
413  // If a compressive convection scheme is used
414  // the interface normal must be cached
415  tmp<surfaceScalarField> nHatf;
416 
417  if (compressionSchemes.found(sScheme))
418  {
419  const surfaceVectorField gradsf(fvc::interpolate(fvc::grad(s_)));
420 
421  nHatf = new surfaceScalarField
422  (
423  IOobject
424  (
425  "nHatf",
426  s_.time().name(),
427  mesh_
428  ),
429  gradsf/(mag(gradsf) + deltaN_) & mesh_.Sf()
430  );
431  }
432 
433  // Set the off-centering coefficient according to ddt scheme
434  scalar ocCoeff = 0;
435  {
436  tmp<fv::ddtScheme<scalar>> tddtS
437  (
439  (
440  mesh_,
441  mesh_.schemes().ddt("ddt(s)")
442  )
443  );
444  const fv::ddtScheme<scalar>& ddtS = tddtS();
445 
446  if
447  (
448  isType<fv::EulerDdtScheme<scalar>>(ddtS)
449  || isType<fv::localEulerDdtScheme<scalar>>(ddtS)
450  )
451  {
452  ocCoeff = 0;
453  }
454  else if (isType<fv::CrankNicolsonDdtScheme<scalar>>(ddtS))
455  {
456  if (nSubCycles > 1)
457  {
459  << "Sub-cycling is not supported "
460  "with the CrankNicolson ddt scheme"
461  << exit(FatalError);
462  }
463 
464  if
465  (
466  sRestart_
467  || time_.timeIndex() > time_.startTimeIndex() + 1
468  )
469  {
470  ocCoeff =
471  refCast<const fv::CrankNicolsonDdtScheme<scalar>>(ddtS)
472  .ocCoeff();
473  }
474  }
475  else
476  {
478  << "Only Euler and CrankNicolson ddt schemes are supported"
479  << exit(FatalError);
480  }
481  }
482 
483  // Set the time blending factor, 1 for Euler
484  scalar cnCoeff = 1.0/(1.0 + ocCoeff);
485 
486  tmp<surfaceScalarField> phiCN(phi);
487 
488  // Calculate the Crank-Nicolson off-centred volumetric flux
489  if (ocCoeff > 0)
490  {
492  (
493  "phiCN",
494  cnCoeff*phi + (1.0 - cnCoeff)*phi.oldTime()
495  );
496  }
497 
498  if (MULESCorr)
499  {
500  fvScalarMatrix sEqn
501  (
502  (
503  LTS
504  ? fv::localEulerDdtScheme<scalar>(mesh_).fvmDdt(s_)
505  : fv::EulerDdtScheme<scalar>(mesh_).fvmDdt(s_)
506  )
507  + fv::gaussConvectionScheme<scalar>
508  (
509  mesh_,
510  phiCN,
511  upwind<scalar>(mesh_, phiCN)
512  ).fvmDiv(phiCN, s_)
513  );
514 
515  sEqn.solve();
516 
517  Info<< fieldName_ << " volume fraction = "
518  << s_.weightedAverage(mesh_.Vsc()).value()
519  << " Min(" << fieldName_ << ") = " << min(s_).value()
520  << " Max(" << fieldName_ << ") = " << max(s_).value()
521  << endl;
522 
523  tmp<surfaceScalarField> tsPhiUD(sEqn.flux());
524  sPhi_ = tsPhiUD();
525 
526  if (applyPrevCorr && tsPhiCorr0_.valid())
527  {
528  Info<< "Applying the previous iteration compression flux" << endl;
530  (
531  geometricOneField(),
532  s_,
533  sPhi_,
534  tsPhiCorr0_.ref(),
535  oneField(),
536  zeroField()
537  );
538 
539  sPhi_ += tsPhiCorr0_();
540  }
541 
542  // Cache the upwind-flux
543  tsPhiCorr0_ = tsPhiUD;
544  }
545 
546  for (int sCorr=0; sCorr<nCorr; sCorr++)
547  {
548  // Split operator
549  tmp<surfaceScalarField> tsPhiUn
550  (
551  fvc::flux
552  (
553  phiCN(),
554  (cnCoeff*s_ + (1.0 - cnCoeff)*s_.oldTime())(),
555  mesh_.schemes().div(divScheme)
556  )
557  );
558 
559  if (MULESCorr)
560  {
561  tmp<surfaceScalarField> tsPhiCorr(tsPhiUn() - sPhi_);
562  volScalarField s0("s0", s_);
563 
565  (
566  geometricOneField(),
567  s_,
568  tsPhiUn(),
569  tsPhiCorr.ref(),
570  oneField(),
571  zeroField()
572  );
573 
574  // Under-relax the correction for all but the 1st corrector
575  if (sCorr == 0)
576  {
577  sPhi_ += tsPhiCorr();
578  }
579  else
580  {
581  s_ = 0.5*s_ + 0.5*s0;
582  sPhi_ += 0.5*tsPhiCorr();
583  }
584  }
585  else
586  {
587  sPhi_ = tsPhiUn;
588 
590  (
591  geometricOneField(),
592  s_,
593  phiCN,
594  sPhi_,
595  oneField(),
596  zeroField()
597  );
598  }
599  }
600 
601  if (applyPrevCorr && MULESCorr)
602  {
603  tsPhiCorr0_ = sPhi_ - tsPhiCorr0_;
604  tsPhiCorr0_.ref().rename("sPhiCorr0");
605  }
606  else
607  {
608  tsPhiCorr0_.clear();
609  }
610 
611  if
612  (
613  word(mesh_.schemes().ddt("ddt(s)"))
614  == fv::CrankNicolsonDdtScheme<scalar>::typeName
615  )
616  {
617  if (ocCoeff > 0)
618  {
619  // Calculate the end-of-time-step s flux
620  sPhi_ = (sPhi_ - (1.0 - cnCoeff)*sPhi_.oldTime())/cnCoeff;
621  }
622  }
623 
624  Info<< fieldName_ << "volume fraction = "
625  << s_.weightedAverage(mesh_.Vsc()).value()
626  << " Min(" << fieldName_ << ") = " << min(s_).value()
627  << " Max(" << fieldName_ << ") = " << max(s_).value()
628  << endl;
629 }
630 
631 
633 {
634  return true;
635 }
636 
637 
638 // ************************************************************************* //
static const Foam::dimensionedScalar D("D", Foam::dimTemperature, 257.14)
CMULES: Multidimensional universal limiter for explicit corrected implicit solution.
Macros for easy insertion into run-time selection tables.
const bool MULESCorr(alphaControls.lookupOrDefault< Switch >("MULESCorr", false))
const dimensionSet & dimensions() const
Return dimensions.
Generic GeometricField class.
static tmp< GeometricField< Type, PatchField, GeoMesh > > New(const word &name, const Internal &, const PtrList< PatchField< Type >> &, const HashPtrTable< Source > &=HashPtrTable< Source >())
Return a temporary field constructed from name,.
bool found(const Key &) const
Return true if hashedEntry is found in table.
Definition: HashTable.C:138
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
static word group(const word &name)
Return group (extension part of name)
Definition: IOobject.C:134
static word groupName(Name name, const word &group)
Initialise the NamedEnum HashTable from the static list of names.
Definition: NamedEnum.H:54
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:61
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:76
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T, if not found return the given default.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:710
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:843
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:509
static autoPtr< dictionary > New(Istream &)
Construct top-level dictionary on freestore from Istream.
Definition: dictionaryIO.C:100
const word & name() const
Return const reference to name.
Abstract base-class for Time/database functionObjects.
Specialisation of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
const fvMesh & mesh_
Reference to the fvMesh.
virtual bool read(const dictionary &)
Read optional controls.
Evolves a passive scalar transport equation.
virtual wordList fields() const
Return the list of fields required.
diffusivityType
Enumeration defining the type of the diffusivity.
virtual bool execute()
Calculate the scalarTransport.
virtual bool read(const dictionary &)
Read the scalarTransport data.
scalarTransport(const word &name, const Time &runTime, const dictionary &dict)
Construct from Time and dictionary.
Finite volume constraints.
Definition: fvConstraints.H:67
bool constrain(fvMatrix< Type > &eqn) const
Apply constraints to an equation.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:118
void relax(const scalar alpha)
Relax matrix (for steady-state solution).
Definition: fvMatrix.C:603
SolverPerformance< Type > solve(const dictionary &)
Solve segregated or coupled returning the solution statistics.
Definition: fvMatrixSolve.C:58
const fvSchemes & schemes() const
Return the fvSchemes.
Definition: fvMesh.C:1751
const fvSolution & solution() const
Return the fvSolution.
Definition: fvMesh.C:1762
Finite volume models.
Definition: fvModels.H:65
tmp< fvMatrix< Type > > source(const VolField< Type > &field) const
Return source for an equation.
void setFluxRequired(const word &name) const
Definition: fvSchemes.C:503
static tmp< ddtScheme< Type > > New(const fvMesh &mesh, Istream &schemeData)
Return a pointer to a new ddtScheme created on freestore.
Definition: ddtScheme.C:45
static tmp< volScalarField > localRSubDeltaT(const fvMesh &mesh, const label nAlphaSubCycles)
Calculate and return the reciprocal of the local sub-cycling.
Definition: localEulerDdt.C:70
static bool enabled(const fvMesh &mesh)
Return true if LTS is enabled.
Definition: localEulerDdt.C:37
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type and name.
const Type & lookupType(const word &group=word::null) const
Lookup and return the object of the given Type.
const dictionary & solversDict() const
Return the solver controls dictionary.
Definition: solution.C:342
const dictionary & solverDict(const word &name) const
Return the solver controls dictionary for the given field.
Definition: solution.C:348
Perform a subCycleTime on a field or list of fields.
Definition: subCycle.H:235
A class for managing temporary objects.
Definition: tmp.H:55
Templated form of IOobject providing type information for file reading and header type checking.
Definition: IOobject.H:531
bool headerOk()
Read header (uses typeGlobalFile to find file) and check.
A class for handling words, derived from string.
Definition: word.H:62
Foam::fvConstraints & fvConstraints(Foam::fvConstraints::New(mesh))
Foam::fvModels & fvModels(Foam::fvModels::New(mesh))
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
Calculate the first temporal derivative.
Calculate the face-flux of the given field.
Calculate the matrix for the first temporal derivative.
Calculate the matrix for the divergence of the given field and flux.
Calculate the matrix for the laplacian of the field.
Calculate the matrix for implicit and explicit sources.
void explicitSolve(const RdeltaTType &rDeltaT, const RhoType &rho, volScalarField &psi, const surfaceScalarField &phiPsi, const SpType &Sp, const SuType &Su)
void correct(const RdeltaTType &rDeltaT, const RhoType &rho, volScalarField &psi, const surfaceScalarField &phiCorr, const SpType &Sp, const SuType &Su)
compressibleMomentumTransportModel momentumTransportModel
defineTypeNameAndDebug(adjustTimeStepToCombustion, 0)
addToRunTimeSelectionTable(functionObject, adjustTimeStepToCombustion, dictionary)
tmp< SurfaceField< typename innerProduct< vector, Type >::type > > flux(const VolField< Type > &vf)
Return the face-flux field obtained from the given volVectorField.
static tmp< SurfaceField< Type > > interpolate(const VolField< Type > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
tmp< VolField< Type > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
Definition: fvcDdt.C:45
tmp< VolField< typename outerProduct< vector, Type >::type > > grad(const SurfaceField< Type > &ssf)
Definition: fvcGrad.C:46
tmp< fvMatrix< Type > > laplacian(const VolField< Type > &vf, const word &name)
Definition: fvmLaplacian.C:47
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const VolField< Type > &vf, const word &name)
Definition: fvmDiv.C:48
tmp< fvMatrix< Type > > ddt(const VolField< Type > &vf)
Definition: fvmDdt.C:46
Namespace for OpenFOAM.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
const doubleScalar e
Definition: doubleScalar.H:106
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
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
const dimensionSet dimKinematicViscosity
bool isType(const Type &t)
Check the typeid.
Definition: typeInfo.H:158
SurfaceField< scalar > surfaceScalarField
messageStream Info
dimensioned< Type > average(const DimensionedField< Type, GeoMesh > &df)
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
const dimensionSet dimTime
const dimensionSet dimVolume
const wordHashSet compressionSchemes
VolField< scalar > volScalarField
Definition: volFieldsFwd.H:64
dimensioned< scalar > mag(const dimensioned< Type > &)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
const dimensionSet dimMass
error FatalError
SurfaceField< vector > surfaceVectorField
static const char nl
Definition: Ostream.H:266
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
dictionary dict
Info<< "Reading field U\n"<< endl;volVectorField U(IOobject("U", runTime.name(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);Info<< "Creating face flux\n"<< endl;surfaceScalarField phi(IOobject("phi", runTime.name(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), mesh, dimensionedScalar(mesh.Sf().dimensions() *U.dimensions(), 0));autoPtr< viscosityModel > viscosity(viscosityModel::New(mesh))
autoPtr< incompressible::momentumTransportModel > turbulence(incompressible::momentumTransportModel::New(U, phi, viscosity))
Foam::surfaceFields.