DAObjFuncVonMisesStressKS.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2 
3  DAFoam : Discrete Adjoint with OpenFOAM
4  Version : v3
5 
6 \*---------------------------------------------------------------------------*/
7 
9 
10 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
11 
12 namespace Foam
13 {
14 
15 defineTypeNameAndDebug(DAObjFuncVonMisesStressKS, 0);
16 addToRunTimeSelectionTable(DAObjFunc, DAObjFuncVonMisesStressKS, dictionary);
17 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
18 
20  const fvMesh& mesh,
21  const DAOption& daOption,
22  const DAModel& daModel,
23  const DAIndex& daIndex,
24  const DAResidual& daResidual,
25  const word objFuncName,
26  const word objFuncPart,
27  const dictionary& objFuncDict)
28  : DAObjFunc(
29  mesh,
30  daOption,
31  daModel,
32  daIndex,
33  daResidual,
34  objFuncName,
35  objFuncPart,
36  objFuncDict)
37 {
38  // Assign type, this is common for all objectives
39  objFuncDict_.readEntry<word>("type", objFuncType_);
40 
41  // setup the connectivity for heat flux, this is needed in Foam::DAJacCondFdW
42  objFuncConInfo_ = {
43  {"D"}, // level 0
44  {"D"}}; // level 1
45 
46  objFuncDict_.readEntry<scalar>("scale", scale_);
47 
48  objFuncDict_.readEntry<scalar>("coeffKS", coeffKS_);
49 }
50 
53  const labelList& objFuncFaceSources,
54  const labelList& objFuncCellSources,
55  scalarList& objFuncFaceValues,
56  scalarList& objFuncCellValues,
57  scalar& objFuncValue)
58 {
59  /*
60  Description:
61  Calculate the maximal von Mises stress aggregated using the KS function
62  von Mises stress = KS( mu*twoSymm(fvc::grad(D)) + lambda*(I*tr(fvc::grad(D))) )
63  where the KS function KS(x) = 1/coeffKS * ln( sum[exp(coeffKS*x_i)] )
64 
65  Input:
66  objFuncFaceSources: List of face source (index) for this objective
67 
68  objFuncCellSources: List of cell source (index) for this objective
69 
70  Output:
71  objFuncFaceValues: the discrete value of objective for each face source (index).
72  This will be used for computing df/dw in the adjoint.
73 
74  objFuncCellValues: the discrete value of objective on each cell source (index).
75  This will be used for computing df/dw in the adjoint.
76 
77  objFuncValue: the sum of objective, reduced across all processors and scaled by "scale"
78  */
79 
80  // initialize faceValues to zero
81  forAll(objFuncCellValues, idxI)
82  {
83  objFuncCellValues[idxI] = 0.0;
84  }
85 
86  objFuncValue = 0.0;
87 
88  const objectRegistry& db = mesh_.thisDb();
89  //const volVectorField& D = db.lookupObject<volVectorField>("D");
90  const volScalarField& lambda = db.lookupObject<volScalarField>("solid:lambda");
91  const volScalarField& mu = db.lookupObject<volScalarField>("solid:mu");
92  const volScalarField& rho = db.lookupObject<volScalarField>("solid:rho");
93  const volTensorField& gradD = db.lookupObject<volTensorField>("gradD");
94 
95  volSymmTensorField sigma = rho * (mu * twoSymm(gradD) + lambda * (I * tr(gradD)));
96 
97  // NOTE: vonMises stress is scaled by scale_ provided in the objFunc dict
98  volScalarField vonMises = scale_* sqrt((3.0 / 2.0) * magSqr(dev(sigma)));
99 
100  scalar objValTmp = 0.0;
101  forAll(objFuncCellSources, idxI)
102  {
103  const label& cellI = objFuncCellSources[idxI];
104 
105  objFuncCellValues[idxI] = exp(coeffKS_ * vonMises[cellI]);
106 
107  objValTmp += objFuncCellValues[idxI];
108 
109  if (objValTmp > 1e200)
110  {
111  FatalErrorIn(" ") << "KS function summation term too large! "
112  << "Reduce coeffKS! " << abort(FatalError);
113  }
114  }
115 
116  // need to reduce the sum of force across all processors
117  reduce(objValTmp, sumOp<scalar>());
118 
119  // expSumKS stores sum[exp(coeffKS*x_i)], it will be used to scale dFdW
120  expSumKS = objValTmp;
121 
122  objFuncValue = log(objValTmp) / coeffKS_;
123 
124  return;
125 }
126 
127 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
128 
129 } // End namespace Foam
130 
131 // ************************************************************************* //
mu
volScalarField & mu
Definition: createRefsSolidDisplacement.H:6
Foam::DAObjFunc::objFuncType_
word objFuncType_
the type of the objective function
Definition: DAObjFunc.H:66
forAll
forAll(pseudoP.boundaryField(), patchI)
Definition: solvePseudoPEqn.H:10
Foam::DAOption
Definition: DAOption.H:29
daOption
DAOption daOption(mesh, pyOptions_)
Foam::DAObjFuncVonMisesStressKS::calcObjFunc
virtual void calcObjFunc(const labelList &objFuncFaceSources, const labelList &objFuncCellSources, scalarList &objFuncFaceValues, scalarList &objFuncCellValues, scalar &objFuncValue)
calculate the value of objective function
Definition: DAObjFuncVonMisesStressKS.C:52
mesh
fvMesh & mesh
Definition: createRefsHeatTransfer.H:4
Foam::DAObjFunc::mesh_
const fvMesh & mesh_
fvMesh
Definition: DAObjFunc.H:45
Foam::DAIndex
Definition: DAIndex.H:32
Foam::DAModel
Definition: DAModel.H:59
sigma
volSymmTensorField sigma(IOobject("sigma", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), rho *sigmaD)
Foam
Definition: multiFreqScalarFvPatchField.C:144
Foam::DAResidual
Definition: DAResidual.H:35
DAObjFuncVonMisesStressKS.H
Foam::defineTypeNameAndDebug
defineTypeNameAndDebug(DAFvSource, 0)
gradD
volTensorField & gradD
Definition: createRefsSolidDisplacement.H:10
Foam::addToRunTimeSelectionTable
addToRunTimeSelectionTable(DAFvSource, DAFvSourceActuatorDisk, dictionary)
lambda
volScalarField & lambda
Definition: createRefsSolidDisplacement.H:7
Foam::DAObjFunc::scale_
scalar scale_
scale of the objective function
Definition: DAObjFunc.H:87
Foam::DAObjFuncVonMisesStressKS::coeffKS_
scalar coeffKS_
coefficient for the KS function
Definition: DAObjFuncVonMisesStressKS.H:33
Foam::DAObjFuncVonMisesStressKS::DAObjFuncVonMisesStressKS
DAObjFuncVonMisesStressKS(const fvMesh &mesh, const DAOption &daOption, const DAModel &daModel, const DAIndex &daIndex, const DAResidual &daResidual, const word objFuncName, const word objFuncPart, const dictionary &objFuncDict)
Definition: DAObjFuncVonMisesStressKS.C:19
daModel
DAModel daModel(mesh, daOption)
Foam::DAObjFunc::objFuncConInfo_
List< List< word > > objFuncConInfo_
the connectivity information for the objective function
Definition: DAObjFunc.H:93
daIndex
DAIndex daIndex(mesh, daOption, daModel)
Foam::DAObjFunc::objFuncDict_
const dictionary & objFuncDict_
dictionary containing the information for the objective function
Definition: DAObjFunc.H:69
Foam::DAObjFunc
Definition: DAObjFunc.H:33
rho
volScalarField & rho
Definition: createRefsRhoSimpleC.H:8
Foam::DAObjFunc::expSumKS
scalar expSumKS
expSumKS stores sum[exp(coeffKS*x_i)] for KS function which will be used to scale dFdW
Definition: DAObjFunc.H:235