waves
Basic FE playground
wAdjoint.h
Go to the documentation of this file.
1 /*
2  * Copyright 2022 University of Liège
3  *
4  * Licensed under the Apache License, Version 2.0 (the "License");
5  * you may not use this file except in compliance with the License.
6  * You may obtain a copy of the License at
7  *
8  * http://www.apache.org/licenses/LICENSE-2.0
9  *
10  * Unless required by applicable law or agreed to in writing, software
11  * distributed under the License is distributed on an "AS IS" BASIS,
12  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13  * See the License for the specific language governing permissions and
14  * limitations under the License.
15  */
16 
17 #ifndef WADJOINT_H
18 #define WADJOINT_H
19 
20 #include "flow.h"
21 #include "wObject.h"
22 #include "wTimers.h"
23 
24 #include <iostream>
25 #include <vector>
26 #include <memory>
27 #include <Eigen/Sparse>
28 
29 namespace flow
30 {
31 
36 class FLOW_API Adjoint : public fwk::wSharedObject
37 {
38 public:
39  std::shared_ptr<Newton> sol;
40  std::shared_ptr<tbox::MshDeform> morph;
41  std::shared_ptr<tbox::LinearSolver> linsol;
42 
43  int nthreads;
44  int verbose;
45 
46  std::vector<double> lamClFlo;
47  std::vector<double> lamCdFlo;
48  std::vector<Eigen::Vector3d> dClMsh;
49  std::vector<Eigen::Vector3d> dCdMsh;
50  double dClAoa;
51  double dCdAoa;
52 
53 private:
54  fwk::Timers tms;
55  std::vector<std::vector<int>> arows;
56 
57 public:
58  Adjoint(std::shared_ptr<Newton> _sol, std::shared_ptr<tbox::MshDeform> _morph, std::shared_ptr<tbox::LinearSolver> _linsol);
59  virtual ~Adjoint() { std::cout << "~Adjoint()\n"; }
60 
61  void run();
62  void save(tbox::MshExport *mshWriter, int n = 0);
63 
64 #ifndef SWIG
65  virtual void write(std::ostream &out) const override;
66 #endif
67 
68 private:
69  // Adjoint solution and total gradients
70  void solve();
71  void computeGradientAoa();
72  void computeGradientMesh();
73  // Partial gradients wrt with respect to flow variables
74  void buildGradientLoadsFlow(Eigen::VectorXd &dL, Eigen::VectorXd &dD);
75  // Partial gradients with respect to angle of attack
76  void buildGradientResidualAoa(Eigen::VectorXd &dR);
77  void buildGradientLoadsAoa(double &dL, double &dD);
78  // Partial gradients with respect to mesh
79  void buildGradientResidualMesh(Eigen::SparseMatrix<double, Eigen::RowMajor> &dR);
80  void buildGradientLoadsMesh(Eigen::RowVectorXd &dL, Eigen::RowVectorXd &dD);
81 
82  // @todo
83  void mapRows();
84 };
85 
86 } // namespace flow
87 #endif //WADJOINT_H
flow::Adjoint::dClMsh
std::vector< Eigen::Vector3d > dClMsh
lift gradient with respect to mesh
Definition: wAdjoint.h:48
flow::Adjoint::linsol
std::shared_ptr< tbox::LinearSolver > linsol
linear solver
Definition: wAdjoint.h:41
flow::Adjoint::lamCdFlo
std::vector< double > lamCdFlo
drag flow adjoint
Definition: wAdjoint.h:47
flow::Adjoint::dCdMsh
std::vector< Eigen::Vector3d > dCdMsh
drag gradient with respect to mesh
Definition: wAdjoint.h:49
FLOW_API
#define FLOW_API
Definition: flow.h:29
flow::Adjoint::lamClFlo
std::vector< double > lamClFlo
lift flow adjoint
Definition: wAdjoint.h:46
flow
Namespace for flow module.
Definition: flow.h:37
flow::Adjoint::sol
std::shared_ptr< Newton > sol
Newton solver.
Definition: wAdjoint.h:39
flow::Adjoint::arows
std::vector< std::vector< int > > arows
rows (unknown index) with duplicated dof on wake
Definition: wAdjoint.h:55
flow::Adjoint::verbose
int verbose
display more info
Definition: wAdjoint.h:44
flow::Adjoint::morph
std::shared_ptr< tbox::MshDeform > morph
mesh morpher
Definition: wAdjoint.h:40
flow::Adjoint::tms
fwk::Timers tms
internal timers
Definition: wAdjoint.h:54
flow::Adjoint
Adjoint solver class.
Definition: wAdjoint.h:36
flow::Adjoint::dCdAoa
double dCdAoa
drag gradient with respect to angle of attack
Definition: wAdjoint.h:51
flow::Adjoint::dClAoa
double dClAoa
lift gradient with respect to angle of attack
Definition: wAdjoint.h:50
flow::Adjoint::~Adjoint
virtual ~Adjoint()
Definition: wAdjoint.h:59
flow::Adjoint::nthreads
int nthreads
number of threads for the assembly
Definition: wAdjoint.h:43
flow.h