All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
UEqns.H
Go to the documentation of this file.
1 Info<< "Constructing face momentum equations" << endl;
2 
3 PtrList<fvVectorMatrix> UEqns(phases.size());
4 
5 {
6  fluid.momentumTransfer(); // !!! Update coefficients shouldn't be necessary
7  // This should be done on demand
8 
9  autoPtr<phaseSystem::momentumTransferTable>
10  momentumTransferPtr(fluid.momentumTransferf());
11 
12  phaseSystem::momentumTransferTable&
14 
15  forAll(fluid.movingPhases(), movingPhasei)
16  {
17  phaseModel& phase = fluid.movingPhases()[movingPhasei];
18 
19  const volScalarField& alpha = phase;
20  const volScalarField& rho = phase.rho();
21  volVectorField& U = phase.URef();
22 
23  UEqns.set
24  (
25  phase.index(),
26  new fvVectorMatrix
27  (
28  phase.UfEqn()
29  ==
30  *momentumTransfer[phase.name()]
31  + fvModels.source(alpha, rho, U)
32  )
33  );
34 
35  UEqns[phase.index()].relax();
36  fvConstraints.constrain(UEqns[phase.index()]);
37  U.correctBoundaryConditions();
39  }
40 }
Info<< "Constructing momentum equations"<< endl;PtrList< fvVectorMatrix > UEqns(phases.size())
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:58
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
Foam::fvConstraints & fvConstraints
phaseSystem & fluid
Definition: createFields.H:11
phaseSystem::momentumTransferTable & momentumTransfer(momentumTransferPtr())
forAll(fluid.movingPhases(), movingPhasei)
Definition: UEqns.H:12
bool constrain(fvMatrix< Type > &eqn) const
Apply constraints to an equation.
Foam::fvModels & fvModels
U
Definition: pEqn.H:72
fvMatrix< vector > fvVectorMatrix
Definition: fvMatricesFwd.H:45
phaseSystem::phaseModelList & phases
Definition: createFields.H:12
messageStream Info
autoPtr< phaseSystem::momentumTransferTable > momentumTransferPtr(fluid.momentumTransferf())