applications/solvers/incompressible/boundaryFoam/boundaryFoam.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | foam-extend: Open Source CFD
   4    \\    /   O peration     | Version:     3.2
   5     \\  /    A nd           | Web:         http://www.foam-extend.org
   6      \\/     M anipulation  | For copyright notice see file Copyright
   7 -------------------------------------------------------------------------------
   8 License
   9     This file is part of foam-extend.
  10
  11     foam-extend is free software: you can redistribute it and/or modify it
  12     under the terms of the GNU General Public License as published by the
  13     Free Software Foundation, either version 3 of the License, or (at your
  14     option) any later version.
  15
  16     foam-extend is distributed in the hope that it will be useful, but
  17     WITHOUT ANY WARRANTY; without even the implied warranty of
  18     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  19     General Public License for more details.
  20
  21     You should have received a copy of the GNU General Public License
  22     along with foam-extend.  If not, see <http://www.gnu.org/licenses/>.
  23
  24 Application
  25     boundaryFoam
  26
  27 Description
  28     Steady-state solver for 1D turbulent flow, typically to generate boundary
  29     layer conditions at an inlet, for use in a simulation.
  30
  31     Boundary layer code to calculate the U, k and epsilon distributions.
  32     Used to create inlet boundary conditions for experimental comparisons
  33     for which U and k have not been measured.
  34     Turbulence model is runtime selectable.
  35
  36 \*---------------------------------------------------------------------------*/
  37
  38 #include "fvCFD.H"
  39 #include "singlePhaseTransportModel.H"
  40 #include "RASModel.H"
  41 #include "wallFvPatch.H"
  42 #include "makeGraph.H"
  43
  44
  45 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  46
  47 int main(int argc, char *argv[])
  48 {
  49 #   include "setRootCase.H"
  50
  51 #   include "createTime.H"
  52 #   include "createMesh.H"
  53 #   include "createFields.H"
  54
  55 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  56
  57     Info<< "\nStarting time loop\n" << endl;
  58
  59     while (runTime.loop())
  60     {
  61         Info<< "Time = " << runTime.timeName() << nl << endl;
  62
  63         fvVectorMatrix divR = turbulence->divDevReff(U);
  64         divR.source() = flowMask & divR.source();
  65
  66         fvVectorMatrix UEqn
  67         (
  68             divR == gradP
  69         );
  70
  71         UEqn.relax();
  72
  73         UEqn.solve();
  74
  75
  76         // Correct driving force for a constant mass flow rate
  77
  78         dimensionedVector UbarStar = flowMask & U.weightedAverage(mesh.V());
  79
  80         U += (Ubar - UbarStar);
  81         gradP += (Ubar - UbarStar)/(1.0/UEqn.A())().weightedAverage(mesh.V());
  82
  83         label id = y.size() - 1;
  84
  85         scalar wallShearStress =
  86             flowDirection & turbulence->R()()[id] & wallNormal;
  87
  88         scalar yplusWall
  89 //            = Foam::sqrt(mag(wallShearStress))*y[id]/laminarTransport.nu()()[id];
  90             = Foam::sqrt(mag(wallShearStress))*y[id]/turbulence->nuEff()()[id];
  91
  92         Info<< "Uncorrected Ubar = " << (flowDirection & UbarStar.value())<< tab
  93             << "pressure gradient = " << (flowDirection & gradP.value()) << tab
  94             << "min y+ = " << yplusWall << endl;
  95
  96
  97         turbulence->correct();
  98
  99
 100         if (runTime.outputTime())
 101         {
 102             volSymmTensorField R
 103             (
 104                 IOobject
 105                 (
 106                     "R",
 107                     runTime.timeName(),
 108                     mesh,
 109                     IOobject::NO_READ,
 110                     IOobject::AUTO_WRITE
 111                 ),
 112                 turbulence->R()
 113             );
 114
 115             runTime.write();
 116
 117             const word& gFormat = runTime.graphFormat();
 118
 119             makeGraph(y, flowDirection & U, "Uf", gFormat);
 120
 121             makeGraph(y, laminarTransport.nu(), gFormat);
 122
 123             makeGraph(y, turbulence->k(), gFormat);
 124             makeGraph(y, turbulence->epsilon(), gFormat);
 125
 126             //makeGraph(y, flowDirection & R & flowDirection, "Rff", gFormat);
 127             //makeGraph(y, wallNormal & R & wallNormal, "Rww", gFormat);
 128             //makeGraph(y, flowDirection & R & wallNormal, "Rfw", gFormat);
 129
 130             //makeGraph(y, sqrt(R.component(tensor::XX)), "u", gFormat);
 131             //makeGraph(y, sqrt(R.component(tensor::YY)), "v", gFormat);
 132             //makeGraph(y, sqrt(R.component(tensor::ZZ)), "w", gFormat);
 133             makeGraph(y, R.component(tensor::XY), "uv", gFormat);
 134
 135             makeGraph(y, mag(fvc::grad(U)), "gammaDot", gFormat);
 136         }
 137
 138         Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
 139             << "  ClockTime = " << runTime.elapsedClockTime() << " s"
 140             << nl << endl;
 141     }
 142
 143     Info<< "End\n" << endl;
 144
 145     return 0;
 146 }
 147
 148
 149 // ************************************************************************* //