applications/solvers/multiphase/bubbleFoam/wallFunctions.H

   1 {
   2     labelList cellBoundaryFaceCount(epsilon.size(), 0);
   3
   4     scalar Cmu25 = ::pow(Cmu.value(), 0.25);
   5     scalar Cmu75 = ::pow(Cmu.value(), 0.75);
   6     scalar kappa_ = kappa.value();
   7     scalar nub_ = nub.value();
   8
   9     const fvPatchList& patches = mesh.boundary();
  10
  11     //- Initialise the near-wall P field to zero
  12     forAll(patches, patchi)
  13     {
  14         const fvPatch& currPatch = patches[patchi];
  15
  16         if (isA<wallFvPatch>(currPatch))
  17         {
  18             forAll(currPatch, facei)
  19             {
  20                 label faceCelli = currPatch.faceCells()[facei];
  21
  22                 epsilon[faceCelli] = 0.0;
  23                 G[faceCelli] = 0.0;
  24             }
  25         }
  26     }
  27
  28     //- Accumulate the wall face contributions to epsilon and G
  29     //  Increment cellBoundaryFaceCount for each face for averaging
  30     forAll(patches, patchi)
  31     {
  32         const fvPatch& currPatch = patches[patchi];
  33
  34         if (isA<wallFvPatch>(currPatch))
  35         {
  36             const scalarField& nutbw = nutb.boundaryField()[patchi];
  37
  38             scalarField magFaceGradU(mag(Ub.boundaryField()[patchi].snGrad()));
  39
  40             forAll(currPatch, facei)
  41             {
  42                 label faceCelli = currPatch.faceCells()[facei];
  43
  44                 scalar yPlus =
  45                     Cmu25*y[patchi][facei]
  46                     *::sqrt(k[faceCelli])
  47                     /nub_;
  48
  49
  50                 // For corner cells (with two boundary or more faces),
  51                 // epsilon and G in the near-wall cell are calculated
  52                 // as an average
  53
  54                 cellBoundaryFaceCount[faceCelli]++;
  55
  56                 epsilon[faceCelli] +=
  57                      Cmu75*::pow(k[faceCelli], 1.5)
  58                     /(kappa_*y[patchi][facei]);
  59
  60                 if (yPlus > 11.6)
  61                 {
  62                     G[faceCelli] +=
  63                         (nutbw[facei] + nub_)*magFaceGradU[facei]
  64                         *Cmu25*::sqrt(k[faceCelli])
  65                         /(kappa_*y[patchi][facei]);
  66                 }
  67             }
  68         }
  69     }
  70
  71
  72     // perform the averaging
  73
  74     forAll(patches, patchi)
  75     {
  76         const fvPatch& curPatch = patches[patchi];
  77
  78         if (isA<wallFvPatch>(curPatch))
  79         {
  80             forAll(curPatch, facei)
  81             {
  82                 label faceCelli = curPatch.faceCells()[facei];
  83
  84                 epsilon[faceCelli] /= cellBoundaryFaceCount[faceCelli];
  85                 G[faceCelli] /= cellBoundaryFaceCount[faceCelli];
  86             }
  87         }
  88     }
  89 }