MAST
complex_assembly_base.h
Go to the documentation of this file.
1 /*
2  * MAST: Multidisciplinary-design Adaptation and Sensitivity Toolkit
3  * Copyright (C) 2013-2019 Manav Bhatia
4  *
5  * This library is free software; you can redistribute it and/or
6  * modify it under the terms of the GNU Lesser General Public
7  * License as published by the Free Software Foundation; either
8  * version 2.1 of the License, or (at your option) any later version.
9  *
10  * This library is distributed in the hope that it will be useful,
11  * but WITHOUT ANY WARRANTY; without even the implied warranty of
12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
13  * Lesser General Public License for more details.
14  *
15  * You should have received a copy of the GNU Lesser General Public
16  * License along with this library; if not, write to the Free Software
17  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
18  */
19 
20 #ifndef __mast__complex_assembly_base_h__
21 #define __mast__complex_assembly_base_h__
22 
23 
24 // MAST includes
25 #include "base/assembly_base.h"
26 
27 // libMesh includes
28 #include "libmesh/nonlinear_implicit_system.h"
29 
30 
31 
32 namespace MAST {
33 
34  // Forward declerations
35  class ComplexSolverBase;
36  class Parameter;
37 
38 
40  public MAST::AssemblyBase {
41  public:
42 
47 
48 
53  virtual ~ComplexAssemblyBase();
54 
55 
59  Real residual_l2_norm(const libMesh::NumericVector<Real>& real,
60  const libMesh::NumericVector<Real>& imag);
61 
62 
70  void set_base_solution(const libMesh::NumericVector<Real>& sol,
71  bool if_sens = false);
72 
73 
79  void clear_base_solution(bool if_sens = false);
80 
81 
87 
88 
94  const libMesh::NumericVector<Real>& base_sol(bool if_sens = false) const;
95 
96 
97  void
98  residual_and_jacobian_field_split (const libMesh::NumericVector<Real>& X_R,
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);
104 
113  void
114  residual_and_jacobian_blocked (const libMesh::NumericVector<Real>& X,
115  libMesh::NumericVector<Real>& R,
116  libMesh::SparseMatrix<Real>& J,
117  MAST::Parameter* p = nullptr);
118 
130  virtual bool
132  libMesh::NumericVector<Real>& sensitivity_rhs);
133 
134  protected:
135 
141  const libMesh::NumericVector<Real> * _base_sol;
142 
148  const libMesh::NumericVector<Real> * _base_sol_sensitivity;
149  };
150 
151 }
152 
153 
154 
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...
Definition: parameter.h:35
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...
libMesh::Real Real
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