v0.14.0
OperatorsTester.hpp
Go to the documentation of this file.
1 /**
2  * @file OperatorsTester.hpp
3  * @brief Used to calculate derivatives and testing tangent matrix, i.e. hessian.
4  * @version 0.13.2
5  * @date 2022-12-04
6  *
7  * @copyright Copyright (c) 2022
8  *
9  */
10 
11 #ifndef __OPERATORS_TESTER_HPP__
12 #define __OPERATORS_TESTER_HPP__
13 
14 namespace MoFEM {
15 
16 /**
17  * @brief Calculate directional derivative of the right hand side and compare it
18  * with tangent matrix derivative.
19  *
20  */
22 
23  MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
24  UnknownInterface **iface) const;
25 
26  OperatorsTester(const MoFEM::Core &core);
27  virtual ~OperatorsTester() = default;
28 
29  using RandomFieldData = std::pair<std::string, std::array<double, 2>>;
30 
31  /**
32  * @brief Generate random fileds
33  *
34  * Example: generate random vector for DM (problem) from simple interface,
35  where FIELD random values of DOFs are in range from -1 to 1, and FIELD2
36  random values are in range from 0 to 1.
37  * @code {.cpp}
38  * auto x = opt->setRandomFields(simple->getDM(),
39  {{"FIELD1", {-1, 1}}, {"FIELD2", {0,1}}});
40  * @endcode
41  *
42  *
43  * TODO: Set random field to specific entities, and potentially order for
44  testing proposes to dissect error in tangent matrix.
45  *
46  * @param dm
47  * @param random_fields look at definition @ref RandomFieldData
48  * @return SmartPetscObj<Vec> smart vector
49  */
51  std::vector<RandomFieldData> random_fields,
52  boost::shared_ptr<Range> ents = nullptr);
53 
54  /**
55  * @brief Assemble the right hand side vector
56  *
57  * @param dm
58  * @param fe_name // fe name
59  * @param pipeline // pipeline, i.e. fe instance
60  * @param x // problem (dm) vector
61  * @param delta_x // vector for x rate, can be null, i.e, SmartPetscObj<Vec>()
62  * @param delta2_x // vector for x second rate, i.e. acceleration
63  * @param time // time
64  * @param delta_t // time increment
65  * @param cache_ptr // finite element data cache, can be null
66  * @return SmartPetscObj<Vec>
67  */
68  SmartPetscObj<Vec> assembleVec(SmartPetscObj<DM> dm, std::string fe_name,
69  boost::shared_ptr<FEMethod> pipeline,
71  SmartPetscObj<Vec> delta_x,
72  SmartPetscObj<Vec> delta2_x, double time,
73  double delta_t, CacheTupleWeakPtr cache_ptr);
74 
75  /**
76  * @brief Assemble the left hand side vector
77  *
78  * @param dm
79  * @param fe_name // fe name
80  * @param pipeline // pipeline, i.e. fe instance
81  * @param x // problem (dm) vector
82  * @param delta_x // vector for x rate, can be null, i.e, SmartPetscObj<Vec>()
83  * @param delta2_x // vector for x second rate, i.e. acceleration
84  * @param time // time
85  * @param delta_t // time increment
86  * @param cache_ptr // finite element data cache, can be null
87  * @return SmartPetscObj<Vec>
88  */
89  SmartPetscObj<Mat> assembleMat(SmartPetscObj<DM> dm, std::string fe_name,
90  boost::shared_ptr<FEMethod> pipeline,
92  SmartPetscObj<Vec> delta_x,
93  SmartPetscObj<Vec> delta2_x, double time,
94  double delta_t, CacheTupleWeakPtr cache_ptr);
95 
96  /**
97  * @brief Calculate directional directive using finite difference
98  *
99  * @param dm
100  * @param fe_name
101  * @param pipeline
102  * @param x
103  * @param delta_x
104  * @param delta2_x
105  * @param diff_x // direction of derivative
106  * @param time
107  * @param delta_t
108  * @param eps
109  * @param cache_ptr
110  * @return SmartPetscObj<Vec>
111  */
113  SmartPetscObj<DM> dm, std::string fe_name,
114  boost::shared_ptr<FEMethod> pipeline, SmartPetscObj<Vec> x,
115  SmartPetscObj<Vec> delta_x, SmartPetscObj<Vec> delta2_x,
116  SmartPetscObj<Vec> diff_x, double time, double delta_t, double eps,
117  CacheTupleWeakPtr cache_ptr = CacheTupleSharedPtr());
118 
119  /**
120  * @brief Check consistency between directional derivative with matrix
121  *
122  * @param dm
123  * @param fe_name
124  * @param pipeline_rhs
125  * @param pipeline_lhs
126  * @param x
127  * @param delta_x
128  * @param delta2_x
129  * @param diff_x
130  * @param time
131  * @param delta_t
132  * @param eps
133  * @param cache_ptr
134  * @return SmartPetscObj<Vec>
135  */
137  SmartPetscObj<DM> dm, std::string fe_name,
138  boost::shared_ptr<FEMethod> pipeline_rhs,
139  boost::shared_ptr<FEMethod> pipeline_lhs, SmartPetscObj<Vec> x,
140  SmartPetscObj<Vec> delta_x, SmartPetscObj<Vec> delta2_x,
141  SmartPetscObj<Vec> diff_x, double time, double delta_t, double eps,
142  CacheTupleWeakPtr cache_ptr = CacheTupleSharedPtr());
143 
144 private:
146 
147  /**
148  * @brief Set vectors, x, x_t, and x_tt to finite element instance
149  *
150  * Finite element instance is a pipeline. x_t and x_tt are evaluated for given
151  * delta_x, delta2_x and delta_t.
152  *
153  * @param pipeline
154  * @param x
155  * @param delta_x
156  * @param delta2_x
157  * @param delta_t
158  * @return std::pair<SmartPetscObj<Vec>, SmartPetscObj<Vec>>
159  */
160  std::pair<SmartPetscObj<Vec>, SmartPetscObj<Vec>>
161  setPipelineX(boost::shared_ptr<FEMethod> pipeline, SmartPetscObj<Vec> x,
162  SmartPetscObj<Vec> delta_x, SmartPetscObj<Vec> delta2_x,
163  double delta_t);
164 };
165 
166 } // namespace MoFEM
167 
168 #endif //__OPERATORS_TESTER_HPP__
MoFEM::OperatorsTester::OperatorsTester
OperatorsTester(const MoFEM::Core &core)
Definition: OperatorsTester.cpp:13
MoFEM::CoreTmp< 0 >
Core (interface) class.
Definition: Core.hpp:82
MoFEM::OperatorsTester::~OperatorsTester
virtual ~OperatorsTester()=default
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::OperatorsTester::query_interface
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Definition: OperatorsTester.cpp:17
MoFEM::OperatorsTester::setPipelineX
std::pair< SmartPetscObj< Vec >, SmartPetscObj< Vec > > setPipelineX(boost::shared_ptr< FEMethod > pipeline, SmartPetscObj< Vec > x, SmartPetscObj< Vec > delta_x, SmartPetscObj< Vec > delta2_x, double delta_t)
Set vectors, x, x_t, and x_tt to finite element instance.
Definition: OperatorsTester.cpp:233
MoFEM::OperatorsTester::assembleMat
SmartPetscObj< Mat > assembleMat(SmartPetscObj< DM > dm, std::string fe_name, boost::shared_ptr< FEMethod > pipeline, SmartPetscObj< Vec > x, SmartPetscObj< Vec > delta_x, SmartPetscObj< Vec > delta2_x, double time, double delta_t, CacheTupleWeakPtr cache_ptr)
Assemble the left hand side vector.
Definition: OperatorsTester.cpp:120
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
MoFEM::OperatorsTester::cOre
MoFEM::Core & cOre
Definition: OperatorsTester.hpp:145
MoFEM::OperatorsTester::directionalCentralFiniteDifference
SmartPetscObj< Vec > directionalCentralFiniteDifference(SmartPetscObj< DM > dm, std::string fe_name, boost::shared_ptr< FEMethod > pipeline, SmartPetscObj< Vec > x, SmartPetscObj< Vec > delta_x, SmartPetscObj< Vec > delta2_x, SmartPetscObj< Vec > diff_x, double time, double delta_t, double eps, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Calculate directional directive using finite difference.
Definition: OperatorsTester.cpp:153
MoFEM::CacheTupleSharedPtr
boost::shared_ptr< CacheTuple > CacheTupleSharedPtr
Definition: FEMultiIndices.hpp:495
MoFEM::OperatorsTester::setRandomFields
SmartPetscObj< Vec > setRandomFields(SmartPetscObj< DM > dm, std::vector< RandomFieldData > random_fields, boost::shared_ptr< Range > ents=nullptr)
Generate random fileds.
Definition: OperatorsTester.cpp:24
MoFEM::OperatorsTester
Calculate directional derivative of the right hand side and compare it with tangent matrix derivative...
Definition: OperatorsTester.hpp:21
MoFEM::UnknownInterface
base class for all interface classes
Definition: UnknownInterface.hpp:34
MoFEM::OperatorsTester::assembleVec
SmartPetscObj< Vec > assembleVec(SmartPetscObj< DM > dm, std::string fe_name, boost::shared_ptr< FEMethod > pipeline, SmartPetscObj< Vec > x, SmartPetscObj< Vec > delta_x, SmartPetscObj< Vec > delta2_x, double time, double delta_t, CacheTupleWeakPtr cache_ptr)
Assemble the right hand side vector.
Definition: OperatorsTester.cpp:83
eps
static const double eps
Definition: check_base_functions_derivatives_on_tet.cpp:11
MoFEM::OperatorsTester::RandomFieldData
std::pair< std::string, std::array< double, 2 > > RandomFieldData
Definition: OperatorsTester.hpp:29
MoFEM::CacheTupleWeakPtr
boost::weak_ptr< CacheTuple > CacheTupleWeakPtr
Definition: FEMultiIndices.hpp:494
MoFEM::SmartPetscObj< Vec >
MoFEM::OperatorsTester::checkCentralFiniteDifference
SmartPetscObj< Vec > checkCentralFiniteDifference(SmartPetscObj< DM > dm, std::string fe_name, boost::shared_ptr< FEMethod > pipeline_rhs, boost::shared_ptr< FEMethod > pipeline_lhs, SmartPetscObj< Vec > x, SmartPetscObj< Vec > delta_x, SmartPetscObj< Vec > delta2_x, SmartPetscObj< Vec > diff_x, double time, double delta_t, double eps, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Check consistency between directional derivative with matrix.
Definition: OperatorsTester.cpp:210