waves
Basic FE playground
flow::Adjoint Class Reference

Adjoint solver class. More...

#include <wAdjoint.h>

Inheritance diagram for flow::Adjoint:
Collaboration diagram for flow::Adjoint:

Public Member Functions

 Adjoint (std::shared_ptr< Newton > _sol, std::shared_ptr< tbox::MshDeform > _morph, std::shared_ptr< tbox::LinearSolver > _linsol)
 Initialize the adjoint solver. More...
 
virtual ~Adjoint ()
 
void run ()
 Run the adjoint method. More...
 
void save (tbox::MshExport *mshWriter, int n=0)
 Write the results. More...
 
virtual void write (std::ostream &out) const override
 

Public Attributes

std::shared_ptr< Newtonsol
 Newton solver. More...
 
std::shared_ptr< tbox::MshDeform > morph
 mesh morpher More...
 
std::shared_ptr< tbox::LinearSolver > linsol
 linear solver More...
 
int nthreads
 number of threads for the assembly More...
 
int verbose
 display more info More...
 
std::vector< double > lamClFlo
 lift flow adjoint More...
 
std::vector< double > lamCdFlo
 drag flow adjoint More...
 
std::vector< Eigen::Vector3d > dClMsh
 lift gradient with respect to mesh More...
 
std::vector< Eigen::Vector3d > dCdMsh
 drag gradient with respect to mesh More...
 
double dClAoa
 lift gradient with respect to angle of attack More...
 
double dCdAoa
 drag gradient with respect to angle of attack More...
 

Private Member Functions

void solve ()
 Solve the adjoint steady transonic Full Potential Equation. More...
 
void computeGradientAoa ()
 Compute total gradients of the loads with respect to angle of attack. More...
 
void computeGradientMesh ()
 Compute total gradients of the loads with respect to mesh nodes. More...
 
void buildGradientLoadsFlow (Eigen::VectorXd &dL, Eigen::VectorXd &dD)
 Build the gradients of the lift and drag with respect to the flow variables. More...
 
void buildGradientResidualAoa (Eigen::VectorXd &dR)
 Compute gradients of the residual with respect to angle of attack. More...
 
void buildGradientLoadsAoa (double &dL, double &dD)
 Compute gradients of the lift and drag with respect to angle of attack. More...
 
void buildGradientResidualMesh (Eigen::SparseMatrix< double, Eigen::RowMajor > &dR)
 Compute gradients of the residual with respect to mesh nodes. More...
 
void buildGradientLoadsMesh (Eigen::RowVectorXd &dL, Eigen::RowVectorXd &dD)
 Compute gradients of the lift and drag with respect to mesh nodes. More...
 
void mapRows ()
 Create the rows vector for the adjoint mesh. More...
 

Private Attributes

fwk::Timers tms
 internal timers More...
 
std::vector< std::vector< int > > arows
 rows (unknown index) with duplicated dof on wake More...
 

Detailed Description

Adjoint solver class.

Authors
Adrien Crovato

Constructor & Destructor Documentation

◆ Adjoint()

Adjoint::Adjoint ( std::shared_ptr< Newton _sol,
std::shared_ptr< tbox::MshDeform >  _morph,
std::shared_ptr< tbox::LinearSolver >  _linsol 
)

Initialize the adjoint solver.

Authors
Adrien Crovato

◆ ~Adjoint()

virtual flow::Adjoint::~Adjoint ( )
inlinevirtual

Member Function Documentation

◆ buildGradientLoadsAoa()

void Adjoint::buildGradientLoadsAoa ( double &  dL,
double &  dD 
)
private

Compute gradients of the lift and drag with respect to angle of attack.

Todo:
one gradient for each surface node

◆ buildGradientLoadsFlow()

void Adjoint::buildGradientLoadsFlow ( Eigen::VectorXd &  dL,
Eigen::VectorXd &  dD 
)
private

Build the gradients of the lift and drag with respect to the flow variables.

Todo:
one gradient for each surface node

◆ buildGradientLoadsMesh()

void Adjoint::buildGradientLoadsMesh ( Eigen::RowVectorXd &  dL,
Eigen::RowVectorXd &  dD 
)
private

Compute gradients of the lift and drag with respect to mesh nodes.

Todo:
one gradient for each surface node

◆ buildGradientResidualAoa()

void Adjoint::buildGradientResidualAoa ( Eigen::VectorXd &  dR)
private

Compute gradients of the residual with respect to angle of attack.

◆ buildGradientResidualMesh()

void Adjoint::buildGradientResidualMesh ( Eigen::SparseMatrix< double, Eigen::RowMajor > &  dR)
private

Compute gradients of the residual with respect to mesh nodes.

◆ computeGradientAoa()

void Adjoint::computeGradientAoa ( )
private

Compute total gradients of the loads with respect to angle of attack.

Todo:
one gradient for each surface node

◆ computeGradientMesh()

void Adjoint::computeGradientMesh ( )
private

Compute total gradients of the loads with respect to mesh nodes.

◆ mapRows()

void Adjoint::mapRows ( )
private

Create the rows vector for the adjoint mesh.

Todo:
duplicated from wake

◆ run()

void Adjoint::run ( )

Run the adjoint method.

Solve the adjoint steady transonic Full Potential Equation and compute lift and drag sensitivites

◆ save()

void Adjoint::save ( tbox::MshExport *  mshWriter,
int  n = 0 
)

Write the results.

◆ solve()

void Adjoint::solve ( )
private

Solve the adjoint steady transonic Full Potential Equation.

◆ write()

void Adjoint::write ( std::ostream &  out) const
overridevirtual

Member Data Documentation

◆ arows

std::vector<std::vector<int> > flow::Adjoint::arows
private

rows (unknown index) with duplicated dof on wake

◆ dCdAoa

double flow::Adjoint::dCdAoa

drag gradient with respect to angle of attack

◆ dCdMsh

std::vector<Eigen::Vector3d> flow::Adjoint::dCdMsh

drag gradient with respect to mesh

◆ dClAoa

double flow::Adjoint::dClAoa

lift gradient with respect to angle of attack

◆ dClMsh

std::vector<Eigen::Vector3d> flow::Adjoint::dClMsh

lift gradient with respect to mesh

◆ lamCdFlo

std::vector<double> flow::Adjoint::lamCdFlo

drag flow adjoint

◆ lamClFlo

std::vector<double> flow::Adjoint::lamClFlo

lift flow adjoint

◆ linsol

std::shared_ptr<tbox::LinearSolver> flow::Adjoint::linsol

linear solver

◆ morph

std::shared_ptr<tbox::MshDeform> flow::Adjoint::morph

mesh morpher

◆ nthreads

int flow::Adjoint::nthreads

number of threads for the assembly

◆ sol

std::shared_ptr<Newton> flow::Adjoint::sol

Newton solver.

◆ tms

fwk::Timers flow::Adjoint::tms
private

internal timers

◆ verbose

int flow::Adjoint::verbose

display more info


The documentation for this class was generated from the following files: