DARhoSimpleFoam.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2 
3  DAFoam : Discrete Adjoint with OpenFOAM
4  Version : v4
5 
6  This class is modified from OpenFOAM's source code
7  applications/solvers/compressible/rhoSimpleFoam
8 
9  OpenFOAM: The Open Source CFD Toolbox
10 
11  Copyright (C): 2011-2016 OpenFOAM Foundation
12 
13  OpenFOAM License:
14 
15  OpenFOAM is free software: you can redistribute it and/or modify it
16  under the terms of the GNU General Public License as published by
17  the Free Software Foundation, either version 3 of the License, or
18  (at your option) any later version.
19 
20  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
21  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
22  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
23  for more details.
24 
25  You should have received a copy of the GNU General Public License
26  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
27 
28 \*---------------------------------------------------------------------------*/
29 
30 #include "DARhoSimpleFoam.H"
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 
37 defineTypeNameAndDebug(DARhoSimpleFoam, 0);
38 addToRunTimeSelectionTable(DASolver, DARhoSimpleFoam, dictionary);
39 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
40 
42  char* argsAll,
43  PyObject* pyOptions)
44  : DASolver(argsAll, pyOptions),
45  simplePtr_(nullptr),
46  pThermoPtr_(nullptr),
47  pPtr_(nullptr),
48  rhoPtr_(nullptr),
49  UPtr_(nullptr),
50  phiPtr_(nullptr),
51  pressureControlPtr_(nullptr),
52  turbulencePtr_(nullptr),
53  daTurbulenceModelPtr_(nullptr),
54  daFvSourcePtr_(nullptr),
55  fvSourcePtr_(nullptr),
56  fvSourceEnergyPtr_(nullptr),
57  initialMass_(dimensionedScalar("initialMass", dimensionSet(1, 0, 0, 0, 0, 0, 0), 0.0)),
58  MRFPtr_(nullptr)
59 {
60 }
61 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
62 
64 {
65  /*
66  Description:
67  Initialize variables for DASolver
68  */
69 
70  Info << "Initializing fields for DARhoSimpleFoam" << endl;
71  Time& runTime = runTimePtr_();
72  fvMesh& mesh = meshPtr_();
73  argList& args = argsPtr_();
75 #include "createFieldsRhoSimple.H"
76 
77  // read the RAS model from constant/turbulenceProperties
78  const word turbModelName(
79  IOdictionary(
80  IOobject(
81  "turbulenceProperties",
82  mesh.time().constant(),
83  mesh,
84  IOobject::MUST_READ,
85  IOobject::NO_WRITE,
86  false))
87  .subDict("RAS")
88  .lookup("RASModel"));
90 
91 #include "createAdjoint.H"
92 
93  // initialize fvSource and compute the source term
94  const dictionary& allOptions = daOptionPtr_->getAllOptions();
95  if (allOptions.subDict("fvSource").toc().size() != 0)
96  {
97  hasFvSource_ = 1;
98  Info << "Initializing DASource" << endl;
99  word sourceName = allOptions.subDict("fvSource").toc()[0];
100  word fvSourceType = allOptions.subDict("fvSource").subDict(sourceName).getWord("type");
102  fvSourceType, mesh, daOptionPtr_(), daModelPtr_(), daIndexPtr_()));
103  }
104 }
105 
107 {
108  /*
109  Description:
110  Call the primal solver to get converged state variables
111  */
112 
113 #include "createRefsRhoSimple.H"
114 #include "createFvOptions.H"
115 
116  // call correctNut, this is equivalent to turbulence->validate();
117  daTurbulenceModelPtr_->updateIntermediateVariables();
118 
119  Info << "\nStarting time loop\n"
120  << endl;
121 
122  while (this->loop(runTime)) // using simple.loop() will have seg fault in parallel
123  {
124 
125  if (printToScreen_)
126  {
127  Info << "Time = " << runTime.timeName() << nl << endl;
128  }
129 
130  p.storePrevIter();
131  rho.storePrevIter();
132 
133  // Pressure-velocity SIMPLE corrector
134 #include "UEqnRhoSimple.H"
135 #include "EEqnRhoSimple.H"
136 #include "pEqnRhoSimple.H"
137 
139 
140  // calculate all functions
142  // calculate yPlus
144  // compute the regression model and print the feature
145  regModelFail_ = daRegressionPtr_->compute();
146  daRegressionPtr_->printInputInfo(printToScreen_);
147  // print run time
149 
150  runTime.write();
151  }
152 
153  // write the mesh to files
154  mesh.write();
155 
156  // write associated fields such as URel
157  this->writeAssociatedFields();
158 
159  Info << "End\n"
160  << endl;
161 
162  return this->checkPrimalFailure();
163 }
164 
165 } // End namespace Foam
166 
167 // ************************************************************************* //
Foam::DASolver::checkPrimalFailure
label checkPrimalFailure()
check whether the primal fails based on residual and regression fail flag
Definition: DASolver.C:2472
Foam::DASolver::loop
label loop(Time &runTime)
return whether to loop the primal solution, similar to runTime::loop() except we don't do file IO
Definition: DASolver.C:147
Foam::DASolver::printToScreen_
label printToScreen_
whether to print primal information to the screen
Definition: DASolver.H:150
Foam::DASolver::daRegressionPtr_
autoPtr< DARegression > daRegressionPtr_
DARegression pointer.
Definition: DASolver.H:111
Foam::DASolver::daModelPtr_
autoPtr< DAModel > daModelPtr_
DAModel pointer.
Definition: DASolver.H:84
UEqnRhoSimple.H
Foam::DASolver::daIndexPtr_
autoPtr< DAIndex > daIndexPtr_
DAIndex pointer.
Definition: DASolver.H:87
Foam::DARhoSimpleFoam::daTurbulenceModelPtr_
autoPtr< DATurbulenceModel > daTurbulenceModelPtr_
DATurbulenceModel pointer.
Definition: DARhoSimpleFoam.H:83
initialMass_
initialMass_
Definition: createRefsRhoSimpleC.H:13
Foam::DASolver
Definition: DASolver.H:54
Foam::addToRunTimeSelectionTable
addToRunTimeSelectionTable(DAFunction, DAFunctionForce, dictionary)
DARhoSimpleFoam.H
Foam::DASolver::daOptionPtr_
autoPtr< DAOption > daOptionPtr_
DAOption pointer.
Definition: DASolver.H:81
createRefsRhoSimple.H
Foam::DATurbulenceModel::New
static autoPtr< DATurbulenceModel > New(const word modelType, const fvMesh &mesh, const DAOption &daOption)
Definition: DATurbulenceModel.C:161
createAdjoint.H
Foam::DARhoSimpleFoam::hasFvSource_
label hasFvSource_
whether to have fvSource term
Definition: DARhoSimpleFoam.H:95
mesh
fvMesh & mesh
Definition: createRefsHeatTransfer.H:4
p
volScalarField & p
Definition: createRefsPimpleDyM.H:6
pEqnRhoSimple.H
createSimpleControlPython.H
Foam::DAFvSource::New
static autoPtr< DAFvSource > New(const word modelType, const fvMesh &mesh, const DAOption &daOption, const DAModel &daModel, const DAIndex &daIndex)
Definition: DAFvSource.C:47
Foam::DARhoSimpleFoam::initSolver
virtual void initSolver()
initialize fields and variables
Definition: DARhoSimpleFoam.C:63
Foam::defineTypeNameAndDebug
defineTypeNameAndDebug(DAFunction, 0)
Foam
Definition: checkGeometry.C:32
createFieldsRhoSimple.H
Foam::DARhoSimpleFoam::daFvSourcePtr_
autoPtr< DAFvSource > daFvSourcePtr_
DASource pointer.
Definition: DARhoSimpleFoam.H:86
Foam::DASolver::writeAssociatedFields
void writeAssociatedFields()
write associated fields such as relative velocity
Definition: DASolver.C:2526
rho
volScalarField & rho
Definition: createRefsRhoPimple.H:8
allOptions
const dictionary & allOptions
Definition: createRefsPimpleDyM.H:14
Foam::DASolver::regModelFail_
label regModelFail_
whether the regModel compute fails
Definition: DASolver.H:153
Foam::DARhoSimpleFoam::solvePrimal
virtual label solvePrimal()
solve the primal equations
Definition: DARhoSimpleFoam.C:106
Foam::DASolver::runTimePtr_
autoPtr< Time > runTimePtr_
runTime pointer
Definition: DASolver.H:75
Foam::DASolver::meshPtr_
autoPtr< fvMesh > meshPtr_
fvMesh pointer
Definition: DASolver.H:78
Foam::DASolver::argsPtr_
autoPtr< argList > argsPtr_
args pointer
Definition: DASolver.H:72
args
Foam::argList & args
Definition: setRootCasePython.H:42
runTime
Time & runTime
Definition: createRefsHeatTransfer.H:1
EEqnRhoSimple.H
Foam::DARhoSimpleFoam::DARhoSimpleFoam
DARhoSimpleFoam(char *argsAll, PyObject *pyOptions)
Definition: DARhoSimpleFoam.C:41
Foam::DASolver::printElapsedTime
void printElapsedTime(const Time &runTime, const label printToScreen)
Definition: DASolver.H:256
Foam::DASolver::calcAllFunctions
void calcAllFunctions(label print=0)
calculate the values of all objective functions and print them to screen
Definition: DASolver.C:204