MAST
Multidisciplinary-design Adaptation and Sensitivity Toolkit (MAST)
fe_base.h
Go to the documentation of this file.
1 /*
2  * MAST: Multidisciplinary-design Adaptation and Sensitivity Toolkit
3  * Copyright (C) 2013-2020 Manav Bhatia and MAST authors
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 
21 #ifndef __mast_fe_base_h__
22 #define __mast_fe_base_h__
23 
24 // MAST includes
25 #include "base/mast_data_types.h"
26 
27 // libMesh includes
28 #include "libmesh/elem.h"
29 #include "libmesh/fe_base.h"
30 #include "libmesh/quadrature.h"
31 
32 
33 namespace MAST {
34 
35  // forward declerations
36  class SystemInitialization;
37  class GeomElem;
38 
39  class FEBase {
40 
41  public:
42 
44 
45  virtual ~FEBase();
46 
47  bool initialized() const { return _initialized;}
48 
53  void set_extra_quadrature_order(int n);
54 
59 
71  virtual void init(const MAST::GeomElem& elem,
72  bool init_grads,
73  const std::vector<libMesh::Point>* pts = nullptr);
74 
75 
83  virtual void init_for_side(const MAST::GeomElem& elem,
84  unsigned int s,
85  bool if_calculate_dphi);
86 
87 
88  libMesh::FEType
89  get_fe_type() const;
90 
91  virtual const std::vector<Real>&
92  get_JxW() const;
93 
98  virtual const std::vector<libMesh::Point>&
99  get_xyz() const;
100 
101  virtual unsigned int
102  n_shape_functions() const;
103 
104  virtual const std::vector<std::vector<Real> >&
105  get_phi() const;
106 
107  virtual const std::vector<std::vector<libMesh::RealVectorValue> >&
108  get_dphi() const;
109 
110  virtual const std::vector<std::vector<libMesh::RealTensorValue>>&
111  get_d2phi() const;
112 
113  virtual const std::vector<Real>&
114  get_dxidx() const;
115 
116  virtual const std::vector<Real>&
117  get_dxidy() const;
118 
119  virtual const std::vector<Real>&
120  get_dxidz() const;
121 
122  virtual const std::vector<Real>&
123  get_detadx() const;
124 
125  virtual const std::vector<Real>&
126  get_detady() const;
127 
128  virtual const std::vector<Real>&
129  get_detadz() const;
130 
131  virtual const std::vector<Real>&
132  get_dzetadx() const;
133 
134  virtual const std::vector<Real>&
135  get_dzetady() const;
136 
137  virtual const std::vector<Real>&
138  get_dzetadz() const;
139 
140  virtual const std::vector<libMesh::RealVectorValue>&
141  get_dxyzdxi() const;
142 
143  virtual const std::vector<libMesh::RealVectorValue>&
144  get_dxyzdeta() const;
145 
146  virtual const std::vector<libMesh::RealVectorValue>&
147  get_dxyzdzeta() const;
148 
149  virtual const std::vector<std::vector<Real> >&
150  get_dphidxi() const;
151 
152  virtual const std::vector<std::vector<Real> >&
153  get_dphideta() const;
154 
155  virtual const std::vector<std::vector<Real> >&
156  get_dphidzeta() const;
157 
162  virtual const std::vector<libMesh::Point>&
164 
169  virtual const std::vector<libMesh::Point>&
171 
172 
173  virtual const std::vector<libMesh::Point>&
174  get_qpoints() const;
175 
176  virtual const libMesh::QBase&
177  get_qrule() const;
178 
179  protected:
180 
186  libMesh::FEBase* _fe;
187  libMesh::QBase* _qrule;
188  std::vector<libMesh::Point> _qpoints;
189  std::vector<libMesh::Point> _global_xyz;
190  std::vector<libMesh::Point> _local_normals;
191  std::vector<libMesh::Point> _global_normals;
192  };
193 }
194 
195 
196 #endif // __mast_fe_base_h__
virtual ~FEBase()
Definition: fe_base.cpp:39
virtual const std::vector< Real > & get_dxidy() const
Definition: fe_base.cpp:279
std::vector< libMesh::Point > _global_normals
Definition: fe_base.h:191
virtual const std::vector< Real > & get_dzetadz() const
Definition: fe_base.cpp:335
virtual const std::vector< Real > & get_dxidz() const
Definition: fe_base.cpp:287
libMesh::QBase * _qrule
Definition: fe_base.h:187
virtual const std::vector< libMesh::Point > & get_qpoints() const
Definition: fe_base.cpp:407
virtual const std::vector< Real > & get_detadz() const
Definition: fe_base.cpp:311
FEBase(const MAST::SystemInitialization &sys)
Definition: fe_base.cpp:27
virtual const std::vector< Real > & get_dzetady() const
Definition: fe_base.cpp:327
bool _init_second_order_derivatives
Definition: fe_base.h:183
virtual const std::vector< libMesh::Point > & get_xyz() const
physical location of the quadrature point in the global coordinate system for the reference element ...
Definition: fe_base.cpp:228
virtual const std::vector< std::vector< libMesh::RealTensorValue > > & get_d2phi() const
Definition: fe_base.cpp:263
virtual unsigned int n_shape_functions() const
Definition: fe_base.cpp:239
std::vector< libMesh::Point > _local_normals
Definition: fe_base.h:190
virtual const std::vector< std::vector< libMesh::RealVectorValue > > & get_dphi() const
Definition: fe_base.cpp:255
std::vector< libMesh::Point > _qpoints
Definition: fe_base.h:188
virtual const std::vector< std::vector< Real > > & get_dphidxi() const
Definition: fe_base.cpp:367
virtual const std::vector< libMesh::RealVectorValue > & get_dxyzdeta() const
Definition: fe_base.cpp:351
bool _initialized
Definition: fe_base.h:184
virtual const std::vector< std::vector< Real > > & get_phi() const
Definition: fe_base.cpp:247
virtual const std::vector< libMesh::Point > & get_normals_for_reference_coordinate() const
normals defined in the global coordinate system for the reference element
Definition: fe_base.cpp:391
virtual void init_for_side(const MAST::GeomElem &elem, unsigned int s, bool if_calculate_dphi)
Initializes the quadrature and finite element for element side integration.
Definition: fe_base.cpp:137
virtual const std::vector< std::vector< Real > > & get_dphideta() const
Definition: fe_base.cpp:375
const MAST::SystemInitialization & _sys
Definition: fe_base.h:181
void set_evaluate_second_order_derivatives(bool f)
sets the flag for evaluation of second order derivative
Definition: fe_base.cpp:57
libMesh::FEBase * _fe
Definition: fe_base.h:186
void set_extra_quadrature_order(int n)
this is used, in addition to libMesh::System::extra_quadrature_order to set the quadrature rule...
Definition: fe_base.cpp:50
virtual const std::vector< libMesh::RealVectorValue > & get_dxyzdzeta() const
Definition: fe_base.cpp:359
std::vector< libMesh::Point > _global_xyz
Definition: fe_base.h:189
bool initialized() const
Definition: fe_base.h:47
virtual const std::vector< libMesh::Point > & get_normals_for_local_coordinate() const
normals defined in the coordinate system for the local reference element.
Definition: fe_base.cpp:399
This class acts as a wrapper around libMesh::Elem for the purpose of providing a uniform interface fo...
Definition: geom_elem.h:59
virtual void init(const MAST::GeomElem &elem, bool init_grads, const std::vector< libMesh::Point > *pts=nullptr)
Initializes the quadrature and finite element for element volume integration.
Definition: fe_base.cpp:64
virtual const std::vector< Real > & get_dzetadx() const
Definition: fe_base.cpp:319
libMesh::FEType get_fe_type() const
Definition: fe_base.cpp:212
virtual const std::vector< Real > & get_JxW() const
Definition: fe_base.cpp:220
virtual const std::vector< Real > & get_detady() const
Definition: fe_base.cpp:303
virtual const std::vector< Real > & get_detadx() const
Definition: fe_base.cpp:295
virtual const std::vector< libMesh::RealVectorValue > & get_dxyzdxi() const
Definition: fe_base.cpp:343
virtual const std::vector< Real > & get_dxidx() const
Definition: fe_base.cpp:271
virtual const std::vector< std::vector< Real > > & get_dphidzeta() const
Definition: fe_base.cpp:383
const MAST::GeomElem * _elem
Definition: fe_base.h:185
virtual const libMesh::QBase & get_qrule() const
Definition: fe_base.cpp:418
unsigned int _extra_quadrature_order
Definition: fe_base.h:182