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_thermal_face.h
Go to the documentation of this file.
1 // KRATOS ___ ___ _ ___ __ ___ ___ ___ ___
2 // / __/ _ \| \| \ \ / /__| \_ _| __| __|
3 // | (_| (_) | .` |\ V /___| |) | || _|| _|
4 // \___\___/|_|\_| \_/ |___/___|_| |_| APPLICATION
5 //
6 // License: BSD License
7 // Kratos default license: kratos/license.txt
8 //
9 // Main authors: Jordi Cotela
10 //
11 
12 #if !defined(KRATOS_ADJOINT_THERMAL_FACE_H_INCLUDED )
13 #define KRATOS_ADJOINT_THERMAL_FACE_H_INCLUDED
14 
15 // System includes
16 
17 // External includes
18 
19 // Project includes
21 
22 namespace Kratos
23 {
24 
27 
30 
34 
38 
42 
46 
49 {
50 public:
53 
56 
64 
65  constexpr static double StefanBoltzmann = ThermalFace::StefanBoltzmann;
66 
70 
71  AdjointThermalFace(IndexType NewId, typename GeometryType::Pointer pGeometry);
72 
73  AdjointThermalFace(IndexType NewId, typename GeometryType::Pointer pGeometry, Properties::Pointer pProperties);
74 
76  ~AdjointThermalFace() override;
77 
81 
82  Condition::Pointer Create(IndexType NewId, NodesArrayType const& ThisNodes, Properties::Pointer pProperties) const override;
83 
84  Condition::Pointer Create(IndexType NewId, typename GeometryType::Pointer pGeom, Properties::Pointer pProperties) const override;
85 
86  void CalculateLocalSystem(MatrixType& rLeftHandSideMatrix, VectorType& rRightHandSideVector, const ProcessInfo& rCurrentProcessInfo) override;
87 
88  void CalculateRightHandSide(VectorType& rRightHandSideVector, const ProcessInfo& rCurrentProcessInfo) override;
89 
90  void GetValuesVector(Vector& rValues, int Step = 0) const override;
91 
93  const ProcessInfo& rCurrentProcessInfo) const override;
94 
95  void GetDofList(DofsVectorType& rConditionDofList, const ProcessInfo& rCurrentProcessInfo) const override;
96 
97  void CalculateSensitivityMatrix(const Variable<array_1d<double, 3>>& rDesignVariable,
98  Matrix& rOutput,
99  const ProcessInfo& rCurrentProcessInfo) override;
100 
104 
105  int Check(const ProcessInfo &rCurrentProcessInfo) const override;
106 
110 
112  std::string Info() const override;
113 
115  void PrintInfo(std::ostream& rOStream) const override;
116 
118 
119 protected:
122 
123 
127 
128 
132 
133 
137 
138 
142 
143 
147 
148 
152 
153 
155 
156 private:
159 
160 
164 
165 
169 
170  friend class Serializer;
171 
172  // A private default constructor necessary for serialization
174  {
175  }
176 
177  void save(Serializer& rSerializer) const override
178  {
179  using BaseType = ThermalFace;
181  }
182 
183  void load(Serializer& rSerializer) override
184  {
185  using BaseType = ThermalFace;
187  }
188 
192 
193 
197 
198  // Note: the Jacobian should come from the geometry, but non-square ones are often wrongly implemented at the moment :(
199  MatrixType GetJacobian(GeometryData::IntegrationMethod QuadratureOrder, unsigned int IntegrationPointIndex) const;
200 
204 
205 
209 
210 
214 
216  AdjointThermalFace& operator=(const AdjointThermalFace& rOther) = delete;
217 
219  AdjointThermalFace(const AdjointThermalFace& rOther) = delete;
220 
222 
223 }; // Class AdjointThermalFace
224 
226 
229 
230 
234 
235 
237 inline std::istream& operator >> (std::istream& rIStream, AdjointThermalFace& rThis)
238 {
239  return rIStream;
240 }
241 
243 inline std::ostream& operator << (std::ostream& rOStream, const AdjointThermalFace& rThis)
244 {
245  rThis.PrintInfo(rOStream);
246  rOStream << std::endl;
247  rThis.PrintData(rOStream);
248 
249  return rOStream;
250 }
251 
253 
255 
256 } // namespace Kratos.
257 
258 #endif // KRATOS_ADJOINT_THERMAL_FACE_H_INCLUDED defined
259 
260 
PeriodicInterfaceProcess & operator=(const PeriodicInterfaceProcess &)=delete
Heat flux Neumann condition for the ajdoint thermal diffusion problem.
Definition: adjoint_thermal_face.h:49
void GetDofList(DofsVectorType &rConditionDofList, const ProcessInfo &rCurrentProcessInfo) const override
Definition: adjoint_thermal_face.cpp:112
typename ThermalFace::MatrixType MatrixType
Definition: adjoint_thermal_face.h:60
void EquationIdVector(EquationIdVectorType &rResult, const ProcessInfo &rCurrentProcessInfo) const override
Definition: adjoint_thermal_face.cpp:94
constexpr static double StefanBoltzmann
Definition: adjoint_thermal_face.h:65
AdjointThermalFace(IndexType NewId, typename GeometryType::Pointer pGeometry)
Definition: adjoint_thermal_face.cpp:24
void CalculateSensitivityMatrix(const Variable< array_1d< double, 3 >> &rDesignVariable, Matrix &rOutput, const ProcessInfo &rCurrentProcessInfo) override
Definition: adjoint_thermal_face.cpp:169
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(AdjointThermalFace)
Counted pointer of AdjointThermalFace.
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: adjoint_thermal_face.cpp:161
void GetValuesVector(Vector &rValues, int Step=0) const override
Definition: adjoint_thermal_face.cpp:78
~AdjointThermalFace() override
Destructor.
Definition: adjoint_thermal_face.cpp:33
std::string Info() const override
Turn back information as a string.
Definition: adjoint_thermal_face.cpp:154
Condition::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, Properties::Pointer pProperties) const override
Definition: adjoint_thermal_face.cpp:35
friend class Serializer
Definition: adjoint_thermal_face.h:170
int Check(const ProcessInfo &rCurrentProcessInfo) const override
Definition: adjoint_thermal_face.cpp:129
void CalculateRightHandSide(VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: adjoint_thermal_face.cpp:63
void CalculateLocalSystem(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: adjoint_thermal_face.cpp:51
GeometricalObject BaseType
base type: an GeometricalObject that automatically has a unique number
Definition: condition.h:71
std::size_t IndexType
Definition: condition.h:92
Geometry< NodeType >::PointsArrayType NodesArrayType
definition of nodes container type, redefined from GeometryType
Definition: condition.h:86
std::vector< std::size_t > EquationIdVectorType
Definition: condition.h:98
std::vector< DofType::Pointer > DofsVectorType
Definition: condition.h:100
Geometry< NodeType > GeometryType
definition of the geometry type with given NodeType
Definition: condition.h:83
std::size_t IndexType
Definition: flags.h:74
IntegrationMethod
Definition: geometry_data.h:76
Geometry base class.
Definition: geometry.h:71
PointerVector is a container like stl vector but using a vector to store pointers to its data.
Definition: pointer_vector.h:72
ProcessInfo holds the current value of different solution parameters.
Definition: process_info.h:59
The serialization consists in storing the state of an object into a storage format like data file or ...
Definition: serializer.h:123
A basic Neumann condition for convection-diffusion problems.
Definition: thermal_face.h:44
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: thermal_face.cpp:311
ThermalFace()
Definition: thermal_face.cpp:390
Condition::MatrixType MatrixType
Definition: thermal_face.h:88
constexpr static double StefanBoltzmann
Stefan Boltzmann constant for radiation in SI units: [W / (m^2 K^4)].
Definition: thermal_face.h:92
Condition::VectorType VectorType
Definition: thermal_face.h:89
Variable class contains all information needed to store and retrive data from a data container.
Definition: variable.h:63
#define KRATOS_SERIALIZE_SAVE_BASE_CLASS(Serializer, BaseType)
Definition: define.h:812
#define KRATOS_SERIALIZE_LOAD_BASE_CLASS(Serializer, BaseType)
Definition: define.h:815
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
std::istream & operator>>(std::istream &rIStream, LinearMasterSlaveConstraint &rThis)
input stream function
std::ostream & operator<<(std::ostream &rOStream, const LinearMasterSlaveConstraint &rThis)
output stream function
Definition: linear_master_slave_constraint.h:432
def load(f)
Definition: ode_solve.py:307