KratosMultiphysics
KRATOS Multiphysics (Kratos) is a framework for building parallel, multi-disciplinary simulation software, aiming at modularity, extensibility, and high performance. Kratos is written in C++, and counts with an extensive Python interface.
adjoint_response_function.h
Go to the documentation of this file.
1 // | / |
2 // ' / __| _` | __| _ \ __|
3 // . \ | ( | | ( |\__ `
4 // _|\_\_| \__,_|\__|\___/ ____/
5 // Multi-Physics
6 //
7 // License: BSD License
8 // Kratos default license: kratos/license.txt
9 //
10 // Main authors:
11 //
12 
13 #ifndef ADJOINT_RESPONSE_FUNCTION_H
14 #define ADJOINT_RESPONSE_FUNCTION_H
15 
16 // System includes
17 
18 // External includes
19 
20 // Project includes
21 #include "includes/define.h"
22 #include "includes/element.h"
23 #include "includes/condition.h"
24 #include "includes/process_info.h"
26 #include "containers/array_1d.h"
27 #include "containers/variable.h"
28 
29 
30 // Application includes
31 
32 namespace Kratos
33 {
36 
39 {
40 public:
43 
45 
49 
52  {
53  }
54 
57  {
58  }
59 
63 
67 
68  virtual void Initialize()
69  {
70  }
71 
72  virtual void InitializeSolutionStep()
73  {
74  }
75 
76  virtual void FinalizeSolutionStep()
77  {
78  }
79 
81 
88  virtual void CalculateGradient(const Element& rAdjointElement,
89  const Matrix& rResidualGradient,
90  Vector& rResponseGradient,
91  const ProcessInfo& rProcessInfo)
92  {
93  KRATOS_ERROR << "Calling base class response function method.\n";
94  }
95 
97 
104  virtual void CalculateGradient(const Condition& rAdjointCondition,
105  const Matrix& rResidualGradient,
106  Vector& rResponseGradient,
107  const ProcessInfo& rProcessInfo)
108  {
109  KRATOS_ERROR << "Calling base class response function method.\n";
110  }
111 
113 
120  virtual void CalculateFirstDerivativesGradient(const Element& rAdjointElement,
121  const Matrix& rResidualGradient,
122  Vector& rResponseGradient,
123  const ProcessInfo& rProcessInfo)
124  {
125  KRATOS_ERROR << "Calling base class response function method.\n";
126  }
127 
129 
136  virtual void CalculateFirstDerivativesGradient(const Condition& rAdjointCondition,
137  const Matrix& rResidualGradient,
138  Vector& rResponseGradient,
139  const ProcessInfo& rProcessInfo)
140  {
141  KRATOS_ERROR << "Calling base class response function method.\n";
142  }
143 
145 
152  virtual void CalculateSecondDerivativesGradient(const Element& rAdjointElement,
153  const Matrix& rResidualGradient,
154  Vector& rResponseGradient,
155  const ProcessInfo& rProcessInfo)
156  {
157  KRATOS_ERROR << "Calling base class response function method.\n";
158  }
159 
161 
168  virtual void CalculateSecondDerivativesGradient(const Condition& rAdjointCondition,
169  const Matrix& rResidualGradient,
170  Vector& rResponseGradient,
171  const ProcessInfo& rProcessInfo)
172  {
173  KRATOS_ERROR << "Calling base class response function method.\n";
174  }
175 
177 
185  virtual void CalculatePartialSensitivity(Element& rAdjointElement,
186  const Variable<double>& rVariable,
187  const Matrix& rSensitivityMatrix,
188  Vector& rSensitivityGradient,
189  const ProcessInfo& rProcessInfo)
190  {
191  KRATOS_ERROR << "Calling base class response function method.\n";
192  }
193 
195 
203  virtual void CalculatePartialSensitivity(Condition& rAdjointCondition,
204  const Variable<double>& rVariable,
205  const Matrix& rSensitivityMatrix,
206  Vector& rSensitivityGradient,
207  const ProcessInfo& rProcessInfo)
208  {
209  KRATOS_ERROR << "Calling base class response function method.\n";
210  }
211 
213 
221  virtual void CalculatePartialSensitivity(Element& rAdjointElement,
222  const Variable<array_1d<double, 3>>& rVariable,
223  const Matrix& rSensitivityMatrix,
224  Vector& rSensitivityGradient,
225  const ProcessInfo& rProcessInfo)
226  {
227  KRATOS_ERROR << "Calling base class response function method.\n";
228  }
229 
231 
239  virtual void CalculatePartialSensitivity(Condition& rAdjointCondition,
240  const Variable<array_1d<double, 3>>& rVariable,
241  const Matrix& rSensitivityMatrix,
242  Vector& rSensitivityGradient,
243  const ProcessInfo& rProcessInfo)
244  {
245  KRATOS_ERROR << "Calling base class response function method.\n";
246  }
247 
249  virtual double CalculateValue(ModelPart& rModelPart) = 0;
250 
252 
253 protected:
256 
260 
264 
266 
267 private:
270 
274 
278 
280 };
281 
283 
284 } /* namespace Kratos.*/
285 
286 #endif /* ADJOINT_RESPONSE_FUNCTION_H defined */
A base class for adjoint response functions.
Definition: adjoint_response_function.h:39
virtual void CalculateFirstDerivativesGradient(const Condition &rAdjointCondition, const Matrix &rResidualGradient, Vector &rResponseGradient, const ProcessInfo &rProcessInfo)
Calculate the local gradient w.r.t. first derivatives of primal solution.
Definition: adjoint_response_function.h:136
virtual void CalculatePartialSensitivity(Element &rAdjointElement, const Variable< double > &rVariable, const Matrix &rSensitivityMatrix, Vector &rSensitivityGradient, const ProcessInfo &rProcessInfo)
Calculate the partial sensitivity w.r.t. design variable.
Definition: adjoint_response_function.h:185
virtual void Initialize()
Definition: adjoint_response_function.h:68
KRATOS_CLASS_POINTER_DEFINITION(AdjointResponseFunction)
virtual void InitializeSolutionStep()
Definition: adjoint_response_function.h:72
virtual void CalculatePartialSensitivity(Condition &rAdjointCondition, const Variable< array_1d< double, 3 >> &rVariable, const Matrix &rSensitivityMatrix, Vector &rSensitivityGradient, const ProcessInfo &rProcessInfo)
Calculate the partial sensitivity w.r.t. design variable.
Definition: adjoint_response_function.h:239
virtual double CalculateValue(ModelPart &rModelPart)=0
Calculate the scalar valued response function.
virtual void CalculatePartialSensitivity(Condition &rAdjointCondition, const Variable< double > &rVariable, const Matrix &rSensitivityMatrix, Vector &rSensitivityGradient, const ProcessInfo &rProcessInfo)
Calculate the partial sensitivity w.r.t. design variable.
Definition: adjoint_response_function.h:203
virtual void CalculatePartialSensitivity(Element &rAdjointElement, const Variable< array_1d< double, 3 >> &rVariable, const Matrix &rSensitivityMatrix, Vector &rSensitivityGradient, const ProcessInfo &rProcessInfo)
Calculate the partial sensitivity w.r.t. design variable.
Definition: adjoint_response_function.h:221
virtual void CalculateGradient(const Condition &rAdjointCondition, const Matrix &rResidualGradient, Vector &rResponseGradient, const ProcessInfo &rProcessInfo)
Calculate the local gradient w.r.t. primal solution.
Definition: adjoint_response_function.h:104
virtual void CalculateFirstDerivativesGradient(const Element &rAdjointElement, const Matrix &rResidualGradient, Vector &rResponseGradient, const ProcessInfo &rProcessInfo)
Calculate the local gradient w.r.t. first derivatives of primal solution.
Definition: adjoint_response_function.h:120
AdjointResponseFunction()
Constructor.
Definition: adjoint_response_function.h:51
virtual void CalculateSecondDerivativesGradient(const Element &rAdjointElement, const Matrix &rResidualGradient, Vector &rResponseGradient, const ProcessInfo &rProcessInfo)
Calculate the local gradient w.r.t. second derivatives of primal solution.
Definition: adjoint_response_function.h:152
virtual void CalculateGradient(const Element &rAdjointElement, const Matrix &rResidualGradient, Vector &rResponseGradient, const ProcessInfo &rProcessInfo)
Calculate the local gradient w.r.t. primal solution.
Definition: adjoint_response_function.h:88
virtual void CalculateSecondDerivativesGradient(const Condition &rAdjointCondition, const Matrix &rResidualGradient, Vector &rResponseGradient, const ProcessInfo &rProcessInfo)
Calculate the local gradient w.r.t. second derivatives of primal solution.
Definition: adjoint_response_function.h:168
virtual void FinalizeSolutionStep()
Definition: adjoint_response_function.h:76
virtual ~AdjointResponseFunction()
Destructor.
Definition: adjoint_response_function.h:56
Base class for all Conditions.
Definition: condition.h:59
Base class for all Elements.
Definition: element.h:60
This class aims to manage meshes for multi-physics simulations.
Definition: model_part.h:77
ProcessInfo holds the current value of different solution parameters.
Definition: process_info.h:59
#define KRATOS_ERROR
Definition: exception.h:161
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21