libMesh
Loading...
Searching...
No Matches
Public Member Functions | Private Attributes | List of all members
LinearElasticity Class Reference
Inheritance diagram for LinearElasticity:
[legend]

Public Member Functions

 LinearElasticity (EquationSystems &es_in)
 
Real kronecker_delta (unsigned int i, unsigned int j)
 Kronecker delta function.
 
Real elasticity_tensor (unsigned int i, unsigned int j, unsigned int k, unsigned int l)
 Evaluate the fourth order tensor (C_ijkl) that relates stress to strain.
 
void assemble ()
 Assemble the system matrix and right-hand side vector.
 
void compute_stresses ()
 
 LinearElasticity (EquationSystems &es_in)
 
Real kronecker_delta (unsigned int i, unsigned int j)
 Kronecker delta function.
 
Real elasticity_tensor (unsigned int i, unsigned int j, unsigned int k, unsigned int l)
 Evaluate the fourth order tensor (C_ijkl) that relates stress to strain.
 
void assemble ()
 Assemble the system matrix and right-hand side vector.
 

Private Attributes

EquationSystemses
 

Detailed Description

Definition at line 113 of file systems_of_equations_ex6.C.

Constructor & Destructor Documentation

◆ LinearElasticity() [1/2]

LinearElasticity::LinearElasticity ( EquationSystems es_in)
inline

Definition at line 120 of file systems_of_equations_ex6.C.

120 :
121 es(es_in)
122 {}

◆ LinearElasticity() [2/2]

LinearElasticity::LinearElasticity ( EquationSystems es_in)
inline

Definition at line 201 of file systems_of_equations_ex9.C.

201 :
202 es(es_in)
203 {}

Member Function Documentation

◆ assemble() [1/2]

void LinearElasticity::assemble ( )
inline

Assemble the system matrix and right-hand side vector.

Definition at line 156 of file systems_of_equations_ex6.C.

157 {
158 const MeshBase & mesh = es.get_mesh();
159
160 const unsigned int dim = mesh.mesh_dimension();
161
162 LinearImplicitSystem & system = es.get_system<LinearImplicitSystem>("Elasticity");
163
164 const unsigned int u_var = system.variable_number ("u");
165
166 const DofMap & dof_map = system.get_dof_map();
167 FEType fe_type = dof_map.variable_type(u_var);
168 std::unique_ptr<FEBase> fe (FEBase::build(dim, fe_type));
169 QGauss qrule (dim, fe_type.default_quadrature_order());
170 fe->attach_quadrature_rule (&qrule);
171
172 std::unique_ptr<FEBase> fe_face (FEBase::build(dim, fe_type));
173 QGauss qface(dim-1, fe_type.default_quadrature_order());
174 fe_face->attach_quadrature_rule (&qface);
175
176 const std::vector<Real> & JxW = fe->get_JxW();
177 const std::vector<std::vector<Real>> & phi = fe->get_phi();
178 const std::vector<std::vector<RealGradient>> & dphi = fe->get_dphi();
179
181 DenseSubMatrix<Number> Ke_var[3][3] =
182 {
186 };
187
189
190 DenseSubVector<Number> Fe_var[3] =
194
195 std::vector<dof_id_type> dof_indices;
196 std::vector<std::vector<dof_id_type>> dof_indices_var(3);
197
198 SparseMatrix<Number> & matrix = system.get_system_matrix();
199
200 for (const auto & elem : mesh.active_local_element_ptr_range())
201 {
202 dof_map.dof_indices (elem, dof_indices);
203 for (unsigned int var=0; var<3; var++)
204 dof_map.dof_indices (elem, dof_indices_var[var], var);
205
206 const unsigned int n_dofs = dof_indices.size();
207 const unsigned int n_var_dofs = dof_indices_var[0].size();
208
209 fe->reinit (elem);
210
211 Ke.resize (n_dofs, n_dofs);
212 for (unsigned int var_i=0; var_i<3; var_i++)
213 for (unsigned int var_j=0; var_j<3; var_j++)
214 Ke_var[var_i][var_j].reposition (var_i*n_var_dofs, var_j*n_var_dofs, n_var_dofs, n_var_dofs);
215
216 Fe.resize (n_dofs);
217 for (unsigned int var=0; var<3; var++)
218 Fe_var[var].reposition (var*n_var_dofs, n_var_dofs);
219
220 for (unsigned int qp=0; qp<qrule.n_points(); qp++)
221 {
222 // assemble \int_Omega C_ijkl u_k,l v_i,j \dx
223 for (unsigned int dof_i=0; dof_i<n_var_dofs; dof_i++)
224 for (unsigned int dof_j=0; dof_j<n_var_dofs; dof_j++)
225 for (unsigned int i=0; i<3; i++)
226 for (unsigned int j=0; j<3; j++)
227 for (unsigned int k=0; k<3; k++)
228 for (unsigned int l=0; l<3; l++)
229 Ke_var[i][k](dof_i,dof_j) +=
230 JxW[qp] * elasticity_tensor(i,j,k,l) * dphi[dof_j][qp](l) * dphi[dof_i][qp](j);
231
232 // assemble \int_Omega f_i v_i \dx
233 VectorValue<Number> f_vec(0., 0., -1.);
234 for (unsigned int dof_i=0; dof_i<n_var_dofs; dof_i++)
235 for (unsigned int i=0; i<3; i++)
236 Fe_var[i](dof_i) += JxW[qp] * (f_vec(i) * phi[dof_i][qp]);
237 }
238
239 // assemble \int_\Gamma g_i v_i \ds
240 VectorValue<Number> g_vec(0., 0., -1.);
241 {
242 for (auto side : elem->side_index_range())
243 if (elem->neighbor_ptr(side) == nullptr)
244 {
245 const std::vector<std::vector<Real>> & phi_face = fe_face->get_phi();
246 const std::vector<Real> & JxW_face = fe_face->get_JxW();
247
248 fe_face->reinit(elem, side);
249
250 // Apply a traction
251 for (unsigned int qp=0; qp<qface.n_points(); qp++)
252 if (mesh.get_boundary_info().has_boundary_id(elem, side, BOUNDARY_ID_MAX_X))
253 for (unsigned int dof_i=0; dof_i<n_var_dofs; dof_i++)
254 for (unsigned int i=0; i<3; i++)
255 Fe_var[i](dof_i) += JxW_face[qp] * (g_vec(i) * phi_face[dof_i][qp]);
256 }
257 }
258
259 dof_map.constrain_element_matrix_and_vector (Ke, Fe, dof_indices);
260
261 matrix.add_matrix (Ke, dof_indices);
262 system.rhs->add_vector (Fe, dof_indices);
263 }
264 }
unsigned int dim
Real elasticity_tensor(unsigned int i, unsigned int j, unsigned int k, unsigned int l)
Evaluate the fourth order tensor (C_ijkl) that relates stress to strain.
bool has_boundary_id(const Node *const node, const boundary_id_type id) const
Defines a dense matrix for use in Finite Element-type computations.
void resize(const unsigned int new_m, const unsigned int new_n)
Resizes the matrix to the specified size and calls zero().
Defines a dense submatrix for use in Finite Element-type computations.
Defines a dense subvector for use in finite element computations.
Defines a dense vector for use in Finite Element-type computations.
void resize(const unsigned int n)
Resize the vector.
This class handles the numbering of degrees of freedom on a mesh.
Definition dof_map.h:181
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
Definition dof_map.C:2201
const FEType & variable_type(const unsigned int i) const
Definition dof_map.h:2388
void constrain_element_matrix_and_vector(DenseMatrix< Number > &matrix, DenseVector< Number > &rhs, std::vector< dof_id_type > &elem_dofs, bool asymmetric_constraint_rows=true) const
Constrains the element matrix and vector.
Definition dof_map.h:2498
const MeshBase & get_mesh() const
const T_sys & get_system(std::string_view name) const
NumericVector< Number > * rhs
The system matrix.
static std::unique_ptr< FEGenericBase > build(const unsigned int dim, const FEType &type)
Builds a specific finite element type.
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition fe_type.h:197
Order default_quadrature_order() const
Definition fe_type.h:415
const SparseMatrix< Number > & get_system_matrix() const
Manages consistently variables, degrees of freedom, coefficient vectors, matrices and linear solvers ...
This is the MeshBase class.
Definition mesh_base.h:81
const BoundaryInfo & get_boundary_info() const
The information about boundary ids on the mesh.
Definition mesh_base.h:170
unsigned int mesh_dimension() const
Definition mesh_base.C:430
virtual void add_vector(const T *v, const std::vector< numeric_index_type > &dof_indices)
Computes , where v is a pointer and each dof_indices[i] specifies where to add value v[i].
This class implements specific orders of Gauss quadrature.
Generic sparse matrix.
virtual void add_matrix(const DenseMatrix< T > &dm, const std::vector< numeric_index_type > &rows, const std::vector< numeric_index_type > &cols)=0
Add the full matrix dm to the SparseMatrix.
unsigned int variable_number(std::string_view var) const
Definition system.C:1398
const DofMap & get_dof_map() const
Definition system.h:2417
This class defines a vector in LIBMESH_DIM dimensional Real or Complex space.
MeshBase & mesh

References libMesh::SparseMatrix< T >::add_matrix(), libMesh::NumericVector< T >::add_vector(), libMesh::FEGenericBase< OutputType >::build(), libMesh::DofMap::constrain_element_matrix_and_vector(), libMesh::FEType::default_quadrature_order(), dim, libMesh::DofMap::dof_indices(), elasticity_tensor(), es, libMesh::MeshBase::get_boundary_info(), libMesh::System::get_dof_map(), libMesh::EquationSystems::get_mesh(), libMesh::EquationSystems::get_system(), libMesh::ImplicitSystem::get_system_matrix(), libMesh::BoundaryInfo::has_boundary_id(), mesh, libMesh::MeshBase::mesh_dimension(), libMesh::QBase::n_points(), libMesh::DenseVector< T >::resize(), libMesh::DenseMatrix< T >::resize(), libMesh::ExplicitSystem::rhs, libMesh::System::variable_number(), and libMesh::DofMap::variable_type().

◆ assemble() [2/2]

void LinearElasticity::assemble ( )
inline

Assemble the system matrix and right-hand side vector.

Definition at line 237 of file systems_of_equations_ex9.C.

238 {
239 const MeshBase & mesh = es.get_mesh();
240
241 const unsigned int dim = mesh.mesh_dimension();
242
243 LinearImplicitSystem & system = es.get_system<LinearImplicitSystem>("Elasticity");
244
245 const unsigned int u_var = system.variable_number ("u");
246
247 const DofMap & dof_map = system.get_dof_map();
248 FEType fe_type = dof_map.variable_type(u_var);
249 std::unique_ptr<FEBase> fe (FEBase::build(dim, fe_type));
250 QGauss qrule (dim, fe_type.default_quadrature_order());
251 fe->attach_quadrature_rule (&qrule);
252
253 std::unique_ptr<FEBase> fe_face (FEBase::build(dim, fe_type));
254 QGauss qface(dim-1, fe_type.default_quadrature_order());
255 fe_face->attach_quadrature_rule (&qface);
256
257 const std::vector<Real> & JxW = fe->get_JxW();
258 const std::vector<std::vector<Real>> & phi = fe->get_phi();
259 const std::vector<std::vector<RealGradient>> & dphi = fe->get_dphi();
260
262 DenseSubMatrix<Number> Ke_var[3][3] =
263 {
267 };
268
270
271 DenseSubVector<Number> Fe_var[3] =
275
276 std::vector<dof_id_type> dof_indices;
277 std::vector<std::vector<dof_id_type>> dof_indices_var(3);
278
279 SparseMatrix<Number> & matrix = system.get_system_matrix();
280
281 for (const auto & elem : mesh.active_local_element_ptr_range())
282 {
283 dof_map.dof_indices (elem, dof_indices);
284 for (unsigned int var=0; var<3; var++)
285 dof_map.dof_indices (elem, dof_indices_var[var], var);
286
287 const unsigned int n_dofs = dof_indices.size();
288 const unsigned int n_var_dofs = dof_indices_var[0].size();
289
290 fe->reinit (elem);
291
292 Ke.resize (n_dofs, n_dofs);
293 for (unsigned int var_i=0; var_i<3; var_i++)
294 for (unsigned int var_j=0; var_j<3; var_j++)
295 Ke_var[var_i][var_j].reposition (var_i*n_var_dofs, var_j*n_var_dofs, n_var_dofs, n_var_dofs);
296
297 Fe.resize (n_dofs);
298 for (unsigned int var=0; var<3; var++)
299 Fe_var[var].reposition (var*n_var_dofs, n_var_dofs);
300
301 for (unsigned int qp=0; qp<qrule.n_points(); qp++)
302 {
303 // assemble \int_Omega C_ijkl u_k,l v_i,j \dx
304 for (unsigned int dof_i=0; dof_i<n_var_dofs; dof_i++)
305 for (unsigned int dof_j=0; dof_j<n_var_dofs; dof_j++)
306 for (unsigned int i=0; i<3; i++)
307 for (unsigned int j=0; j<3; j++)
308 for (unsigned int k=0; k<3; k++)
309 for (unsigned int l=0; l<3; l++)
310 Ke_var[i][k](dof_i,dof_j) +=
311 JxW[qp] * elasticity_tensor(i,j,k,l) * dphi[dof_j][qp](l) * dphi[dof_i][qp](j);
312
313 // assemble \int_Omega f_i v_i \dx
314 // The mesh for this example has two subdomains with ids 1
315 // and 101, and the forcing function is different on each.
316 auto f_vec = (elem->subdomain_id() == 101 ?
317 VectorValue<Number>(1., 1., 0.) :
318 VectorValue<Number>(0.36603, 1.36603, 0.));
319
320 for (unsigned int dof_i=0; dof_i<n_var_dofs; dof_i++)
321 for (unsigned int i=0; i<3; i++)
322 Fe_var[i](dof_i) += JxW[qp] * (f_vec(i) * phi[dof_i][qp]);
323 }
324
325 dof_map.constrain_element_matrix_and_vector (Ke, Fe, dof_indices);
326
327 matrix.add_matrix (Ke, dof_indices);
328 system.rhs->add_vector (Fe, dof_indices);
329 }
330 }

References libMesh::SparseMatrix< T >::add_matrix(), libMesh::NumericVector< T >::add_vector(), libMesh::FEGenericBase< OutputType >::build(), libMesh::DofMap::constrain_element_matrix_and_vector(), libMesh::FEType::default_quadrature_order(), dim, libMesh::DofMap::dof_indices(), elasticity_tensor(), es, libMesh::System::get_dof_map(), libMesh::EquationSystems::get_mesh(), libMesh::EquationSystems::get_system(), libMesh::ImplicitSystem::get_system_matrix(), mesh, libMesh::MeshBase::mesh_dimension(), libMesh::QBase::n_points(), libMesh::DenseVector< T >::resize(), libMesh::DenseMatrix< T >::resize(), libMesh::ExplicitSystem::rhs, libMesh::System::variable_number(), and libMesh::DofMap::variable_type().

◆ compute_stresses()

void LinearElasticity::compute_stresses ( )
inline

Definition at line 267 of file systems_of_equations_ex6.C.

268 {
269 const MeshBase & mesh = es.get_mesh();
270 const unsigned int dim = mesh.mesh_dimension();
271
272 LinearImplicitSystem & system = es.get_system<LinearImplicitSystem>("Elasticity");
273
274 unsigned int displacement_vars[3];
275 displacement_vars[0] = system.variable_number ("u");
276 displacement_vars[1] = system.variable_number ("v");
277 displacement_vars[2] = system.variable_number ("w");
278 const unsigned int u_var = system.variable_number ("u");
279
280 const DofMap & dof_map = system.get_dof_map();
281 FEType fe_type = dof_map.variable_type(u_var);
282 std::unique_ptr<FEBase> fe (FEBase::build(dim, fe_type));
283 QGauss qrule (dim, fe_type.default_quadrature_order());
284 fe->attach_quadrature_rule (&qrule);
285
286 const std::vector<Real> & JxW = fe->get_JxW();
287 const std::vector<std::vector<Real>> & phi = fe->get_phi();
288 const std::vector<std::vector<RealGradient>> & dphi = fe->get_dphi();
289
290 // Also, get a reference to the ExplicitSystem
291 ExplicitSystem & stress_system = es.get_system<ExplicitSystem>("StressSystem");
292 const DofMap & stress_dof_map = stress_system.get_dof_map();
293 unsigned int sigma_vars[6];
294 sigma_vars[0] = stress_system.variable_number ("sigma_00");
295 sigma_vars[1] = stress_system.variable_number ("sigma_01");
296 sigma_vars[2] = stress_system.variable_number ("sigma_02");
297 sigma_vars[3] = stress_system.variable_number ("sigma_11");
298 sigma_vars[4] = stress_system.variable_number ("sigma_12");
299 sigma_vars[5] = stress_system.variable_number ("sigma_22");
300 unsigned int vonMises_var = stress_system.variable_number ("vonMises");
301
302 // Storage for the stress dof indices on each element
303 std::vector<std::vector<dof_id_type>> dof_indices_var(system.n_vars());
304 std::vector<dof_id_type> stress_dof_indices_var;
305 std::vector<dof_id_type> vonmises_dof_indices_var;
306
307 for (const auto & elem : mesh.active_local_element_ptr_range())
308 {
309 for (unsigned int var=0; var<3; var++)
310 dof_map.dof_indices (elem, dof_indices_var[var], displacement_vars[var]);
311
312 const unsigned int n_var_dofs = dof_indices_var[0].size();
313
314 fe->reinit (elem);
315
316 std::vector<TensorValue<Number>> stress_tensor_qp(qrule.n_points());
317 for (unsigned int qp=0; qp<qrule.n_points(); qp++)
318 {
319 // Row is variable u1, u2, or u3, column is x, y, or z
320 TensorValue<Number> grad_u;
321 for (unsigned int var_i=0; var_i<3; var_i++)
322 for (unsigned int var_j=0; var_j<3; var_j++)
323 for (unsigned int j=0; j<n_var_dofs; j++)
324 grad_u(var_i,var_j) += dphi[j][qp](var_j) * system.current_solution(dof_indices_var[var_i][j]);
325
326 for (unsigned int var_i=0; var_i<3; var_i++)
327 for (unsigned int var_j=0; var_j<3; var_j++)
328 for (unsigned int k=0; k<3; k++)
329 for (unsigned int l=0; l<3; l++)
330 stress_tensor_qp[qp](var_i,var_j) += elasticity_tensor(var_i,var_j,k,l) * grad_u(k,l);
331 }
332
333 stress_dof_map.dof_indices (elem, vonmises_dof_indices_var, vonMises_var);
334 std::vector<TensorValue<Number>> elem_sigma_vec(vonmises_dof_indices_var.size());
335
336 // Below we project each component of the stress tensor onto a L2_LAGRANGE discretization.
337 // Note that this gives a discontinuous stress plot on element boundaries, which is
338 // appropriate. We then also get the von Mises stress from the projected stress tensor.
339 unsigned int stress_var_index = 0;
340 for (unsigned int var_i=0; var_i<3; var_i++)
341 for (unsigned int var_j=var_i; var_j<3; var_j++)
342 {
343 stress_dof_map.dof_indices (elem, stress_dof_indices_var, sigma_vars[stress_var_index]);
344
345 const unsigned int n_proj_dofs = stress_dof_indices_var.size();
346
347 DenseMatrix<Real> Me(n_proj_dofs, n_proj_dofs);
348 for (unsigned int qp=0; qp<qrule.n_points(); qp++)
349 {
350 for(unsigned int i=0; i<n_proj_dofs; i++)
351 for(unsigned int j=0; j<n_proj_dofs; j++)
352 {
353 Me(i,j) += JxW[qp]*(phi[i][qp]*phi[j][qp]);
354 }
355 }
356
357 DenseVector<Number> Fe(n_proj_dofs);
358 for (unsigned int qp=0; qp<qrule.n_points(); qp++)
359 for(unsigned int i=0; i<n_proj_dofs; i++)
360 {
361 Fe(i) += JxW[qp] * stress_tensor_qp[qp](var_i,var_j) * phi[i][qp];
362 }
363
364 DenseVector<Number> projected_data;
365 Me.cholesky_solve(Fe, projected_data);
366
367 for(unsigned int index=0; index<n_proj_dofs; index++)
368 {
369 dof_id_type dof_index = stress_dof_indices_var[index];
370 if ((stress_system.solution->first_local_index() <= dof_index) &&
371 (dof_index < stress_system.solution->last_local_index()))
372 stress_system.solution->set(dof_index, projected_data(index));
373
374 elem_sigma_vec[index](var_i,var_j) = projected_data(index);
375 }
376
377 stress_var_index++;
378 }
379
380 for (std::size_t index=0; index<elem_sigma_vec.size(); index++)
381 {
382 elem_sigma_vec[index](1,0) = elem_sigma_vec[index](0,1);
383 elem_sigma_vec[index](2,0) = elem_sigma_vec[index](0,2);
384 elem_sigma_vec[index](2,1) = elem_sigma_vec[index](1,2);
385
386 // Get the von Mises stress from the projected stress tensor
387 Number vonMises_value = std::sqrt(0.5*(Utility::pow<2>(elem_sigma_vec[index](0,0) - elem_sigma_vec[index](1,1)) +
388 Utility::pow<2>(elem_sigma_vec[index](1,1) - elem_sigma_vec[index](2,2)) +
389 Utility::pow<2>(elem_sigma_vec[index](2,2) - elem_sigma_vec[index](0,0)) +
390 6.*(Utility::pow<2>(elem_sigma_vec[index](0,1)) +
391 Utility::pow<2>(elem_sigma_vec[index](1,2)) +
392 Utility::pow<2>(elem_sigma_vec[index](2,0)))));
393
394 dof_id_type dof_index = vonmises_dof_indices_var[index];
395
396 if ((stress_system.solution->first_local_index() <= dof_index) &&
397 (dof_index < stress_system.solution->last_local_index()))
398 stress_system.solution->set(dof_index, vonMises_value);
399 }
400 }
401
402 // Should call close and update when we set vector entries directly
403 stress_system.solution->close();
404 stress_system.update();
405 }
Manages consistently variables, degrees of freedom, and coefficient vectors for explicit systems.
Number current_solution(const dof_id_type global_dof_number) const
Definition system.C:162
std::unique_ptr< NumericVector< Number > > solution
Data structure to hold solution values.
Definition system.h:1655
virtual void update()
Update the local values to reflect the solution on neighboring processors.
Definition system.C:498
unsigned int n_vars() const
Definition system.C:2674
This class defines a tensor in LIBMESH_DIM dimensional Real or Complex space.
uint8_t dof_id_type
Definition id_types.h:67

References libMesh::FEGenericBase< OutputType >::build(), libMesh::DenseMatrix< T >::cholesky_solve(), libMesh::System::current_solution(), libMesh::FEType::default_quadrature_order(), dim, libMesh::DofMap::dof_indices(), elasticity_tensor(), es, libMesh::System::get_dof_map(), libMesh::EquationSystems::get_mesh(), libMesh::EquationSystems::get_system(), mesh, libMesh::MeshBase::mesh_dimension(), libMesh::QBase::n_points(), libMesh::System::n_vars(), libMesh::System::solution, libMesh::System::update(), libMesh::System::variable_number(), and libMesh::DofMap::variable_type().

Referenced by main().

◆ elasticity_tensor() [1/2]

Real LinearElasticity::elasticity_tensor ( unsigned int  i,
unsigned int  j,
unsigned int  k,
unsigned int  l 
)
inline

Evaluate the fourth order tensor (C_ijkl) that relates stress to strain.

Definition at line 136 of file systems_of_equations_ex6.C.

140 {
141 // Hard code material parameters for the sake of simplicity
142 const Real poisson_ratio = 0.3;
143 const Real young_modulus = 1.;
144
145 // Define the Lame constants
146 const Real lambda_1 = (young_modulus*poisson_ratio)/((1.+poisson_ratio)*(1.-2.*poisson_ratio));
147 const Real lambda_2 = young_modulus/(2.*(1.+poisson_ratio));
148
149 return lambda_1 * kronecker_delta(i, j) * kronecker_delta(k, l) +
150 lambda_2 * (kronecker_delta(i, k) * kronecker_delta(j, l) + kronecker_delta(i, l) * kronecker_delta(j, k));
151 }
Real kronecker_delta(unsigned int i, unsigned int j)
Kronecker delta function.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real

References kronecker_delta(), and libMesh::Real.

Referenced by assemble(), and compute_stresses().

◆ elasticity_tensor() [2/2]

Real LinearElasticity::elasticity_tensor ( unsigned int  i,
unsigned int  j,
unsigned int  k,
unsigned int  l 
)
inline

Evaluate the fourth order tensor (C_ijkl) that relates stress to strain.

Definition at line 217 of file systems_of_equations_ex9.C.

221 {
222 // Hard code material parameters for the sake of simplicity
223 const Real poisson_ratio = 0.3;
224 const Real young_modulus = 1.;
225
226 // Define the Lame constants
227 const Real lambda_1 = (young_modulus*poisson_ratio)/((1.+poisson_ratio)*(1.-2.*poisson_ratio));
228 const Real lambda_2 = young_modulus/(2.*(1.+poisson_ratio));
229
230 return lambda_1 * kronecker_delta(i, j) * kronecker_delta(k, l) +
231 lambda_2 * (kronecker_delta(i, k) * kronecker_delta(j, l) + kronecker_delta(i, l) * kronecker_delta(j, k));
232 }

References kronecker_delta(), and libMesh::Real.

◆ kronecker_delta() [1/2]

Real LinearElasticity::kronecker_delta ( unsigned int  i,
unsigned int  j 
)
inline

Kronecker delta function.

Definition at line 127 of file systems_of_equations_ex6.C.

129 {
130 return i == j ? 1. : 0.;
131 }

Referenced by elasticity_tensor().

◆ kronecker_delta() [2/2]

Real LinearElasticity::kronecker_delta ( unsigned int  i,
unsigned int  j 
)
inline

Kronecker delta function.

Definition at line 208 of file systems_of_equations_ex9.C.

210 {
211 return i == j ? 1. : 0.;
212 }

Member Data Documentation

◆ es

EquationSystems & LinearElasticity::es
private

Definition at line 116 of file systems_of_equations_ex6.C.

Referenced by assemble(), and compute_stresses().


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