applications/solvers/multiphase/settlingFoam/createFields.H

   1     Info<< "Reading field p\n" << endl;
   2     volScalarField p
   3     (
   4         IOobject
   5         (
   6             "p",
   7             runTime.timeName(),
   8             mesh,
   9             IOobject::MUST_READ,
  10             IOobject::AUTO_WRITE
  11         ),
  12         mesh
  13     );
  14
  15     Info<< "Reading field alpha\n" << endl;
  16     volScalarField alpha
  17     (
  18         IOobject
  19         (
  20             "alpha",
  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
  44     Info<< "Reading transportProperties\n" << endl;
  45
  46     IOdictionary transportProperties
  47     (
  48         IOobject
  49         (
  50             "transportProperties",
  51             runTime.constant(),
  52             mesh,
  53             IOobject::MUST_READ,
  54             IOobject::NO_WRITE
  55         )
  56     );
  57
  58
  59     dimensionedScalar rhoc
  60     (
  61         transportProperties.lookup("rhoc")
  62     );
  63
  64     dimensionedScalar rhod
  65     (
  66         transportProperties.lookup("rhod")
  67     );
  68
  69     dimensionedScalar muc
  70     (
  71         transportProperties.lookup("muc")
  72     );
  73
  74     dimensionedScalar plasticViscosityCoeff
  75     (
  76         transportProperties.lookup("plasticViscosityCoeff")
  77     );
  78
  79     dimensionedScalar plasticViscosityExponent
  80     (
  81         transportProperties.lookup("plasticViscosityExponent")
  82     );
  83
  84     dimensionedScalar yieldStressCoeff
  85     (
  86         transportProperties.lookup("yieldStressCoeff")
  87     );
  88
  89     dimensionedScalar yieldStressExponent
  90     (
  91         transportProperties.lookup("yieldStressExponent")
  92     );
  93
  94     dimensionedScalar yieldStressOffset
  95     (
  96         transportProperties.lookup("yieldStressOffset")
  97     );
  98
  99     Switch BinghamPlastic
 100     (
 101         transportProperties.lookup("BinghamPlastic")
 102     );
 103
 104     volScalarField rho
 105     (
 106         IOobject
 107         (
 108             "rho",
 109             runTime.timeName(),
 110             mesh,
 111             IOobject::NO_READ,
 112             IOobject::NO_WRITE
 113         ),
 114         (scalar(1) - alpha)*rhoc + alpha*rhod
 115     );
 116
 117     volScalarField Alpha
 118     (
 119         IOobject
 120         (
 121             "Alpha",
 122             runTime.timeName(),
 123             mesh,
 124             IOobject::NO_READ,
 125             IOobject::AUTO_WRITE
 126         ),
 127         alpha*rhod/rho,
 128         alpha.boundaryField().types()
 129     );
 130
 131 #   include "compressibleCreatePhi.H"
 132
 133
 134     label pRefCell = 0;
 135     scalar pRefValue = 0.0;
 136     setRefCell(p, mesh.solutionDict().subDict("PISO"), pRefCell, pRefValue);
 137
 138
 139     Info<< "Calculating field mul\n" << endl;
 140     volScalarField mul
 141     (
 142         IOobject
 143         (
 144             "mul",
 145             runTime.timeName(),
 146             mesh,
 147             IOobject::NO_READ,
 148             IOobject::AUTO_WRITE
 149         ),
 150         muc +
 151         plasticViscosity
 152         (
 153             plasticViscosityCoeff,
 154             plasticViscosityExponent,
 155             Alpha
 156         )
 157     );
 158
 159
 160     Info<< "Initialising field Vdj\n" << endl;
 161     volVectorField Vdj
 162     (
 163         IOobject
 164         (
 165             "Vdj",
 166             runTime.timeName(),
 167             mesh,
 168             IOobject::NO_READ,
 169             IOobject::AUTO_WRITE
 170         ),
 171         mesh,
 172         dimensionedVector("0.0", U.dimensions(), vector::zero),
 173         U.boundaryField().types()
 174     );
 175
 176
 177     Info<< "Selecting Drift-Flux model " << endl;
 178
 179     word VdjModel
 180     (
 181         transportProperties.lookup("VdjModel")
 182     );
 183
 184     Info<< tab << VdjModel << " selected\n" << endl;
 185
 186     const dictionary& VdjModelCoeffs
 187     (
 188         transportProperties.subDict(VdjModel + "Coeffs")
 189     );
 190
 191     dimensionedVector V0
 192     (
 193         VdjModelCoeffs.lookup("V0")
 194     );
 195
 196     dimensionedScalar a
 197     (
 198         VdjModelCoeffs.lookup("a")
 199     );
 200
 201     dimensionedScalar a1
 202     (
 203         VdjModelCoeffs.lookup("a1")
 204     );
 205
 206     dimensionedScalar alphaMin
 207     (
 208         VdjModelCoeffs.lookup("alphaMin")
 209     );
 210
 211
 212
 213     IOdictionary RASProperties
 214     (
 215         IOobject
 216         (
 217             "RASProperties",
 218             runTime.constant(),
 219             mesh,
 220             IOobject::MUST_READ,
 221             IOobject::NO_WRITE
 222         )
 223     );
 224
 225
 226     Switch turbulence
 227     (
 228         RASProperties.lookup("turbulence")
 229     );
 230
 231     const dictionary& kEpsilonCoeffs
 232     (
 233         RASProperties.subDict("kEpsilonCoeffs")
 234     );
 235
 236
 237     scalar Cmu
 238     (
 239         readScalar(kEpsilonCoeffs.lookup("Cmu"))
 240     );
 241
 242     scalar C1
 243     (
 244         readScalar(kEpsilonCoeffs.lookup("C1"))
 245     );
 246
 247     scalar C2
 248     (
 249         readScalar(kEpsilonCoeffs.lookup("C2"))
 250     );
 251
 252     scalar C3
 253     (
 254         readScalar(kEpsilonCoeffs.lookup("C3"))
 255     );
 256
 257     scalar alphak
 258     (
 259         readScalar(kEpsilonCoeffs.lookup("alphak"))
 260     );
 261
 262     scalar alphaEps
 263     (
 264         readScalar(kEpsilonCoeffs.lookup("alphaEps"))
 265     );
 266
 267     const dictionary& wallFunctionCoeffs
 268     (
 269         RASProperties.subDict("wallFunctionCoeffs")
 270     );
 271
 272     scalar kappa
 273     (
 274         readScalar(wallFunctionCoeffs.lookup("kappa"))
 275     );
 276
 277     scalar E
 278     (
 279         readScalar(wallFunctionCoeffs.lookup("E"))
 280     );
 281
 282     nearWallDist y(mesh);
 283
 284     Info<< "Reading field k\n" << endl;
 285     volScalarField k
 286     (
 287         IOobject
 288         (
 289             "k",
 290             runTime.timeName(),
 291             mesh,
 292             IOobject::MUST_READ,
 293             IOobject::AUTO_WRITE
 294         ),
 295         mesh
 296     );
 297
 298     Info<< "Reading field epsilon\n" << endl;
 299     volScalarField epsilon
 300     (
 301         IOobject
 302         (
 303             "epsilon",
 304             runTime.timeName(),
 305             mesh,
 306             IOobject::MUST_READ,
 307             IOobject::AUTO_WRITE
 308         ),
 309         mesh
 310     );
 311
 312     Info<< "Calculating field mut\n" << endl;
 313     volScalarField mut
 314     (
 315         IOobject
 316         (
 317             "mut",
 318             runTime.timeName(),
 319             mesh,
 320             IOobject::NO_READ,
 321             IOobject::AUTO_WRITE
 322         ),
 323         Cmu*rho*sqr(k)/epsilon
 324     );
 325
 326
 327     Info<< "Calculating field mu\n" << endl;
 328     volScalarField mu
 329     (
 330         IOobject
 331         (
 332             "mu",
 333             runTime.timeName(),
 334             mesh,
 335             IOobject::NO_READ,
 336             IOobject::AUTO_WRITE
 337         ),
 338         mut + mul
 339     );
 340
 341     Info<< "Calculating field (g.h)f\n" << endl;
 342     surfaceScalarField ghf = surfaceScalarField("gh", g & mesh.Cf());