MAST
structural_buckling_eigenproblem_assembly.cpp
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 // MAST includes
27 #include "base/assembly_base.h"
28 #include "base/nonlinear_system.h"
29 #include "base/parameter.h"
30 #include "numerics/utility.h"
31 #include "mesh/geom_elem.h"
32 
33 // libMesh includes
34 #include "libmesh/dof_map.h"
35 #include "libmesh/numeric_vector.h"
36 
37 
41 _use_linearized_formulation(true),
42 _load_param(nullptr),
43 _lambda1(0.),
44 _lambda2(0.),
45 _sol1(nullptr),
46 _sol2(nullptr) {
47 
48 }
49 
50 
51 
52 
55 { }
56 
57 
58 
59 
60 void
62 set_buckling_data (bool use_linearized_form,
63  MAST::Parameter& p,
64  const Real lambda1,
65  const Real lambda2,
66  libMesh::NumericVector<Real>& x1,
67  libMesh::NumericVector<Real>& x2) {
68 
69  _use_linearized_formulation = use_linearized_form;
70  _load_param = &p;
71  _lambda1 = lambda1;
72  _lambda2 = lambda2;
73  _sol1 = &x1;
74  _sol2 = &x2;
75 }
76 
77 
78 void
80 
82  _load_param = nullptr;
83  _lambda1 = 0.;
84  _lambda2 = 0.;
85  _sol1 = nullptr;
86  _sol2 = nullptr;
87 
89 }
90 
91 
92 
93 void
95 eigenproblem_assemble(libMesh::SparseMatrix<Real> *A,
96  libMesh::SparseMatrix<Real> *B) {
97 
98  libmesh_assert(_system);
99  libmesh_assert(_discipline);
100  libmesh_assert(_elem_ops);
101 
102  MAST::NonlinearSystem& eigen_sys =
103  dynamic_cast<MAST::NonlinearSystem&>(_system->system());
104 
105  // create the localized solution vectors
106  std::unique_ptr<libMesh::NumericVector<Real> >
107  localized_solution1,
108  localized_solution2;
109 
110  localized_solution1.reset(build_localized_vector(eigen_sys,
111  *_sol1).release());
112  localized_solution2.reset(build_localized_vector(eigen_sys,
113  *_sol2).release());
114 
115  libMesh::SparseMatrix<Real>
116  &matrix_A = *A,
117  &matrix_B = *B;
118 
119  matrix_A.zero();
120  matrix_B.zero();
121 
122  // iterate over each element, initialize it and get the relevant
123  // analysis quantities
124  RealVectorX sol;
125  RealMatrixX dummy, mat_A, mat_B;
126  std::vector<libMesh::dof_id_type> dof_indices;
127  const libMesh::DofMap& dof_map = eigen_sys.get_dof_map();
128 
129 
130  libMesh::MeshBase::const_element_iterator el =
131  eigen_sys.get_mesh().active_local_elements_begin();
132  const libMesh::MeshBase::const_element_iterator end_el =
133  eigen_sys.get_mesh().active_local_elements_end();
134 
136  &ops = dynamic_cast<MAST::EigenproblemAssemblyElemOperations&>(*_elem_ops);
137 
138  for ( ; el != end_el; ++el) {
139 
140  const libMesh::Elem* elem = *el;
141 
142  dof_map.dof_indices (elem, dof_indices);
143 
144  // get the solution
145  unsigned int ndofs = (unsigned int)dof_indices.size();
146  sol.setZero(ndofs);
147  dummy.setZero(ndofs, ndofs);
148  mat_A.setZero(ndofs, ndofs);
149  mat_B.setZero(ndofs, ndofs);
150 
152  // first calculate the tangent stiffness matrix about the
153  // first load factor and solution
155  (*_load_param) = _lambda1;
156  for (unsigned int i=0; i<dof_indices.size(); i++)
157  sol(i) = (*localized_solution1)(dof_indices[i]);
158 
159  MAST::GeomElem geom_elem;
160  ops.set_elem_data(elem->dim(), *elem, geom_elem);
161  geom_elem.init(*elem, *_system);
162 
163  ops.init(geom_elem);
164  ops.set_elem_solution(sol);
165  ops.elem_calculations(mat_A, dummy);
166 
167  DenseRealMatrix AA, BB;
168  MAST::copy(AA, mat_A); // copy to the libMesh matrix for further processing
169  dof_map.constrain_element_matrix(AA, dof_indices); // constrain the element matrices.
170  matrix_A.add_matrix (AA, dof_indices); // add to the global matrices
171 
172 
174  // next, calculate the tangent stiffness matrix about the
175  // second load factor and solution
177  (*_load_param) = _lambda2;
178  for (unsigned int i=0; i<dof_indices.size(); i++)
179  sol(i) = (*localized_solution2)(dof_indices[i]);
180 
181  ops.set_elem_solution(sol);
182  ops.elem_calculations(mat_B, dummy);
183  mat_B *= -1.;
185  mat_B += mat_A;
186 
187  MAST::copy(BB, mat_B); // copy to the libMesh matrix for further processing
188  dof_map.constrain_element_matrix(BB, dof_indices); // constrain the element matrices.
189  matrix_B.add_matrix (BB, dof_indices); // add to the global matrices
190  }
191 
192  // finalize the matrices for futher use.
193  A->close();
194  B->close();
195 }
196 
197 
198 
199 
200 
201 
202 
203 
204 Real
207 
208  Real rval = 0.;
209 
211  rval = _lambda1 + v * (_lambda2 - _lambda1);
212  else
213  rval = _lambda1 + v/(1.+v) * (_lambda2 - _lambda1);
214 
215  return rval;
216 }
217 
218 
219 
220 
MAST::AssemblyElemOperations * _elem_ops
provides assembly elem operations for use by this class
MAST::NonlinearSystem & system()
libMesh::DenseMatrix< Real > DenseRealMatrix
StructuralBucklingEigenproblemAssembly()
constructor associates the eigen system with this assembly object
Real critical_point_estimate_from_eigenproblem(Real v) const
calculates the critical load factor based on the eigensolution
This class implements a system for solution of nonlinear systems.
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
virtual void elem_calculations(RealMatrixX &mat_A, RealMatrixX &mat_B)=0
performs the element calculations over elem, and returns the element matrices for the eigenproblem ...
virtual void set_elem_solution(const RealVectorX &sol)
sets the element solution
libMesh::Real Real
MAST::SystemInitialization * _system
System for which this assembly is performed.
virtual void clear_discipline_and_system()
clears the states and load factors for buckling eigenproblem
MAST::PhysicsDisciplineBase * _discipline
PhysicsDisciplineBase object for which this class is assembling.
bool _use_linearized_formulation
whether or not to use the linearized formulation
virtual void set_elem_data(unsigned int dim, const libMesh::Elem &ref_elem, MAST::GeomElem &elem) const =0
some analyses may want to set additional element data before initialization of the GeomElem...
virtual void clear_discipline_and_system()
clears association with a system to this discipline
Assembles the system of equations for an eigenproblem of type .
Real _lambda1
values of load factors for which the two stiffness matrices are calculated for the buckling eigenvalu...
Matrix< Real, Dynamic, Dynamic > RealMatrixX
virtual ~StructuralBucklingEigenproblemAssembly()
destructor resets the association with the eigen system from this assembly object ...
void copy(DenseRealMatrix &m1, const RealMatrixX &m2)
Definition: utility.h:167
Matrix< Real, Dynamic, 1 > RealVectorX
virtual void init(const MAST::GeomElem &elem)=0
initializes the object for calculation of element quantities for the specified elem.
std::unique_ptr< libMesh::NumericVector< Real > > build_localized_vector(const libMesh::System &sys, const libMesh::NumericVector< Real > &global) const
localizes the parallel vector so that the local copy stores all values necessary for calculation of t...
This class acts as a wrapper around libMesh::Elem for the purpose of providing a uniform interface fo...
Definition: geom_elem.h:59
void set_buckling_data(bool use_linearized_approach, MAST::Parameter &p, const Real lambda1, const Real lambda2, libMesh::NumericVector< Real > &x1, libMesh::NumericVector< Real > &x2)
set the states and load factors for buckling eigenproblem.
virtual void eigenproblem_assemble(libMesh::SparseMatrix< Real > *A, libMesh::SparseMatrix< Real > *B)
assembles the global A and B matrices for the modal eigenvalue problem
MAST::Parameter * _load_param
load parameter used to define the two stiffness matrices
libMesh::NumericVector< Real > * _sol1
the equilibrium solutions associated with _lambda1 and _lambda2 load factors.