adjointkOmegaSST.H
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) 2014-2022 PCOpt/NTUA
9 Copyright (C) 2014-2022 FOSS GP
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27
28Class
29 Foam::incompressibleAdjoint::adjointRASModels::adjointkOmegaSST
30
31Description
32 Continuous adjoint to the kOmegaSST turbulence model for incompressible
33 flows.
34
35 Reference:
36 \verbatim
37 The code is based on the following reference, with a number of
38 changes in the numerical implementation
39
40 Kavvadias, I., Papoutsis-Kiachagias, E., Dimitrakopoulos, G., &
41 Giannakoglou, K. (2014).
42 The continuous adjoint approach to the k–ω SST turbulence model
43 with applications in shape optimization
44 Engineering Optimization, 47(11), 1523-1542.
45 https://doi.org/10.1080/0305215X.2014.979816
46 \endverbatim
47
48SourceFiles
49 adjointkOmegaSST.C
50
51\*---------------------------------------------------------------------------*/
52
53#ifndef adjointkOmegaSST_H
54#define adjointkOmegaSST_H
55
56#include "adjointRASModel.H"
57#include "wallDist.H"
58
59// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
60
61namespace Foam
62{
63namespace incompressibleAdjoint
64{
65namespace adjointRASModels
66{
67
68/*---------------------------------------------------------------------------*\
69 Class adjointkOmegaSST Declaration
70\*---------------------------------------------------------------------------*/
73:
74 public adjointRASModel
75{
76 // Private Member Functions
77
78 //- No copy construct
79 adjointkOmegaSST(const adjointkOmegaSST&) = delete;
80
81 //- No copy assignment
82 void operator=(const adjointkOmegaSST&) = delete;
83
84
85protected:
86
87 // Protected data
88
89 // Primal Model coefficients
109
110 //- Flag to include the F3 term
111 Switch F3_;
112
113
114 // Fields
115
116 // Primal Fields
117
118 //- Wall distance
119 // Note: reference to the distance known by the primal model
120 const volScalarField& y_;
121
122 //- Cached primal gradient fields
126
127 //- Primal cached fields involved in the solution of the
128 // adjoint equations
129 // Cached to reduce the computational cost
137
138 // Model Field coefficients
143
144
145 // Switch fields
146
147 // Switch fields for the differentiation of F1
152
153 //- Switch fields for the production in the k Eqn
157
158 // Switch fields for the differentiation of nut
159 // Holds also for the differentiation of the second branch of
160 // GbyNu
164
165 // Switch fields for GPrime
168
169 // Zero first cell field and IDs
170 // Since the omega equation is a two-zonal one, some
171 // of the terms in the adjoint equations need to ba canceled
172 // at the cells next to omegaWallFunction BCs
175
176
177 // Turbulence model multipliers
178
179 //- Nut Jacobian w.r.t. omega
181
182 //- Nut Jacobian w.r.t. k
184
185 //- Diffusivity of the omega equation
187
188 //- Diffusivity of the k equation
190
191
192 // Protected Member Functions
193
194 // Primal functions
195
196 virtual tmp<volScalarField> F1() const;
197 virtual tmp<volScalarField> F2() const;
199 (
200 const volScalarField& GbyNu0,
201 const volScalarField& F2,
202 const volScalarField& S2
203 ) const;
204
205 //- Return G/nu
207 (
208 const volScalarField::Internal& GbyNu0,
211 ) const;
214 (
215 const volScalarField& F1,
216 const dimensionedScalar& psi1,
218 ) const
219 {
220 return F1*(psi1 - psi2) + psi2;
221 }
224 (
226 const dimensionedScalar& psi1,
228 ) const
229 {
230 return F1*(psi1 - psi2) + psi2;
231 }
234 {
235 return blend(F1, alphaK1_, alphaK2_);
236 }
239 {
241 }
244 (
246 ) const
247 {
249 (
250 this->type() + ":beta",
251 blend(F1, beta1_, beta2_)
252 );
253 }
256 (
257 const volScalarField& F1
258 ) const
259 {
261 (
262 this->type() + ":beta",
263 blend(F1, beta1_, beta2_)
264 );
265 }
268 (
270 ) const
271 {
273 (
274 this->type() + ":gamma",
275 blend(F1, gamma1_, gamma2_)
276 );
277 }
280 (
281 const volScalarField& F1
282 ) const
283 {
285 (
286 this->type() + ":gamma",
287 blend(F1, gamma1_, gamma2_)
288 );
289 }
290
292
293
294 // References to the primal turbulence model variables
296 inline const volScalarField& k() const
297 {
298 return primalVars_.RASModelVariables()().TMVar1();
299 }
301 inline volScalarField& k()
302 {
303 return primalVars_.RASModelVariables()().TMVar1();
304 }
306 inline const volScalarField& omega() const
307 {
308 return primalVars_.RASModelVariables()().TMVar2();
309 }
311 inline volScalarField& omega()
312 {
313 return primalVars_.RASModelVariables()().TMVar2();
314 }
316 inline const volScalarField& nutRef() const
317 {
319 }
321 inline volScalarField& nutRef()
322 {
324 }
325
326
327 // Adjoint related protected member functions
328
329 //- Derivative of the primal equations wrt nut
331
332 //- Nut Jacobian wrt omega
334
335 //- Nut Jacobian wrt k
337
338 //- F2 Jacobian wrt omega
340
341 //- F2 Jacobian wrt k
343
344 //- GbyNu Jacobian wrt omega
346
347 //- GbyNu Jacobian wrt k
349
350 //- Derivative of the primal equations wrt F1
352
353 //- F1 Jacobian wrt omega (no contributions from grad(omega))
355
356 //- F1 Jacobian wrt grad(omega)
358 (
359 const volScalarField& arg1
360 ) const;
361
362 //- Source to waEqn from the differentiation of F1
364
365 //- Source to waEqn from the differentiation of CDkOmega
367
368 //- F1 Jacobian wrt k (no contributions from grad(k))
369 tmp<volScalarField> dF1_dk(const volScalarField& arg1) const;
370
371 //- F1 Jacobian wrt grad(k)
373
374 //- Source to kaEqn from the differentiation of F1
376
377 //- Source to kaEqn from the differentiation of CDkOmega
379
380 //- Differentiation of the turbulence model diffusion coefficients
382 (
383 const volScalarField& primalField,
384 const volScalarField& adjointField,
385 const word& schemeName
386 ) const;
387
388 //- Term multiplying dnut/db, coming from the turbulence model
390 (
391 const volScalarField& primalField,
392 const volScalarField& adjointField,
393 const volScalarField& coeffField,
394 const volScalarField& bcField,
395 const word& schemeName
396 ) const;
397
398 //- Term multiplying dnut/db, coming from the momentum equations
400 (
401 const volVectorField& primalField,
402 const volVectorField& adjointField,
403 const volScalarField& bcField,
404 const word& schemeName
405 ) const;
406
407 // Functions computing the adjoint mean flow source
408
409 //- Contributions from the turbulence model convection terms
411 (
412 const volScalarField& primalField,
413 const volScalarField& adjointField
414 ) const;
415
416 //- Contributions from the G
418 (
419 tmp<volSymmTensorField>& GbyNuMult
420 ) const;
421
422 //- Contributions from the divU
424 (
425 tmp<volScalarField>& divUMult
426 ) const;
427
428 //- Contributions from nut(U), in the diffusion coefficients
429 //- of the turbulence model
431 (
432 const volScalarField& primalField,
433 const volScalarField& adjointField,
434 const volScalarField& coeffField
435 ) const;
436
437 //- Contributions from nut(U)
439 (
441 ) const;
442
443
444 //- Contributions from the differentiation of k existing in
445 //- nutkWallFunction.
446 // This could also be implemented in kaqRWallFunction but all
447 // the fields required for the computation already exist here,
448 // hence the code complexity is reduced
450 (
451 fvScalarMatrix& kaEqn,
453 );
454
455 // References to the adjoint turbulence model fields
457 inline volScalarField& ka()
458 {
459 return adjointTMVariable1Ptr_();
460 };
462 inline const volScalarField& ka() const
463 {
464 return adjointTMVariable1Ptr_();
465 };
467 inline volScalarField& wa()
468 {
469 return adjointTMVariable2Ptr_();
470 };
472 inline const volScalarField& wa() const
473 {
474 return adjointTMVariable2Ptr_();
475 };
476
477
478 //- Update of the primal cached fields
480
481 //- Return the requested interpolation scheme if it exists,
482 //- otherwise return a reverseLinear scheme
483 template<class Type>
485 (
486 const word& schemeName
487 ) const;
488
489 //- Return the interpolation scheme used by the primal convection
490 //- term of the equation corresponding to the argument
492 (
493 const word& varName
494 ) const;
495
496
497public:
498
499 //- Runtime type information
500 TypeName("adjointkOmegaSST");
501
502
503 // Constructors
504
505 //- Construct from components
507 (
508 incompressibleVars& primalVars,
510 objectiveManager& objManager,
511 const word& adjointTurbulenceModelName
513 const word& modelName = typeName
514 );
515
516
517 //- Destructor
518 virtual ~adjointkOmegaSST() = default;
519
520
521 // Member Functions
522
523 //- Return the effective diffusivity for k
525 {
527 (
528 new volScalarField("DkEff", alphaK(F1)*this->nut() + this->nu())
529 );
530 }
531
532 //- Return the effective diffusivity for omega
534 {
536 (
538 (
539 "DomegaEff",
540 alphaOmega(F1)*this->nut() + this->nu()
541 )
542 );
543 }
544
545 virtual tmp<volSymmTensorField> devReff() const;
546
547 virtual tmp<volSymmTensorField> devReff(const volVectorField& U) const;
548
549 //- Return the transpose part of the adjoint momentum stresses
551
552 //- Non-conservative part of the terms added to the mean flow equations
554
555 //- Source term added to the adjoint mean flow due to the
556 //- differentiation of the turbulence model
558
559 //- Jacobian of nut wrt to k
560 // Needs to be implemented for objectives related to nut, defined in
561 // the internal field
563
564 //- Jacobian of nut wrt to omega
565 // Needs to be implemented for objectives related to nut, defined in
566 // the internal field
568
569 //- Diffusion coeff at the boundary for k
570 virtual tmp<scalarField> diffusionCoeffVar1(label patchI) const;
571
572 //- Diffusion coeff at the boundary for omega
573 virtual tmp<scalarField> diffusionCoeffVar2(label patchI) const;
574
575 virtual const boundaryVectorField& adjointMomentumBCSource() const;
576
577 //- Sensitivity derivative contributions when using the (E)SI approach
579
581
582 //- Contributions to the adjoint eikonal equation (zero for now)
584
585 //- Sensitivity derivative contributions when using the FI approach
587
589 (
590 const word& designVarsName
591 ) const;
592
593 //- Nullify all adjoint turbulence model fields and their old times
594 virtual void nullify();
595
596 //- Solve the adjoint turbulence equations
597 virtual void correct();
598
599 //- Read adjointRASProperties dictionary
600 virtual bool read();
601};
602
603
604// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
605
606} // End namespace adjointRASModels
607} // End namespace incompressibleAdjoint
608} // End namespace Foam
609
610// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
611
612#ifdef NoRepository
614#endif
615
616// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
617
618#endif
619
620// ************************************************************************* //
const volScalarField & psi2
const volScalarField & psi1
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Generic GeometricBoundaryField class.
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:78
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:121
Manages the adjoint mean flow fields and their mean values.
Abstract base class for incompressible turbulence models.
autoPtr< volScalarField > adjointTMVariable1Ptr_
Adjoint turbulence model variable 1.
autoPtr< volScalarField > adjointTMVariable2Ptr_
Adjoint turbulence model variable 2.
Continuous adjoint to the kOmegaSST turbulence model for incompressible flows.
virtual tmp< fvVectorMatrix > divDevReff(volVectorField &U) const
Return the transpose part of the adjoint momentum stresses.
virtual tmp< scalarField > diffusionCoeffVar1(label patchI) const
Diffusion coeff at the boundary for k.
virtual tmp< volScalarField > nutJacobianTMVar2() const
Jacobian of nut wrt to omega.
tmp< volScalarField > dGPrime_domega() const
GbyNu Jacobian wrt omega.
tmp< volVectorField > GMeanFlowSource(tmp< volSymmTensorField > &GbyNuMult) const
Contributions from the G.
virtual tmp< volTensorField > FISensitivityTerm()
Sensitivity derivative contributions when using the FI approach.
tmp< volVectorField > dF1_dGradK(const volScalarField &arg1) const
F1 Jacobian wrt grad(k)
tmp< volScalarField::Internal > blend(const volScalarField::Internal &F1, const dimensionedScalar &psi1, const dimensionedScalar &psi2) const
tmp< volScalarField > dR_dF1() const
Derivative of the primal equations wrt F1.
tmp< volScalarField > dF1_dk(const volScalarField &arg1) const
F1 Jacobian wrt k (no contributions from grad(k))
tmp< volScalarField > dnut_dk() const
Nut Jacobian wrt k.
tmp< volScalarField > dGPrime_dk() const
GbyNu Jacobian wrt k.
tmp< volVectorField > convectionMeanFlowSource(const volScalarField &primalField, const volScalarField &adjointField) const
Contributions from the turbulence model convection terms.
tmp< volScalarField::Internal > gamma(const volScalarField::Internal &F1) const
tmp< surfaceInterpolationScheme< Type > > interpolationScheme(const word &schemeName) const
virtual void correct()
Solve the adjoint turbulence equations.
volScalarField DkEff_
Diffusivity of the k equation.
volScalarField S2_
Primal cached fields involved in the solution of the.
tmp< volScalarField > dNutdbMult(const volScalarField &primalField, const volScalarField &adjointField, const volScalarField &coeffField, const volScalarField &bcField, const word &schemeName) const
Term multiplying dnut/db, coming from the turbulence model.
tmp< fvScalarMatrix > waEqnSourceFromCDkOmega() const
Source to waEqn from the differentiation of CDkOmega.
tmp< volScalarField > kaEqnSourceFromF1() const
Source to kaEqn from the differentiation of F1.
tmp< volScalarField > dR_dnut()
Derivative of the primal equations wrt nut.
void updatePrimalRelatedFields()
Update of the primal cached fields.
tmp< volScalarField > dF2_dk() const
F2 Jacobian wrt k.
virtual const boundaryVectorField & adjointMomentumBCSource() const
tmp< surfaceInterpolationScheme< scalar > > convectionScheme(const word &varName) const
virtual tmp< volSymmTensorField > devReff() const
Return the effective stress tensor including the laminar stress.
tmp< volScalarField > beta(const volScalarField &F1) const
virtual void nullify()
Nullify all adjoint turbulence model fields and their old times.
tmp< volScalarField > blend(const volScalarField &F1, const dimensionedScalar &psi1, const dimensionedScalar &psi2) const
virtual tmp< volScalarField > nutJacobianTMVar1() const
Jacobian of nut wrt to k.
tmp< volVectorField > nutMeanFlowSource(tmp< volScalarField > &mult) const
Contributions from nut(U)
volTensorField gradU_
Cached primal gradient fields.
volScalarField DOmegaEff_
Diffusivity of the omega equation.
virtual tmp< volScalarField > GbyNu(const volScalarField &GbyNu0, const volScalarField &F2, const volScalarField &S2) const
tmp< volScalarField > alphaK(const volScalarField &F1) const
tmp< volScalarField > diffusionNutMeanFlowMult(const volScalarField &primalField, const volScalarField &adjointField, const volScalarField &coeffField) const
virtual const boundaryVectorField & wallShapeSensitivities()
Sensitivity derivative contributions when using the (E)SI approach.
tmp< volScalarField > DomegaEff(const volScalarField &F1) const
Return the effective diffusivity for omega.
tmp< volScalarField > alphaOmega(const volScalarField &F1) const
void addWallFunctionTerms(fvScalarMatrix &kaEqn, const volScalarField &dR_dnut)
tmp< volScalarField > waEqnSourceFromF1() const
Source to waEqn from the differentiation of F1.
volScalarField case_1_Pk_
Switch fields for the production in the k Eqn.
tmp< volScalarField::Internal > beta(const volScalarField::Internal &F1) const
virtual const boundaryVectorField & wallFloCoSensitivities()
Sensitivity terms for flow control, emerging from the.
virtual tmp< volVectorField > nonConservativeMomentumSource() const
Non-conservative part of the terms added to the mean flow equations.
TypeName("adjointkOmegaSST")
Runtime type information.
tmp< volScalarField > kaEqnSourceFromCDkOmega() const
Source to kaEqn from the differentiation of CDkOmega.
tmp< volScalarField > coeffsDifferentiation(const volScalarField &primalField, const volScalarField &adjointField, const word &schemeName) const
Differentiation of the turbulence model diffusion coefficients.
virtual tmp< volScalarField > distanceSensitivities()
Contributions to the adjoint eikonal equation (zero for now)
virtual tmp< scalarField > diffusionCoeffVar2(label patchI) const
Diffusion coeff at the boundary for omega.
tmp< volVectorField > divUMeanFlowSource(tmp< volScalarField > &divUMult) const
Contributions from the divU.
tmp< volVectorField > dF1_dGradOmega(const volScalarField &arg1) const
F1 Jacobian wrt grad(omega)
virtual tmp< scalarField > topologySensitivities(const word &designVarsName) const
tmp< volScalarField > dnut_domega() const
Nut Jacobian wrt omega.
tmp< volScalarField > DkEff(const volScalarField &F1) const
Return the effective diffusivity for k.
tmp< volScalarField > dF1_domega(const volScalarField &arg1) const
F1 Jacobian wrt omega (no contributions from grad(omega))
virtual bool read()
Read adjointRASProperties dictionary.
tmp< volScalarField > gamma(const volScalarField &F1) const
tmp< volScalarField > dF2_domega() const
F2 Jacobian wrt omega.
tmp< volScalarField > nu() const
Return the laminar viscosity.
virtual const volScalarField & nut() const
Return the turbulence viscosity.
Base class for solution control classes.
const autoPtr< incompressible::RASModelVariables > & RASModelVariables() const
Return const reference to the turbulence model variables.
class for managing incompressible objective functions.
A class for managing temporary objects.
Definition: tmp.H:65
A class for handling words, derived from Foam::string.
Definition: word.H:68
U
Definition: pEqn.H:72
const scalar gamma
Definition: EEqn.H:9
Namespace for OpenFOAM.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:598
static const char *const typeName
The type name used in ensight case files.
#define TypeName(TypeNameString)
Declare a ClassName() with extra virtual type info.
Definition: typeInfo.H:73