applications/solvers/multiphase/twoPhaseEulerFoam/createFields.H

   1     Info<< "Reading transportProperties\n" << endl;
   2
   3     IOdictionary transportProperties
   4     (
   5         IOobject
   6         (
   7             "transportProperties",
   8             runTime.constant(),
   9             mesh,
  10             IOobject::MUST_READ,
  11             IOobject::NO_WRITE
  12         )
  13     );
  14
  15     autoPtr<phaseModel> phasea = phaseModel::New
  16     (
  17         mesh,
  18         transportProperties,
  19         "a"
  20     );
  21
  22     autoPtr<phaseModel> phaseb = phaseModel::New
  23     (
  24         mesh,
  25         transportProperties,
  26         "b"
  27     );
  28
  29     volVectorField& Ua = phasea->U();
  30     surfaceScalarField& phia = phasea->phi();
  31     const dimensionedScalar& rhoa = phasea->rho();
  32     const dimensionedScalar& nua = phasea->nu();
  33
  34     volVectorField& Ub = phaseb->U();
  35     surfaceScalarField& phib = phaseb->phi();
  36     const dimensionedScalar& rhob = phaseb->rho();
  37     const dimensionedScalar& nub = phaseb->nu();
  38
  39     Info<< "Reading field alpha\n" << endl;
  40     volScalarField alpha
  41     (
  42         IOobject
  43         (
  44             "alpha",
  45             runTime.timeName(),
  46             mesh,
  47             IOobject::MUST_READ,
  48             IOobject::AUTO_WRITE
  49         ),
  50         mesh
  51     );
  52
  53     volScalarField beta
  54     (
  55         IOobject
  56         (
  57             "beta",
  58             runTime.timeName(),
  59             mesh,
  60             IOobject::NO_READ,
  61             IOobject::NO_WRITE
  62         ),
  63         scalar(1) - alpha
  64         //,alpha.boundaryField().types()
  65     );
  66
  67     Info<< "Reading field p\n" << endl;
  68     volScalarField p
  69     (
  70         IOobject
  71         (
  72             "p",
  73             runTime.timeName(),
  74             mesh,
  75             IOobject::MUST_READ,
  76             IOobject::AUTO_WRITE
  77         ),
  78         mesh
  79     );
  80
  81     volVectorField U
  82     (
  83         IOobject
  84         (
  85             "U",
  86             runTime.timeName(),
  87             mesh,
  88             IOobject::NO_READ,
  89             IOobject::AUTO_WRITE
  90         ),
  91         alpha*Ua + beta*Ub
  92     );
  93
  94     dimensionedScalar Cvm
  95     (
  96         transportProperties.lookup("Cvm")
  97     );
  98
  99     dimensionedScalar Cl
 100     (
 101         transportProperties.lookup("Cl")
 102     );
 103
 104     dimensionedScalar Ct
 105     (
 106         transportProperties.lookup("Ct")
 107     );
 108
 109     surfaceScalarField phi
 110     (
 111         IOobject
 112         (
 113             "phi",
 114             runTime.timeName(),
 115             mesh
 116         ),
 117         fvc::interpolate(alpha)*phia + fvc::interpolate(beta)*phib
 118     );
 119
 120     volScalarField rho
 121     (
 122         IOobject
 123         (
 124             "rho",
 125             runTime.timeName(),
 126             mesh
 127         ),
 128         alpha*rhoa + beta*rhob
 129     );
 130
 131     #include "createRASTurbulence.H"
 132
 133     Info<< "Calculating field DDtUa and DDtUb\n" << endl;
 134
 135     volVectorField DDtUa =
 136         fvc::ddt(Ua)
 137       + fvc::div(phia, Ua)
 138       - fvc::div(phia)*Ua;
 139
 140     volVectorField DDtUb =
 141         fvc::ddt(Ub)
 142       + fvc::div(phib, Ub)
 143       - fvc::div(phib)*Ub;
 144
 145
 146     Info<< "Calculating field g.h\n" << endl;
 147     volScalarField gh("gh", g & mesh.C());
 148
 149     IOdictionary interfacialProperties
 150     (
 151         IOobject
 152         (
 153             "interfacialProperties",
 154             runTime.constant(),
 155             mesh,
 156             IOobject::MUST_READ,
 157             IOobject::NO_WRITE
 158         )
 159     );
 160
 161     autoPtr<dragModel> draga = dragModel::New
 162     (
 163         interfacialProperties,
 164         alpha,
 165         phasea,
 166         phaseb
 167     );
 168
 169     autoPtr<dragModel> dragb = dragModel::New
 170     (
 171         interfacialProperties,
 172         beta,
 173         phaseb,
 174         phasea
 175     );
 176
 177     word dragPhase("blended");
 178     if (interfacialProperties.found("dragPhase"))
 179     {
 180         dragPhase = word(interfacialProperties.lookup("dragPhase"));
 181
 182         bool validDrag =
 183             dragPhase == "a" || dragPhase == "b" || dragPhase == "blended";
 184
 185         if (!validDrag)
 186         {
 187             FatalErrorIn(args.executable())
 188                 << "invalid dragPhase " << dragPhase
 189                 << exit(FatalError);
 190         }
 191     }
 192
 193     Info << "dragPhase is " << dragPhase << endl;
 194     kineticTheoryModel kineticTheory
 195     (
 196         phasea,
 197         Ub,
 198         alpha,
 199         draga
 200     );
 201
 202     surfaceScalarField rUaAf
 203     (
 204         IOobject
 205         (
 206             "rUaAf",
 207             runTime.timeName(),
 208             mesh,
 209             IOobject::NO_READ,
 210             IOobject::NO_WRITE
 211         ),
 212         mesh,
 213         dimensionedScalar("zero", dimensionSet(0, 0, 1, 0, 0), 0.0)
 214     );
 215
 216     surfaceScalarField ppMagf
 217     (
 218         IOobject
 219         (
 220             "ppMagf",
 221             runTime.timeName(),
 222             mesh,
 223             IOobject::NO_READ,
 224             IOobject::NO_WRITE
 225         ),
 226         mesh,
 227         dimensionedScalar("zero", dimensionSet(0, 2, -1, 0, 0), 0.0)
 228     );
 229
 230
 231     label pRefCell = 0;
 232     scalar pRefValue = 0.0;
 233     setRefCell(p, mesh.solutionDict().subDict("PISO"), pRefCell, pRefValue);