incompressibleDenseParticleFluid.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) 2023-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 
27 #include "surfaceInterpolate.H"
28 #include "fvMatrix.H"
29 #include "fvcFlux.H"
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37 namespace solvers
38 {
41  (
42  solver,
44  fvMesh
45  );
46 }
47 }
48 
49 
50 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
51 
53 {
54  fluidSolver::correctCoNum(phic);
55 }
56 
57 
59 {
61 }
62 
63 
64 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
65 
68 (
69  fvMesh& mesh
70 )
71 :
72  fluidSolver(mesh),
73 
74  continuousPhaseName
75  (
77  (
78  IOobject
79  (
80  "physicalProperties",
81  runTime.constant(),
82  mesh,
83  IOobject::MUST_READ
84  )
85  ).lookup("continuousPhaseName")
86  ),
87 
88  p_
89  (
90  IOobject
91  (
92  "p",
93  runTime.name(),
94  mesh,
95  IOobject::MUST_READ,
96  IOobject::AUTO_WRITE
97  ),
98  mesh
99  ),
100 
102 
103  g
104  (
105  IOobject
106  (
107  "g",
108  runTime.constant(),
109  mesh,
110  IOobject::MUST_READ,
111  IOobject::NO_WRITE
112  )
113  ),
114 
115  Uc_
116  (
117  IOobject
118  (
119  IOobject::groupName("U", continuousPhaseName),
120  runTime.name(),
121  mesh,
122  IOobject::MUST_READ,
123  IOobject::AUTO_WRITE
124  ),
125  mesh
126  ),
127 
128  phic_
129  (
130  IOobject
131  (
132  IOobject::groupName("phi", continuousPhaseName),
133  runTime.name(),
134  mesh,
135  IOobject::READ_IF_PRESENT,
136  IOobject::AUTO_WRITE
137  ),
138  linearInterpolate(Uc_) & mesh.Sf()
139  ),
140 
141  viscosity(viscosityModel::New(mesh)),
142 
143  rhoc
144  (
145  IOobject
146  (
147  IOobject::groupName("rho", continuousPhaseName),
148  runTime.name(),
149  mesh
150  ),
151  mesh,
153  (
154  IOobject::groupName("rho", continuousPhaseName),
155  dimDensity,
156  viscosity->lookup
157  (
158  IOobject::groupName("rho", continuousPhaseName)
159  )
160  )
161  ),
162 
163  muc
164  (
165  IOobject
166  (
167  IOobject::groupName("mu", continuousPhaseName),
168  runTime.name(),
169  mesh
170  ),
171  rhoc*viscosity->nu()
172  ),
173 
174  alphac_
175  (
176  IOobject
177  (
178  IOobject::groupName("alpha", continuousPhaseName),
179  runTime.name(),
180  mesh,
181  IOobject::READ_IF_PRESENT,
182  IOobject::AUTO_WRITE
183  ),
184  mesh,
186  ),
187 
188  alphacMin
189  (
190  1 - mesh.solution().solverDict(alphac_.name()).lookup<scalar>("max")
191  ),
192 
193  alphacf("alphacf", fvc::interpolate(alphac_)),
194 
195  alphaPhic
196  (
197  IOobject::groupName("alphaPhi", continuousPhaseName),
198  alphacf*phic_
199  ),
200 
201  momentumTransport
202  (
203  phaseIncompressible::momentumTransportModel::New
204  (
205  alphac_,
206  Uc_,
207  alphaPhic,
208  phic_,
209  viscosity
210  )
211  ),
212 
213  clouds
214  (
215  parcelClouds::New(mesh, rhoc, Uc_, muc, g)
216  ),
217 
218  p(p_),
219  Uc(Uc_),
220  phic(phic_),
221  alphac(alphac_)
222 {
224 
225  momentumTransport->validate();
226 
227  // Update alphac from the particle locations
228  alphac_ = max(1 - clouds.alpha(), alphacMin);
232 
233  correctCoNum();
234 }
235 
236 
237 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
238 
241 {}
242 
243 
244 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
245 
247 {
248  if (mesh.dynamic() && !Ucf.valid())
249  {
250  Info<< "Constructing face momentum Ucf" << endl;
251 
252  // Ensure the U BCs are up-to-date before constructing Ucf
253  Uc_.correctBoundaryConditions();
254 
255  Ucf = new surfaceVectorField
256  (
257  IOobject
258  (
259  "Ucf",
260  runTime.name(),
261  mesh,
264  ),
265  fvc::interpolate(Uc)
266  );
267  }
268 
269  // Store the particle positions
270  if (mesh.topoChanging() || mesh.distributing())
271  {
272  clouds.storeGlobalPositions();
273  }
274 
276 
277  correctCoNum();
278 
279  // Update the mesh for topology change, mesh to mesh mapping
280  mesh_.update();
281 }
282 
283 
285 {
286  if (pimple.firstIter())
287  {
288  clouds.evolve();
289 
290  // Update continuous phase volume fraction field
291  alphac_ = max(1 - clouds.alpha(), alphacMin);
292  alphac_.correctBoundaryConditions();
293  alphacf = fvc::interpolate(alphac);
294 
295  // ... and continuous phase volumetric flux
296  alphaPhic = alphacf*phic;
297 
298  // Update the continuous phase dynamic viscosity
299  muc = rhoc*viscosity->nu();
300 
301  Fd = new volVectorField
302  (
303  IOobject
304  (
305  "Fd",
306  runTime.name(),
307  mesh
308  ),
309  mesh,
311  zeroGradientFvPatchVectorField::typeName
312  );
313 
314  Dc = new volScalarField
315  (
316  IOobject
317  (
318  "Dc",
319  runTime.name(),
320  mesh
321  ),
322  mesh,
324  zeroGradientFvPatchVectorField::typeName
325  );
326 
327  const fvVectorMatrix cloudSU(clouds.SU(Uc));
328 
329  Fd().primitiveFieldRef() = -cloudSU.source()/mesh.V()/rhoc;
330  Fd().correctBoundaryConditions();
331 
332  Dc().primitiveFieldRef() = -cloudSU.diag()/mesh.V()/rhoc;
333  Dc().correctBoundaryConditions();
334 
335  Dcf = fvc::interpolate(Dc()).ptr();
336 
337  phid =
338  (
339  fvc::flux(Fd())
340  /(Dcf() + dimensionedScalar(Dc().dimensions(), small))
341  ).ptr();
342  }
343 
344  if (pimple.predictTransport())
345  {
346  momentumTransport->predict();
347  }
348 }
349 
350 
352 {}
353 
354 
356 {
357  while (pimple.correct())
358  {
359  correctPressure();
360  }
361 
362  tUcEqn.clear();
363 }
364 
365 
367 {
368  if (pimple.correctTransport())
369  {
370  viscosity->correct();
371  momentumTransport->correct();
372  }
373 }
374 
375 
377 {
378  Fd.clear();
379  Dc.clear();
380  Dcf.clear();
381  phid.clear();
382 }
383 
384 
385 // ************************************************************************* //
Macros for easy insertion into run-time selection tables.
void correctBoundaryConditions()
Correct boundary field.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:57
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
const word & name() const
Return name.
Definition: IOobject.H:310
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:118
Field< Type > & source()
Definition: fvMatrix.H:307
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
const fvSchemes & schemes() const
Return the fvSchemes.
Definition: fvMesh.C:1751
virtual void preUpdateMesh()
Prepare for mesh update.
Definition: fvModels.C:260
void setFluxRequired(const word &name) const
Definition: fvSchemes.C:503
scalarField & diag()
Definition: lduMatrix.C:186
Abstract base class for turbulence models (RAS, LES and laminar).
const tmp< volScalarField > alpha() const
Return the particle volume fraction field.
List of parcel clouds, with the same interface as an individual parcel cloud. Is a mesh object,...
Definition: parcelClouds.H:65
Provides controls for the pressure reference in closed-volume simulations.
Selector class for relaxation factors, solver type and solution.
Definition: solution.H:51
Abstract base class for run-time selectable region solvers.
Definition: solver.H:56
const fvMesh & mesh
Region mesh.
Definition: solver.H:101
Base solver module for fluid solvers.
Definition: fluidSolver.H:62
void continuityErrors(const surfaceScalarField &phi)
Calculate and print the continuity errors.
Definition: fluidSolver.C:146
Solver module for transient flow of incompressible isothermal fluids coupled with particle clouds inc...
virtual void thermophysicalPredictor()
Construct and solve the energy equation,.
autoPtr< phaseIncompressible::momentumTransportModel > momentumTransport
Pointer to the momentum transport model.
surfaceScalarField alphacf
Interpolated continuous phase-fraction.
virtual void prePredictor()
Called at the start of the PIMPLE loop.
const volScalarField & alphac
Reference continuous phase-fraction.
virtual void postSolve()
Called after the PIMPLE loop at the end of the time-step.
surfaceScalarField alphaPhic
Continuous phase volumetric-flux field.
virtual void pressureCorrector()
Construct and solve the pressure equation in the PISO loop.
virtual void postCorrector()
Correct the momentum and thermophysical transport modelling.
void continuityErrors()
Calculate and print the continuity errors.
void correctCoNum()
Correct the cached Courant numbers.
incompressibleDenseParticleFluid(fvMesh &mesh)
Construct from region mesh.
virtual void preSolve()
Called at the start of the time-step, before the PIMPLE loop.
const surfaceScalarField & phic
Reference to the continuous phase volumetric-flux field.
const volScalarField & p
Reference to the pressure field.
An abstract base class for Newtonian viscosity models.
Abstract base class for all fluid physical properties.
Definition: viscosity.H:50
virtual tmp< volScalarField > nu() const =0
Return the laminar viscosity.
Foam::fvModels & fvModels(Foam::fvModels::New(mesh))
pimpleControl pimple(mesh)
Calculate the face-flux of the given field.
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.
addToRunTimeSelectionTable(solver, compressibleMultiphaseVoF, fvMesh)
defineTypeNameAndDebug(compressibleMultiphaseVoF, 0)
Namespace for OpenFOAM.
static const zero Zero
Definition: zero.H:97
VolField< vector > volVectorField
Definition: volFieldsFwd.H:65
const HashTable< dimensionSet > & dimensions()
Get the table of dimension sets.
Definition: dimensionSets.C:96
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
const dimensionSet dimless
messageStream Info
const dimensionSet dimAcceleration
tmp< SurfaceField< Type > > linearInterpolate(const VolField< Type > &vf)
Definition: linear.H:108
const dimensionSet dimTime
const dimensionSet dimDensity
VolField< scalar > volScalarField
Definition: volFieldsFwd.H:64
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
SurfaceField< vector > surfaceVectorField
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dictionary dict
volScalarField & p