adjointSensitivityIncompressible.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2007-2020 PCOpt/NTUA
9  Copyright (C) 2013-2020 FOSS GP
10  Copyright (C) 2019-2021 OpenCFD Ltd.
11 -------------------------------------------------------------------------------
12 License
13  This file is part of OpenFOAM.
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 "runTimeSelectionTables.H"
34 #include "wallFvPatch.H"
35 #include "fvOptions.H"
36 
37 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
38 
39 namespace Foam
40 {
41 
42 namespace incompressible
43 {
44 
45 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
46 
47 defineTypeNameAndDebug(adjointSensitivity, 0);
48 defineRunTimeSelectionTable(adjointSensitivity, dictionary);
49 
50 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
51 
52 adjointSensitivity::adjointSensitivity
53 (
54  const fvMesh& mesh,
55  const dictionary& dict,
56  incompressibleVars& primalVars,
57  incompressibleAdjointVars& adjointVars,
59 )
60 :
62  derivatives_(0),
63  primalVars_(primalVars),
64  adjointVars_(adjointVars),
65  objectiveManager_(objectiveManager)
66 {}
67 
68 
69 // * * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * //
70 
72 (
73  const fvMesh& mesh,
74  const dictionary& dict,
75  incompressibleVars& primalVars,
76  incompressibleAdjointVars& adjointVars,
78 )
79 {
80  const word modelType(dict.get<word>("type"));
81 
82  Info<< "adjointSensitivity type : " << modelType << endl;
83 
84  auto* ctorPtr = dictionaryConstructorTable(modelType);
85 
86  if (!ctorPtr)
87  {
89  (
90  dict,
91  "adjointSensitivity",
92  modelType,
93  *dictionaryConstructorTablePtr_
94  ) << exit(FatalIOError);
95  }
96 
98  (
99  ctorPtr
100  (
101  mesh,
102  dict,
103  primalVars,
104  adjointVars,
106  )
107  );
108 }
109 
110 
111 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * //
112 
114 {
116  write(type());
117  return derivatives_;
118 }
119 
120 
122 {
123  return derivatives_;
124 }
125 
126 
128 {
129  derivatives_ = scalar(0);
130  if (fieldSensPtr_)
131  {
132  fieldSensPtr_().primitiveFieldRef() = scalar(0);
133  }
134 }
135 
136 
137 void adjointSensitivity::write(const word& baseName)
138 {
139  sensitivity::write(baseName);
140 }
141 
142 
144 {
145  // Term depending on the adjoint turbulence model
147  (
149  );
150  tmp<volTensorField> tturbulenceTerm(adjointRAS->FISensitivityTerm());
151  volTensorField& turbulenceTerm = tturbulenceTerm.ref();
152 
153  // nu effective
154  tmp<volScalarField> tnuEff(adjointRAS->nuEff());
155  const volScalarField& nuEff = tnuEff();
156 
157  tmp<volTensorField> tflowTerm
158  (
159  new volTensorField
160  (
161  IOobject
162  (
163  "flowTerm",
164  mesh_.time().timeName(),
165  mesh_,
168  ),
169  mesh_,
171  )
172  );
173  volTensorField& flowTerm = tflowTerm.ref();
174 
175  const volScalarField& p = primalVars_.p();
176  const volVectorField& U = primalVars_.U();
177  const volScalarField& pa = adjointVars_.pa();
178  const volVectorField& Ua = adjointVars_.Ua();
179  volTensorField gradU(fvc::grad(U));
180  volTensorField gradUa(fvc::grad(Ua));
181 
182  // Explicitly correct the boundary gradient to get rid of
183  // the tangential component
184  forAll(mesh_.boundary(), patchI)
185  {
186  const fvPatch& patch = mesh_.boundary()[patchI];
187  if (isA<wallFvPatch>(patch))
188  {
189  tmp<vectorField> tnf = mesh_.boundary()[patchI].nf();
190  const vectorField& nf = tnf();
191  gradU.boundaryFieldRef()[patchI] =
192  nf*U.boundaryField()[patchI].snGrad();
193  //gradUa.boundaryField()[patchI] =
194  // nf*Ua.boundaryField()[patchI].snGrad();
195  }
196  }
197 
198  volTensorField stress(nuEff*(gradU + T(gradU)));
199  autoPtr<volVectorField> stressXPtr
200  (
201  createZeroFieldPtr<vector>(mesh_, "stressX", stress.dimensions())
202  );
203  autoPtr<volVectorField> stressYPtr
204  (
205  createZeroFieldPtr<vector>(mesh_, "stressY", stress.dimensions())
206  );
207  autoPtr<volVectorField> stressZPtr
208  (
209  createZeroFieldPtr<vector>(mesh_, "stressZ", stress.dimensions())
210  );
211 
212  stressXPtr().replace(0, stress.component(0));
213  stressXPtr().replace(1, stress.component(1));
214  stressXPtr().replace(2, stress.component(2));
215 
216  stressYPtr().replace(0, stress.component(3));
217  stressYPtr().replace(1, stress.component(4));
218  stressYPtr().replace(2, stress.component(5));
219 
220  stressZPtr().replace(0, stress.component(6));
221  stressZPtr().replace(1, stress.component(7));
222  stressZPtr().replace(2, stress.component(8));
223 
224  volTensorField gradStressX(fvc::grad(stressXPtr()));
225  volTensorField gradStressY(fvc::grad(stressYPtr()));
226  volTensorField gradStressZ(fvc::grad(stressZPtr()));
227 
228  // Contribution from objective functions and constraints
229  volTensorField objectiveContributions
230  (
231  IOobject
232  (
233  "objectiveContributions",
234  mesh_.time().timeName(),
235  mesh_,
238  ),
239  mesh_,
241  );
243  forAll(functions, funcI)
244  {
245  objectiveContributions +=
246  functions[funcI].weight()
247  *functions[funcI].gradDxDbMultiplier();
248  }
249 
250  // Note:
251  // term4 (Ua & grad(stress)) is numerically tricky. Its div leads to third
252  // order spatial derivs in E-SI based computations Applying the product
253  // derivative rule (putting Ua inside the grad) gives better results in
254  // NACA0012, SA, WF. However, the original formulation should be kept at
255  // the boundary in order to respect the Ua boundary conditions (necessary
256  // for E-SI to give the same sens as FI). A mixed approach is hence
257  // followed
258  volTensorField term4
259  (
260  - nuEff*(gradUa & (gradU + T(gradU)))
261  + fvc::grad(nuEff * Ua & (gradU + T(gradU)))
262  );
263 
264  forAll(mesh_.boundary(), pI)
265  {
266  if (!isA<coupledFvPatch>(mesh_.boundary()[pI]))
267  {
268  term4.boundaryFieldRef()[pI] =
269  Ua.component(0)().boundaryField()[pI]
270  *gradStressX.boundaryField()[pI]
271  + Ua.component(1)().boundaryField()[pI]
272  *gradStressY.boundaryField()[pI]
273  + Ua.component(2)().boundaryField()[pI]
274  *gradStressZ.boundaryField()[pI];
275  }
276  }
277 
278  const autoPtr<ATCModel>& ATCModel =
280  (
282  ).getATCModel();
283 
284  // Compute dxdb multiplier
285  flowTerm =
286  // Term 1, ATC
288  // Term 2
289  - fvc::grad(p) * Ua
290  // Term 3
291  - nuEff*(gradU & (gradUa + T(gradUa)))
292  // Term 4
293  + term4
294  // Term 5
295  + (pa * gradU)
296  // Term 6, from the adjoint turbulence model
297  + turbulenceTerm.T()
298  // Term 7, term from objective functions
299  + objectiveContributions;
300 
301  flowTerm.correctBoundaryConditions();
302 
303  return (tflowTerm);
304 }
305 
306 
308 {
310  volTensorField& gradDxDbMult = tgradDxDbMult.ref();
311 
312  tmp<volVectorField> tadjointMeshMovementSource
313  (
314  new volVectorField
315  (
316  IOobject
317  (
318  "adjointMeshMovementSource",
319  mesh_.time().timeName(),
320  mesh_,
323  ),
324  mesh_,
325  dimensionedVector(gradDxDbMult.dimensions()/dimLength, Zero)
326  )
327  );
328 
329  volVectorField& source = tadjointMeshMovementSource.ref();
330 
331  source -= fvc::div(gradDxDbMult.T());
332 
333  // Terms from fvOptions
335  (
336  source.primitiveFieldRef(), adjointVars_.solverName()
337  );
338 
339  return (tadjointMeshMovementSource);
340 }
341 
342 
343 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
344 
345 } // End namespace incompressible
346 } // End namespace Foam
347 
348 // ************************************************************************* //
Foam::IOobject::NO_WRITE
Definition: IOobject.H:195
Foam::incompressible::adjointSensitivity::derivatives_
scalarField derivatives_
Definition: adjointSensitivityIncompressible.H:83
Foam::objectiveManager
class for managing incompressible objective functions.
Definition: objectiveManager.H:54
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
Foam::GeometricField::component
tmp< GeometricField< cmptType, PatchField, GeoMesh > > component(const direction) const
Return a component of the field.
Foam::variablesSet::solverName
const word & solverName() const
Return solver name.
Definition: variablesSet.C:84
p
volScalarField & p
Definition: createFieldRefs.H:8
Foam::fvc::grad
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:54
Foam::incompressible::adjointSensitivity::calculateSensitivities
virtual const scalarField & calculateSensitivities()
Calculates and returns sensitivity fields.
Definition: adjointSensitivityIncompressible.C:113
Foam::incompressibleAdjointSolver
Base class for incompressibleAdjoint solvers.
Definition: incompressibleAdjointSolver.H:53
fvOptions.H
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::incompressible::defineTypeNameAndDebug
defineTypeNameAndDebug(adjointEikonalSolver, 0)
Foam::dimLength
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:52
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::fv::options::New
static options & New(const fvMesh &mesh)
Construct fvOptions and register to database if not present.
Definition: fvOptions.C:103
boundaryAdjointContribution.H
Foam::Time::timeName
static word timeName(const scalar t, const int precision=precision_)
Definition: Time.C:780
Foam::incompressible::adjointSensitivity::assembleSensitivities
virtual void assembleSensitivities()=0
Assemble sensitivities.
wallFvPatch.H
Foam::fvc::div
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:49
Foam::FatalIOError
IOerror FatalIOError
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
Foam::incompressibleVars::p
const volScalarField & p() const
Return const reference to pressure.
Definition: incompressibleVars.C:305
Foam::incompressible::defineRunTimeSelectionTable
defineRunTimeSelectionTable(adjointSensitivity, dictionary)
Foam::ATCModel
Base class for selecting the adjoint transpose convection model. Inherits from regIOobject to add loo...
Definition: ATCModel.H:60
Foam::incompressibleAdjointVars
Class including all adjoint fields for incompressible flows.
Definition: incompressibleAdjointVars.H:52
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
incompressibleAdjointSolver.H
Foam::dimensionedVector
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
Definition: dimensionedVector.H:50
Foam::dimTime
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
Foam::tmp::ref
T & ref() const
Definition: tmpI.H:227
FatalIOErrorInLookup
#define FatalIOErrorInLookup(ios, lookupTag, lookupName, lookupTable)
Report an error message using Foam::FatalIOError.
Definition: error.H:478
Foam::sensitivity::mesh_
const fvMesh & mesh_
Definition: sensitivity.H:69
Foam::Field< scalar >
Foam::pow3
dimensionedScalar pow3(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:89
Foam::Info
messageStream Info
Information stream (stdout output on master, null elsewhere)
Foam::incompressible::adjointSensitivity::New
static autoPtr< adjointSensitivity > New(const fvMesh &mesh, const dictionary &dict, incompressibleVars &primalVars, incompressibleAdjointVars &adjointVars, objectiveManager &objectiveManager)
Return a reference to the selected turbulence model.
Definition: adjointSensitivityIncompressible.C:72
Foam::fvPatch
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:65
Foam::ATCModel::getFISensitivityTerm
virtual tmp< volTensorField > getFISensitivityTerm() const =0
Get the FI sensitivity derivatives term coming from the ATC.
Foam::T
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
Definition: FieldFieldFunctions.C:58
Foam::incompressibleAdjointMeanFlowVars::pa
const volScalarField & pa() const
Return const reference to pressure.
Definition: incompressibleAdjointMeanFlowVars.C:147
Foam::incompressibleVars::U
const volVectorField & U() const
Return const reference to velocity.
Definition: incompressibleVars.C:331
Foam::PtrList
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: List.H:59
Foam::incompressible::adjointSensitivity::clearSensitivities
virtual void clearSensitivities()
Zero sensitivity fields and their constituents.
Definition: adjointSensitivityIncompressible.C:127
Foam::objectRegistry::lookupObject
const Type & lookupObject(const word &name, const bool recursive=false) const
Definition: objectRegistryTemplates.C:434
Foam::incompressibleAdjointMeanFlowVars::Ua
const volVectorField & Ua() const
Return const reference to velocity.
Definition: incompressibleAdjointMeanFlowVars.C:173
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::sensitivity::fieldSensPtr_
autoPtr< volScalarField > fieldSensPtr_
Definition: sensitivity.H:74
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
Foam::incompressible::adjointSensitivity::adjointVars_
incompressibleAdjointVars & adjointVars_
Definition: adjointSensitivityIncompressible.H:85
Foam::fvMesh
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:85
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::fv::optionList::postProcessSens
void postProcessSens(Field< Type > &sensField, const word &fieldName=word::null, const word &designVariablesName=word::null)
Post process sensitivity field related to the fvOption.
Definition: fvOptionListTemplates.C:398
Foam::incompressible::adjointSensitivity::primalVars_
incompressibleVars & primalVars_
Definition: adjointSensitivityIncompressible.H:84
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
Foam::GeometricField::T
tmp< GeometricField< Type, PatchField, GeoMesh > > T() const
Return transpose (only if it is a tensor field)
Definition: GeometricField.C:1046
Foam::fvMesh::boundary
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:685
Foam::autoPtr
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: HashPtrTable.H:53
U
U
Definition: pEqn.H:72
Foam::sqr
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Definition: dimensionedSymmTensor.C:51
adjointSensitivityIncompressible.H
Foam::foamVersion::patch
const std::string patch
OpenFOAM patch number as a std::string.
runTimeSelectionTables.H
Macros to ease declaration of run-time selection tables.
Foam::GeometricField::ref
Internal & ref(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
Definition: GeometricField.C:749
Foam::GeometricField::boundaryFieldRef
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
Definition: GeometricField.C:783
Foam::incompressibleAdjointVars::adjointTurbulence
const autoPtr< incompressibleAdjoint::adjointRASModel > & adjointTurbulence() const
Return const reference to the adjointRASModel.
Definition: incompressibleAdjointVars.C:70
Foam::incompressible::adjointSensitivity::adjointMeshMovementSource
tmp< volVectorField > adjointMeshMovementSource()
Compute source term for adjoint mesh movement equation.
Definition: adjointSensitivityIncompressible.C:307
Foam::incompressible::adjointSensitivity::getSensitivities
const scalarField & getSensitivities() const
Returns the sensitivity fields.
Definition: adjointSensitivityIncompressible.C:121
Foam::type
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:590
Foam::fvMesh::time
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:280
Foam::sensitivity::write
virtual void write(const word &baseName=word::null)
Write sensitivity fields.
Definition: sensitivity.C:77
Foam::objectiveManager::adjointSolverName
const word & adjointSolverName() const
Return name of the adjointSolver.
Definition: objectiveManager.C:308
Foam::GeometricField< tensor, fvPatchField, volMesh >
Foam::IOobject::NO_READ
Definition: IOobject.H:188
Foam::sensitivity
Abstract base class for adjoint sensitivities.
Definition: sensitivity.H:63
Foam::dimensionedTensor
dimensioned< tensor > dimensionedTensor
Dimensioned tensor obtained from generic dimensioned type.
Definition: dimensionedTensor.H:52
Foam::incompressible::adjointSensitivity::computeGradDxDbMultiplier
tmp< volTensorField > computeGradDxDbMultiplier()
Definition: adjointSensitivityIncompressible.C:143
Foam::incompressibleVars
Base class for solution control classes.
Definition: incompressibleVars.H:54
Foam::GeometricField::boundaryField
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Definition: GeometricFieldI.H:62
Foam::incompressible::adjointSensitivity::write
virtual void write(const word &baseName=word::null)
Write sensitivity fields.
Definition: adjointSensitivityIncompressible.C:137
Foam::objectiveManager::getObjectiveFunctions
PtrList< objective > & getObjectiveFunctions()
Return reference to objective functions.
Definition: objectiveManager.C:296
Foam::incompressible::adjointSensitivity::objectiveManager_
objectiveManager & objectiveManager_
Definition: adjointSensitivityIncompressible.H:86