applications/solvers/multiphase/bubbleFoam/UEqns.H

   1 fvVectorMatrix UaEqn(Ua, Ua.dimensions()*dimVol/dimTime);
   2 fvVectorMatrix UbEqn(Ub, Ub.dimensions()*dimVol/dimTime);
   3
   4 {
   5     volTensorField Rca(-nuEffa*(T(fvc::grad(Ua))));
   6     Rca = Rca + (2.0/3.0)*sqr(Ct)*I*k - (2.0/3.0)*I*tr(Rca);
   7
   8     surfaceScalarField phiRa
   9     (
  10       - fvc::interpolate(nuEffa)
  11        *mesh.magSf()*fvc::snGrad(alpha)/fvc::interpolate(alpha + scalar(0.001))
  12     );
  13
  14     UaEqn =
  15     (
  16         (scalar(1) + Cvm*rhob*beta/rhoa)*
  17         (
  18             fvm::ddt(Ua)
  19           + fvm::div(phia, Ua, "div(phia,Ua)")
  20           - fvm::Sp(fvc::div(phia), Ua)
  21         )
  22
  23       - fvm::laplacian(nuEffa, Ua)
  24       + fvc::div(Rca)
  25
  26       + fvm::div(phiRa, Ua, "div(phia,Ua)")
  27       - fvm::Sp(fvc::div(phiRa), Ua)
  28       + (fvc::grad(alpha)/(fvc::average(alpha) + scalar(0.001)) & Rca)
  29      ==
  30     //  g                        // Buoyancy term transfered to p-equation
  31       - fvm::Sp(beta/rhoa*dragCoef, Ua)
  32     //+ beta/rhoa*dragCoef*Ub    // Explicit drag transfered to p-equation
  33       - beta/rhoa*(liftCoeff - Cvm*rhob*DDtUb)
  34     );
  35
  36     UaEqn.relax();
  37
  38
  39     volTensorField Rcb(-nuEffb*T(fvc::grad(Ub)));
  40     Rcb = Rcb + (2.0/3.0)*I*k - (2.0/3.0)*I*tr(Rcb);
  41
  42     surfaceScalarField phiRb
  43     (
  44        - fvc::interpolate(nuEffb)
  45         *mesh.magSf()*fvc::snGrad(beta)/fvc::interpolate(beta + scalar(0.001))
  46     );
  47
  48     UbEqn =
  49     (
  50         (scalar(1) + Cvm*rhob*alpha/rhob)*
  51         (
  52             fvm::ddt(Ub)
  53           + fvm::div(phib, Ub, "div(phib,Ub)")
  54           - fvm::Sp(fvc::div(phib), Ub)
  55         )
  56
  57       - fvm::laplacian(nuEffb, Ub)
  58       + fvc::div(Rcb)
  59
  60       + fvm::div(phiRb, Ub, "div(phib,Ub)")
  61       - fvm::Sp(fvc::div(phiRb), Ub)
  62
  63       + (fvc::grad(beta)/(fvc::average(beta) + scalar(0.001)) & Rcb)
  64      ==
  65     //  g                        // Buoyancy term transfered to p-equation
  66       - fvm::Sp(alpha/rhob*dragCoef, Ub)
  67     //+ alpha/rhob*dragCoef*Ua   // Explicit drag transfered to p-equation
  68       + alpha/rhob*(liftCoeff + Cvm*rhob*DDtUa)
  69     );
  70
  71     UbEqn.relax();
  72 }