20 #ifndef __mast__complex_assembly_base_h__ 21 #define __mast__complex_assembly_base_h__ 28 #include "libmesh/nonlinear_implicit_system.h" 35 class ComplexSolverBase;
60 const libMesh::NumericVector<Real>& imag);
71 bool if_sens =
false);
94 const libMesh::NumericVector<Real>&
base_sol(
bool if_sens =
false)
const;
99 const libMesh::NumericVector<Real>& X_I,
100 libMesh::NumericVector<Real>& R_R,
101 libMesh::NumericVector<Real>& R_I,
102 libMesh::SparseMatrix<Real>& J_R,
103 libMesh::SparseMatrix<Real>& J_I);
115 libMesh::NumericVector<Real>& R,
116 libMesh::SparseMatrix<Real>& J,
132 libMesh::NumericVector<Real>& sensitivity_rhs);
155 #endif //__mast__complex_assembly_base_h__
void residual_and_jacobian_blocked(const libMesh::NumericVector< Real > &X, libMesh::NumericVector< Real > &R, libMesh::SparseMatrix< Real > &J, MAST::Parameter *p=nullptr)
Assembles the residual and Jacobian of the N complex system of equations split into 2N real system of...
void residual_and_jacobian_field_split(const libMesh::NumericVector< Real > &X_R, const libMesh::NumericVector< Real > &X_I, libMesh::NumericVector< Real > &R_R, libMesh::NumericVector< Real > &R_I, libMesh::SparseMatrix< Real > &J_R, libMesh::SparseMatrix< Real > &J_I)
virtual ~ComplexAssemblyBase()
destructor resets the association of this assembly object with the system
const libMesh::NumericVector< Real > * _base_sol_sensitivity
sensitivity of base solution may be needed for sensitivity analysis.
This is a scalar function whose value can be changed and one that can be used as a design variable in...
Real residual_l2_norm(const libMesh::NumericVector< Real > &real, const libMesh::NumericVector< Real > &imag)
calculates the L2 norm of the residual complex system of equations.
void set_base_solution(const libMesh::NumericVector< Real > &sol, bool if_sens=false)
if the problem is defined about a non-zero base solution, then this method provides the object with t...
bool if_linearized_about_nonzero_solution() const
ComplexAssemblyBase()
constructor associates this assembly object with the system
void clear_base_solution(bool if_sens=false)
Clears pointer to the solution vector The flag if_sens tells the method to clear pointer of the solut...
const libMesh::NumericVector< Real > * _base_sol
base solution about which this problem is defined.
virtual bool sensitivity_assemble(const MAST::FunctionBase &f, libMesh::NumericVector< Real > &sensitivity_rhs)
Assembly function.
const libMesh::NumericVector< Real > & base_sol(bool if_sens=false) const