electrostaticDepositionFvPatchScalarField.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) 2021 OpenCFD Ltd.
9 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 \*---------------------------------------------------------------------------*/
27 
30 #include "fvPatchFieldMapper.H"
31 #include "volFields.H"
32 #include "surfaceFields.H"
33 
34 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
35 
37 Foam::electrostaticDepositionFvPatchScalarField::eVPatch
38 (
39  const label patchi
40 ) const
41 {
42  const auto& eV =
43  db().lookupObject<volScalarField>(this->internalField().name());
44 
45  const volScalarField::Boundary& bf = eV.boundaryField();
46 
47  const auto& eVpf =
48  refCast<const electrostaticDepositionFvPatchScalarField>(bf[patchi]);
49 
50  return const_cast<electrostaticDepositionFvPatchScalarField&>(eVpf);
51 }
52 
53 
54 void Foam::electrostaticDepositionFvPatchScalarField::setMaster() const
55 {
56  if (master_ != -1)
57  {
58  return;
59  }
60 
61  const auto& eV =
62  db().lookupObject<volScalarField>(this->internalField().name());
63 
64  const volScalarField::Boundary& bf = eV.boundaryField();
65 
66  label master = -1;
67  forAll(bf, patchi)
68  {
69  if (isA<electrostaticDepositionFvPatchScalarField>(bf[patchi]))
70  {
71  electrostaticDepositionFvPatchScalarField& eVpf = eVPatch(patchi);
72 
73  if (master == -1)
74  {
75  master = patchi;
76  }
77 
78  eVpf.master() = master;
79  }
80  }
81 }
82 
83 
84 void Foam::electrostaticDepositionFvPatchScalarField::round
85 (
87  const scalar dcml
88 ) const
89 {
90  for (auto& f : fld)
91  {
92  f = std::round(f*dcml)/dcml;
93  }
94 }
95 
96 
97 void Foam::electrostaticDepositionFvPatchScalarField::writeFilmFields() const
98 {
99  const auto& eV =
100  db().lookupObject<volScalarField>(this->internalField().name());
101 
102  const volScalarField::Boundary& bf = eV.boundaryField();
103 
104  const fvMesh& mesh = eV.mesh();
105 
107  (
108  IOobject
109  (
110  IOobject::scopedName("electrostaticDeposition", "h"),
111  mesh.time().timeName(),
112  mesh,
115  false // do not register
116  ),
117  mesh,
119  );
120 
121  forAll(bf, patchi)
122  {
123  if (isA<electrostaticDepositionFvPatchScalarField>(bf[patchi]))
124  {
125  electrostaticDepositionFvPatchScalarField& eVpf = eVPatch(patchi);
126 
127  auto& hp = h.boundaryFieldRef()[patchi];
128 
129  hp = eVpf.h();
130  }
131  }
132 
133  h.write();
134 }
135 
136 
137 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
138 
141 (
142  const fvPatch& p,
144 )
145 :
146  fixedValueFvPatchScalarField(p, iF),
147  h_(p.size(), 0),
148  qcum_(p.size(), 0),
149  Vfilm_(p.size(), 0),
150  Ceffptr_(nullptr),
151  rptr_(nullptr),
152  jMin_(0),
153  qMin_(0),
154  Rbody_(0),
155  Vi_(0),
156  Vanode_(GREAT),
157  phasesDict_(),
158  phaseNames_(),
159  phases_(),
160  sigmas_(),
161  sigma_(sqr(dimCurrent)*pow3(dimTime)/(dimMass*pow3(dimLength)), scalar(1)),
162  timei_(-1),
163  master_(-1)
164 {}
165 
166 
169 (
170  const fvPatch& p,
172  const dictionary& dict
173 )
174 :
175  fixedValueFvPatchScalarField(p, iF, dict, false),
176  h_("h", dict, p.size()),
177  qcum_
178  (
179  dict.found("qCumulative")
180  ? scalarField("qCumulative", dict, p.size())
181  : scalarField(p.size(), 0)
182  ),
183  Vfilm_
184  (
185  dict.found("Vfilm")
186  ? scalarField("Vfilm", dict, p.size())
187  : scalarField(p.size(), 0)
188  ),
189  Ceffptr_
190  (
191  PatchFunction1<scalar>::New(p.patch(), "CoulombicEfficiency", dict)
192  ),
193  rptr_(PatchFunction1<scalar>::New(p.patch(), "resistivity", dict)),
194  jMin_(dict.getCheckOrDefault<scalar>("jMin", 0, scalarMinMax::ge(0))),
195  qMin_(dict.getCheckOrDefault<scalar>("qMin", 0, scalarMinMax::ge(0))),
196  Rbody_(dict.getCheckOrDefault<scalar>("Rbody", 0, scalarMinMax::ge(0))),
197  Vi_(dict.getOrDefault<scalar>("Vi", 0)),
198  Vanode_(dict.getOrDefault<scalar>("Vanode", GREAT)),
199  phasesDict_(dict.subOrEmptyDict("phases")),
200  phaseNames_(),
201  phases_(),
202  sigmas_(),
203  sigma_
204  (
206  (
208  dict.getCheckOrDefault<scalar>
209  (
210  "sigma",
211  scalar(1),
212  scalarMinMax::ge(SMALL)
213  )
214  )
215  ),
216  timei_(-1),
217  master_(-1)
218 {
219  if (dict.found("value"))
220  {
222  (
223  scalarField("value", dict, p.size())
224  );
225  }
226  else
227  {
228  fvPatchScalarField::operator=(patchInternalField());
229  }
230 
231  // If flow is multiphase
232  if (!phasesDict_.empty())
233  {
234  phaseNames_.setSize(phasesDict_.size());
235  phases_.setSize(phasesDict_.size());
236  sigmas_.setSize(phasesDict_.size());
237 
238  label phasei = 0;
239  for (const entry& dEntry : phasesDict_)
240  {
241  const word& key = dEntry.keyword();
242 
243  if (!dEntry.isDict())
244  {
245  FatalIOErrorInFunction(phasesDict_)
246  << "Entry " << key << " is not a dictionary" << nl
247  << exit(FatalIOError);
248  }
249 
250  const dictionary& subDict = dEntry.dict();
251 
252  phaseNames_[phasei] = key;
253 
254  sigmas_.set
255  (
256  phasei,
258  (
260  subDict.getCheck<scalar>
261  (
262  "sigma",
263  scalarMinMax::ge(SMALL)
264  )
265  )
266  );
267 
268  ++phasei;
269  }
270 
271  forAll(phaseNames_, i)
272  {
273  phases_.set
274  (
275  i,
276  db().getObjectPtr<volScalarField>(phaseNames_[i])
277  );
278  }
279  }
280 }
281 
282 
285 (
287  const fvPatch& p,
289  const fvPatchFieldMapper& mapper
290 )
291 :
292  fixedValueFvPatchScalarField(ptf, p, iF, mapper),
293  h_(ptf.h_, mapper),
294  qcum_(ptf.qcum_, mapper),
295  Vfilm_(ptf.Vfilm_, mapper),
296  Ceffptr_(ptf.Ceffptr_.clone(p.patch())),
297  rptr_(ptf.rptr_.clone(p.patch())),
298  jMin_(ptf.jMin_),
299  qMin_(ptf.qMin_),
300  Rbody_(ptf.Rbody_),
301  Vi_(ptf.Vi_),
302  Vanode_(ptf.Vanode_),
303  phasesDict_(ptf.phasesDict_),
304  phaseNames_(ptf.phaseNames_),
305  phases_(ptf.phases_),
306  sigmas_(),
307  sigma_(ptf.sigma_),
308  timei_(ptf.timei_),
309  master_(-1)
310 {}
311 
312 
315 (
317 )
318 :
319  fixedValueFvPatchScalarField(ptf),
320  h_(ptf.h_),
321  qcum_(ptf.qcum_),
322  Vfilm_(ptf.Vfilm_),
323  Ceffptr_(ptf.Ceffptr_.clone(patch().patch())),
324  rptr_(ptf.rptr_.clone(patch().patch())),
325  jMin_(ptf.jMin_),
326  qMin_(ptf.qMin_),
327  Rbody_(ptf.Rbody_),
328  Vi_(ptf.Vi_),
329  Vanode_(ptf.Vanode_),
330  phasesDict_(ptf.phasesDict_),
331  phaseNames_(ptf.phaseNames_),
332  phases_(ptf.phases_),
333  sigmas_(),
334  sigma_(ptf.sigma_),
335  timei_(ptf.timei_),
336  master_(-1)
337 {}
338 
339 
342 (
345 )
346 :
347  fixedValueFvPatchScalarField(ptf, iF),
348  h_(ptf.h_),
349  qcum_(ptf.qcum_),
350  Vfilm_(ptf.Vfilm_),
351  Ceffptr_(ptf.Ceffptr_.clone(patch().patch())),
352  rptr_(ptf.rptr_.clone(patch().patch())),
353  jMin_(ptf.jMin_),
354  qMin_(ptf.qMin_),
355  Rbody_(ptf.Rbody_),
356  Vi_(ptf.Vi_),
357  Vanode_(ptf.Vanode_),
358  phasesDict_(ptf.phasesDict_),
359  phaseNames_(ptf.phaseNames_),
360  phases_(ptf.phases_),
361  sigmas_(),
362  sigma_(ptf.sigma_),
363  timei_(ptf.timei_),
364  master_(-1)
365 {}
366 
367 
368 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
369 
371 (
372  const fvPatchFieldMapper& m
373 )
374 {
375  fixedValueFvPatchScalarField::autoMap(m);
376 
377  h_.autoMap(m);
378  qcum_.autoMap(m);
379  Vfilm_.autoMap(m);
380 
381  if (Ceffptr_)
382  {
383  Ceffptr_->autoMap(m);
384  }
385 
386  if (rptr_)
387  {
388  rptr_->autoMap(m);
389  }
390 }
391 
392 
394 (
395  const fvPatchScalarField& ptf,
396  const labelList& addr
397 )
398 {
399  fixedValueFvPatchScalarField::rmap(ptf, addr);
400 
401  const auto& tiptf =
402  refCast<const electrostaticDepositionFvPatchScalarField>(ptf);
403 
404  h_.rmap(tiptf.h_, addr);
405  qcum_.rmap(tiptf.qcum_, addr);
406  Vfilm_.rmap(tiptf.Vfilm_, addr);
407 
408  if (Ceffptr_)
409  {
410  Ceffptr_->rmap(tiptf.Ceffptr_(), addr);
411  }
412 
413  if (rptr_)
414  {
415  rptr_->rmap(tiptf.rptr_(), addr);
416  }
417 }
418 
419 
422 {
423  const label patchi = patch().index();
424 
425  if (phases_.size())
426  {
427  tmp<scalarField> tsigma =
428  phases_[0].boundaryField()[patchi]*sigmas_[0].value();
429 
430  for (label i = 1; i < phases_.size(); ++i)
431  {
432  tsigma.ref() +=
433  phases_[i].boundaryField()[patchi]*sigmas_[i].value();
434  }
435 
436  return tsigma;
437  }
438 
439  return tmp<scalarField>::New(patch().size(), sigma_.value());
440 }
441 
442 
444 {
445  if (updated())
446  {
447  return;
448  }
449 
450  if (timei_ == db().time().timeIndex())
451  {
452  return;
453  }
454 
455  const scalar t = db().time().timeOutputValue();
456  const scalar dt = db().time().deltaTValue();
457  const label patchi = patch().index();
458 
459  const auto& eV =
460  db().lookupObject<volScalarField>(this->internalField().name());
461 
462  // Current density on film interface
463  tmp<scalarField> tjnp = -this->sigma()*eV.boundaryField()[patchi].snGrad();
464  scalarField& jnp = tjnp.ref();
465  jnp = max(jnp, scalar(0)); // experimental - do not allow any negative jnp
466  // experimental - avoid micro/nano currents/volts
467  // to reduce snowballing effects of lateral gradients on the patch
468  round(jnp);
469 
470 
471  // Calculate film-thickness finite increments
472  tmp<scalarField> tCoulombicEfficiency = Ceffptr_->value(t);
473  tmp<scalarField> tdh = tCoulombicEfficiency*(jnp - jMin_)*dt;
474  scalarField& dh = tdh.ref();
475 
476  // Do not allow any depletion or abrasion of deposition
477  dh = max(dh, scalar(0));
478 
479  // Do not allow any deposition when accumulative specific
480  // charge is less than minimum accumulative specific charge
481  qcum_ += jnp*dt;
482 
483  forAll(dh, i)
484  {
485  if (qcum_[i] < qMin_)
486  {
487  dh[i] = 0;
488  }
489  }
490 
491  // Add finite increments of film thickness to total film thickness
492  h_ += dh;
493 
494 
495  // Calculate incremental electric potential due to film resistance
496  tmp<scalarField> tresistivity = rptr_->value(t);
497  tmp<scalarField> tRfilm = tresistivity*tdh;
498  tmp<scalarField> tdV = jnp*tRfilm;
499  Vfilm_ += tdV;
500  Vfilm_ = min(Vfilm_, Vanode_);
501 
502 
503  // Calculate electric potential due to body resistance
504  tmp<scalarField> tVbody = tjnp*Rbody_;
505 
506 
507  // Add all electric potential contributions
508  operator==(min(Vi_ + Vfilm_ + tVbody, Vanode_));
509 
510 
511  fixedValueFvPatchScalarField::updateCoeffs();
512 
513  timei_ = db().time().timeIndex();
514 
515  {
516  const scalar hMin = gMin(h_);
517  const scalar hMax = gMax(h_);
518  const scalar hAvg = gAverage(h_);
519 
520  if (Pstream::master())
521  {
522  Info<< " patch: " << patch().name()
523  << ", h: min = " << hMin
524  << ", max = " << hMax
525  << ", average = " << hAvg << nl
526  << endl;
527  }
528  }
529 
530  // Write here to avoid any upset to redistributePar-decompose
531  if (db().time().writeTime())
532  {
533  // Write film thickness fields as patch fields of a volScalarField
534  setMaster();
535 
536  if (patch().index() == master_)
537  {
538  writeFilmFields();
539  }
540  }
541 }
542 
543 
545 {
547 
548  h_.writeEntry("h", os);
549 
550  if (Ceffptr_)
551  {
552  Ceffptr_->writeData(os);
553  }
554 
555  if (rptr_)
556  {
557  rptr_->writeData(os);
558  }
559 
560  if (!phasesDict_.empty())
561  {
562  phasesDict_.writeEntry(phasesDict_.dictName(), os);
563  }
564  else
565  {
566  sigma_.writeEntry("sigma", os);
567  }
568 
569  os.writeEntryIfDifferent<scalar>("jMin", 0, jMin_);
570  os.writeEntryIfDifferent<scalar>("qMin", 0, qMin_);
571  os.writeEntryIfDifferent<scalar>("Rbody", 0, Rbody_);
572  os.writeEntryIfDifferent<scalar>("Vi", 0, Vi_);
573  os.writeEntryIfDifferent<scalar>("Vanode", GREAT, Vanode_);
574  qcum_.writeEntry("qCumulative", os);
575  Vfilm_.writeEntry("Vfilm", os);
576 
577  writeEntry("value", os);
578 }
579 
580 
581 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
582 
583 namespace Foam
584 {
586  (
589  );
590 }
591 
592 // ************************************************************************* //
Foam::fvPatchField< scalar >
Foam::entry
A keyword and a list of tokens is an 'entry'.
Definition: entry.H:67
Foam::IOobject::NO_WRITE
Definition: IOobject.H:195
volFields.H
Foam::fvPatchField< scalar >::write
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:384
Foam::Ostream::writeEntryIfDifferent
Ostream & writeEntryIfDifferent(const word &key, const T &value1, const T &value2)
Write a keyword/value entry only when the two values differ.
Definition: Ostream.H:248
Foam::scalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Definition: primitiveFieldsFwd.H:52
p
volScalarField & p
Definition: createFieldRefs.H:8
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::electrostaticDepositionFvPatchScalarField::autoMap
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
Definition: electrostaticDepositionFvPatchScalarField.C:371
Foam::dimLength
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:52
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::dictionary::found
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search for an entry (const access) with the given keyword.
Definition: dictionaryI.H:87
Foam::gAverage
Type gAverage(const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:604
Foam::Time::timeName
static word timeName(const scalar t, const int precision=precision_)
Definition: Time.C:780
Foam::glTF::key
auto key(const Type &t) -> typename std::enable_if< std::is_enum< Type >::value, typename std::underlying_type< Type >::type >::type
Definition: foamGltfBase.H:108
Foam::UPstream::master
static bool master(const label communicator=worldComm)
Am I the master process.
Definition: UPstream.H:457
Foam::fvPatchField< scalar >::operator
friend Ostream & operator(Ostream &, const fvPatchField< scalar > &)
phasei
label phasei
Definition: pEqn.H:27
Foam::FatalIOError
IOerror FatalIOError
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
surfaceFields.H
Foam::surfaceFields.
fvPatchFieldMapper.H
Foam::min
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::dimCurrent
const dimensionSet dimCurrent(0, 0, 0, 0, 0, 1, 0)
Definition: dimensionSets.H:56
Foam::dimTime
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
Foam::tmp::ref
T & ref() const
Definition: tmpI.H:227
Foam::Field< scalar >
Foam::operator==
tmp< faMatrix< Type > > operator==(const faMatrix< Type > &, const faMatrix< Type > &)
Foam::pow3
dimensionedScalar pow3(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:89
Foam::constant::universal::h
const dimensionedScalar h
Planck constant.
Definition: setRegionSolidFields.H:33
Foam::Info
messageStream Info
Information stream (stdout output on master, null elsewhere)
Foam::fvPatch
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:65
Foam::dimensionedScalar
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Definition: dimensionedScalarFwd.H:42
Foam::volScalarField
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
Foam::dictionary::getCheckOrDefault
T getCheckOrDefault(const word &keyword, const T &deflt, const Predicate &pred, enum keyType::option matchOpt=keyType::REGEX) const
Definition: dictionaryTemplates.C:209
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
fld
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< ' ';}gmvFile<< nl;for(const word &name :lagrangianScalarNames){ IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Definition: gmvOutputLagrangian.H:23
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
os
OBJstream os(runTime.globalPath()/outputName)
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam::PatchFunction1
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
Definition: PatchFunction1.H:60
Foam::dimensioned< scalar >
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::dictionary::subOrEmptyDict
dictionary subOrEmptyDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX, const bool mandatory=false) const
Definition: dictionary.C:540
Foam::electrostaticDepositionFvPatchScalarField::rmap
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
Definition: electrostaticDepositionFvPatchScalarField.C:394
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
Foam::electrostaticDepositionFvPatchScalarField::sigma
tmp< scalarField > sigma() const
Return the isotropic electrical conductivity field of mixture.
Definition: electrostaticDepositionFvPatchScalarField.C:421
Foam::sqr
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Definition: dimensionedSymmTensor.C:51
Foam::nl
constexpr char nl
Definition: Ostream.H:404
f
labelList f(nPoints)
Foam::foamVersion::patch
const std::string patch
OpenFOAM patch number as a std::string.
Foam::MinMax::ge
static MinMax< T > ge(const T &minVal)
A semi-infinite range from minVal to the type max.
Definition: MinMaxI.H:31
Foam::List< label >
electrostaticDepositionFvPatchScalarField.H
Foam::tmp::New
static tmp< T > New(Args &&... args)
Construct tmp of T with forwarding arguments.
Foam::dictionary::getCheck
T getCheck(const word &keyword, const Predicate &pred, enum keyType::option matchOpt=keyType::REGEX) const
Definition: dictionaryTemplates.C:120
Foam::Ostream::writeEntry
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Definition: Ostream.H:236
Foam::electrostaticDepositionFvPatchScalarField
The electrostaticDeposition is a boundary condition to calculate electric potential (V) on a given bo...
Definition: electrostaticDepositionFvPatchScalarField.H:319
Foam::fvMesh::time
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:280
Foam::electrostaticDepositionFvPatchScalarField::updateCoeffs
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Definition: electrostaticDepositionFvPatchScalarField.C:443
FatalIOErrorInFunction
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:473
Foam::fvPatchFieldMapper
Foam::fvPatchFieldMapper.
Definition: fvPatchFieldMapper.H:47
timeIndex
label timeIndex
Definition: getTimeIndex.H:30
Foam::Ostream
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:56
sigma
dimensionedScalar sigma("sigma", dimMass/sqr(dimTime), transportProperties)
Foam::makePatchTypeField
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
Foam::gMin
Type gMin(const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:593
Foam::IOobject::scopedName
static word scopedName(const std::string &scope, const word &name)
Create scope:name or scope_name string.
Definition: IOobjectI.H:47
Foam::dictionary::getOrDefault
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
Definition: dictionaryTemplates.C:148
Foam::GeometricField< scalar, fvPatchField, volMesh >
Foam::IOobject::NO_READ
Definition: IOobject.H:188
Foam::fvPatchField< scalar >::operator=
virtual void operator=(const UList< scalar > &)
Definition: fvPatchField.C:404
Foam::electrostaticDepositionFvPatchScalarField::write
virtual void write(Ostream &) const
Write.
Definition: electrostaticDepositionFvPatchScalarField.C:544
Foam::gMax
Type gMax(const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:592
Foam::GeometricField::boundaryField
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Definition: GeometricFieldI.H:62
Foam::DimensionedField
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Definition: DimensionedField.H:54
Foam::electrostaticDepositionFvPatchScalarField::electrostaticDepositionFvPatchScalarField
electrostaticDepositionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
Definition: electrostaticDepositionFvPatchScalarField.C:141