applications/solvers/multiphase/interDyMFoam/createFields.H

   1     Info<< "Reading field pd\n" << endl;
   2     volScalarField pd
   3     (
   4         IOobject
   5         (
   6             "pd",
   7             runTime.timeName(),
   8             mesh,
   9             IOobject::MUST_READ,
  10             IOobject::AUTO_WRITE
  11         ),
  12         mesh
  13     );
  14
  15     Info<< "Reading field alpha1\n" << endl;
  16     volScalarField alpha1
  17     (
  18         IOobject
  19         (
  20             "alpha1",
  21             runTime.timeName(),
  22             mesh,
  23             IOobject::MUST_READ,
  24             IOobject::AUTO_WRITE
  25         ),
  26         mesh
  27     );
  28
  29     Info<< "Reading field U\n" << endl;
  30     volVectorField U
  31     (
  32         IOobject
  33         (
  34             "U",
  35             runTime.timeName(),
  36             mesh,
  37             IOobject::MUST_READ,
  38             IOobject::AUTO_WRITE
  39         ),
  40         mesh
  41     );
  42
  43 #   include "createPhi.H"
  44
  45
  46     Info<< "Reading transportProperties\n" << endl;
  47     twoPhaseMixture twoPhaseProperties(U, phi, "alpha1");
  48
  49     const dimensionedScalar& rho1 = twoPhaseProperties.rho1();
  50     const dimensionedScalar& rho2 = twoPhaseProperties.rho2();
  51
  52
  53     // Need to store rho for ddt(rho, U)
  54     volScalarField rho
  55     (
  56         IOobject
  57         (
  58             "rho",
  59             runTime.timeName(),
  60             mesh,
  61             IOobject::READ_IF_PRESENT
  62         ),
  63         alpha1*rho1 + (scalar(1) - alpha1)*rho2,
  64         alpha1.boundaryField().types()
  65     );
  66     rho.oldTime();
  67
  68
  69     // Mass flux
  70     // Initialisation does not matter because rhoPhi is reset after the
  71     // alpha1 solution before it is used in the U equation.
  72     surfaceScalarField rhoPhi
  73     (
  74         IOobject
  75         (
  76             "rho*phi",
  77             runTime.timeName(),
  78             mesh,
  79             IOobject::NO_READ,
  80             IOobject::NO_WRITE
  81         ),
  82         rho1*phi
  83     );
  84
  85
  86     // Construct interface from alpha1 distribution
  87     interfaceProperties interface(alpha1, U, twoPhaseProperties);
  88
  89     // Construct incompressible turbulence model
  90     autoPtr<incompressible::turbulenceModel> turbulence
  91     (
  92         incompressible::turbulenceModel::New(U, phi, twoPhaseProperties)
  93     );
  94
  95     wordList pcorrTypes
  96     (
  97         pd.boundaryField().size(),
  98         zeroGradientFvPatchScalarField::typeName
  99     );
 100
 101     for (label i = 0; i < pd.boundaryField().size(); i++)
 102     {
 103         if (pd.boundaryField()[i].fixesValue())
 104         {
 105             pcorrTypes[i] = fixedValueFvPatchScalarField::typeName;
 106         }
 107     }
 108
 109
 110     volScalarField p
 111     (
 112         IOobject
 113         (
 114             "p",
 115             runTime.timeName(),
 116             mesh,
 117             IOobject::NO_READ,
 118             IOobject::AUTO_WRITE
 119         ),
 120         pd + rho*(g & mesh.C())
 121     );
 122
 123
 124     label pdRefCell = 0;
 125     scalar pdRefValue = 0.0;
 126     setRefCell
 127     (
 128         pd,
 129         mesh.solutionDict().subDict("PIMPLE"),
 130         pdRefCell,
 131         pdRefValue
 132     );
 133
 134     scalar pRefValue = 0.0;
 135
 136     if (pd.needReference())
 137     {
 138         pRefValue = readScalar
 139         (
 140             mesh.solutionDict().subDict("PIMPLE").lookup("pRefValue")
 141         );
 142
 143         p += dimensionedScalar
 144         (
 145             "p",
 146             p.dimensions(),
 147             pRefValue - getRefCellValue(p, pdRefCell)
 148         );
 149     }