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

Public Member Functions

 LargeDeformationElasticity (EquationSystems &es_in)
 
Real kronecker_delta (unsigned int i, unsigned int j)
 Kronecker delta function.
 
Real elasticity_tensor (Real young_modulus, Real poisson_ratio, unsigned int i, unsigned int j, unsigned int k, unsigned int l)
 Evaluate the fourth order tensor (C_ijkl) that relates stress to strain.
 
virtual void jacobian (const NumericVector< Number > &soln, SparseMatrix< Number > &jacobian, NonlinearImplicitSystem &)
 Evaluate the Jacobian of the nonlinear system.
 
virtual void residual (const NumericVector< Number > &soln, NumericVector< Number > &residual, NonlinearImplicitSystem &)
 Evaluate the residual of the nonlinear system.
 
void compute_stresses ()
 Compute the Cauchy stress for the current solution.
 

Private Attributes

EquationSystems & es
 

Detailed Description

Definition at line 87 of file systems_of_equations_ex7.C.

Constructor & Destructor Documentation

◆ LargeDeformationElasticity()

LargeDeformationElasticity::LargeDeformationElasticity ( EquationSystems &  es_in)
inline

Definition at line 95 of file systems_of_equations_ex7.C.

95 :
96 es(es_in)
97 {}

Member Function Documentation

◆ compute_stresses()

void LargeDeformationElasticity::compute_stresses ( )
inline

Compute the Cauchy stress for the current solution.

Definition at line 381 of file systems_of_equations_ex7.C.

382 {
383 const Real young_modulus = es.parameters.get<Real>("young_modulus");
384 const Real poisson_ratio = es.parameters.get<Real>("poisson_ratio");
385
386 const MeshBase & mesh = es.get_mesh();
387 const unsigned int dim = mesh.mesh_dimension();
388
389 NonlinearImplicitSystem & system =
390 es.get_system<NonlinearImplicitSystem>("NonlinearElasticity");
391
392 unsigned int displacement_vars[] = {
393 system.variable_number ("u"),
394 system.variable_number ("v"),
395 system.variable_number ("w")};
396 const unsigned int u_var = system.variable_number ("u");
397
398 const DofMap & dof_map = system.get_dof_map();
399 FEType fe_type = dof_map.variable_type(u_var);
400 std::unique_ptr<FEBase> fe (FEBase::build(dim, fe_type));
401 QGauss qrule (dim, fe_type.default_quadrature_order());
402 fe->attach_quadrature_rule (&qrule);
403
404 const std::vector<Real> & JxW = fe->get_JxW();
405 const std::vector<std::vector<RealGradient>> & dphi = fe->get_dphi();
406
407 // Also, get a reference to the ExplicitSystem
408 ExplicitSystem & stress_system = es.get_system<ExplicitSystem>("StressSystem");
409 const DofMap & stress_dof_map = stress_system.get_dof_map();
410 unsigned int sigma_vars[] = {
411 stress_system.variable_number ("sigma_00"),
412 stress_system.variable_number ("sigma_01"),
413 stress_system.variable_number ("sigma_02"),
414 stress_system.variable_number ("sigma_11"),
415 stress_system.variable_number ("sigma_12"),
416 stress_system.variable_number ("sigma_22")};
417
418 // Storage for the stress dof indices on each element
419 std::vector<std::vector<dof_id_type>> dof_indices_var(system.n_vars());
420 std::vector<dof_id_type> stress_dof_indices_var;
421
422 // To store the stress tensor on each element
423 TensorValue<Number> elem_avg_stress_tensor;
424
425 for (const auto & elem : mesh.active_local_element_ptr_range())
426 {
427 for (unsigned int var=0; var<3; var++)
428 dof_map.dof_indices (elem, dof_indices_var[var], displacement_vars[var]);
429
430 const unsigned int n_var_dofs = dof_indices_var[0].size();
431
432 fe->reinit (elem);
433
434 // clear the stress tensor
435 elem_avg_stress_tensor.zero();
436
437 for (unsigned int qp=0; qp<qrule.n_points(); qp++)
438 {
439 TensorValue<Number> grad_u;
440 // Row is variable u1, u2, or u3, column is x, y, or z
441 for (unsigned int var_i=0; var_i<3; var_i++)
442 for (unsigned int var_j=0; var_j<3; var_j++)
443 for (unsigned int j=0; j<n_var_dofs; j++)
444 grad_u(var_i,var_j) += dphi[j][qp](var_j) * system.current_solution(dof_indices_var[var_i][j]);
445
446 TensorValue<Number> strain_tensor;
447 for (unsigned int i=0; i<3; i++)
448 for (unsigned int j=0; j<3; j++)
449 {
450 strain_tensor(i,j) += 0.5 * (grad_u(i,j) + grad_u(j,i));
451
452 for (unsigned int k=0; k<3; k++)
453 strain_tensor(i,j) += 0.5 * grad_u(k,i)*grad_u(k,j);
454 }
455
456 // Define the deformation gradient
457 auto F = grad_u;
458 for (unsigned int var=0; var<3; var++)
459 F(var, var) += 1.;
460
461 TensorValue<Number> stress_tensor;
462 for (unsigned int i=0; i<3; i++)
463 for (unsigned int j=0; j<3; j++)
464 for (unsigned int k=0; k<3; k++)
465 for (unsigned int l=0; l<3; l++)
466 stress_tensor(i,j) +=
467 elasticity_tensor(young_modulus, poisson_ratio, i, j, k, l) * strain_tensor(k, l);
468
469 // stress_tensor now holds the second Piola-Kirchoff stress (PK2) at point qp.
470 // However, in this example we want to compute the Cauchy stress which is given by
471 // 1/det(F) * F * PK2 * F^T, hence we now apply this transformation.
472 stress_tensor = 1. / F.det() * F * stress_tensor * F.transpose();
473
474 // We want to plot the average Cauchy stress on each element, hence
475 // we integrate stress_tensor
476 elem_avg_stress_tensor.add_scaled(stress_tensor, JxW[qp]);
477 }
478
479 // Get the average stress per element by dividing by volume
480 elem_avg_stress_tensor /= elem->volume();
481
482 // load elem_sigma data into stress_system
483 unsigned int stress_var_index = 0;
484 for (unsigned int i=0; i<3; i++)
485 for (unsigned int j=i; j<3; j++)
486 {
487 stress_dof_map.dof_indices (elem, stress_dof_indices_var, sigma_vars[stress_var_index]);
488
489 // We are using CONSTANT MONOMIAL basis functions, hence we only need to get
490 // one dof index per variable
491 dof_id_type dof_index = stress_dof_indices_var[0];
492
493 if ((stress_system.solution->first_local_index() <= dof_index) &&
494 (dof_index < stress_system.solution->last_local_index()))
495 stress_system.solution->set(dof_index, elem_avg_stress_tensor(i,j));
496
497 stress_var_index++;
498 }
499 }
500
501 // Should call close and update when we set vector entries directly
502 stress_system.solution->close();
503 stress_system.update();
504 }
dim
unsigned int dim
Definition adaptivity_ex3.C:114
LargeDeformationElasticity::elasticity_tensor
Real elasticity_tensor(Real young_modulus, Real poisson_ratio, unsigned int i, unsigned int j, unsigned int k, unsigned int l)
Evaluate the fourth order tensor (C_ijkl) that relates stress to strain.
Definition systems_of_equations_ex7.C:111
libMesh::DofMap
This class handles the numbering of degrees of freedom on a mesh.
Definition dof_map.h:181
libMesh::DofMap::dof_indices
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
Definition dof_map.C:2201
libMesh::DofMap::variable_type
const FEType & variable_type(const unsigned int i) const
Definition dof_map.h:2388
libMesh::EquationSystems::get_mesh
const MeshBase & get_mesh() const
Definition equation_systems.h:665
libMesh::EquationSystems::parameters
Parameters parameters
Data structure holding arbitrary parameters.
Definition equation_systems.h:608
libMesh::EquationSystems::get_system
const T_sys & get_system(std::string_view name) const
Definition equation_systems.h:776
libMesh::ExplicitSystem
Manages consistently variables, degrees of freedom, and coefficient vectors for explicit systems.
Definition explicit_system.h:49
libMesh::FEGenericBase::build
static std::unique_ptr< FEGenericBase > build(const unsigned int dim, const FEType &type)
Builds a specific finite element type.
libMesh::FEType
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition fe_type.h:197
libMesh::FEType::default_quadrature_order
Order default_quadrature_order() const
Definition fe_type.h:415
libMesh::MeshBase
This is the MeshBase class.
Definition mesh_base.h:81
libMesh::MeshBase::mesh_dimension
unsigned int mesh_dimension() const
Definition mesh_base.C:430
libMesh::NonlinearImplicitSystem
Manages consistently variables, degrees of freedom, coefficient vectors, matrices and non-linear solv...
Definition nonlinear_implicit_system.h:55
libMesh::Parameters::get
const T & get(std::string_view) const
Definition parameters.h:451
libMesh::QGauss
This class implements specific orders of Gauss quadrature.
Definition quadrature_gauss.h:40
libMesh::System::current_solution
Number current_solution(const dof_id_type global_dof_number) const
Definition system.C:162
libMesh::System::solution
std::unique_ptr< NumericVector< Number > > solution
Data structure to hold solution values.
Definition system.h:1655
libMesh::System::update
virtual void update()
Update the local values to reflect the solution on neighboring processors.
Definition system.C:498
libMesh::System::variable_number
unsigned int variable_number(std::string_view var) const
Definition system.C:1398
libMesh::System::n_vars
unsigned int n_vars() const
Definition system.C:2674
libMesh::System::get_dof_map
const DofMap & get_dof_map() const
Definition system.h:2417
libMesh::TensorValue
This class defines a tensor in LIBMESH_DIM dimensional Real or Complex space.
Definition tensor_value.h:45
libMesh::TypeTensor::add_scaled
void add_scaled(const TypeTensor< T2 > &, const T &)
Add a scaled tensor to this tensor without creating a temporary.
Definition type_tensor.h:858
libMesh::TypeTensor::zero
void zero()
Set all entries of the tensor to 0.
Definition type_tensor.h:1346
libMesh::TypeTensor::transpose
TypeTensor< T > transpose() const
Definition type_tensor.h:1057
mesh
MeshBase & mesh
Definition mesh_communication.C:1460
libMesh::dof_id_type
uint8_t dof_id_type
Definition id_types.h:67
libMesh::Real
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Definition libmesh_common.h:144

References libMesh::TypeTensor< T >::add_scaled(), libMesh::FEGenericBase< OutputType >::build(), libMesh::System::current_solution(), libMesh::FEType::default_quadrature_order(), libMesh::TypeTensor< T >::det(), dim, libMesh::DofMap::dof_indices(), elasticity_tensor(), es, libMesh::Parameters::get(), 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::EquationSystems::parameters, libMesh::Real, libMesh::System::solution, libMesh::TypeTensor< T >::transpose(), libMesh::System::update(), libMesh::System::variable_number(), libMesh::DofMap::variable_type(), and libMesh::TypeTensor< T >::zero().

Referenced by main().

◆ elasticity_tensor()

Real LargeDeformationElasticity::elasticity_tensor ( Real  young_modulus,
Real  poisson_ratio,
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 111 of file systems_of_equations_ex7.C.

117 {
118 // Define the Lame constants
119 const Real lambda_1 = (young_modulus*poisson_ratio)/((1.+poisson_ratio)*(1.-2.*poisson_ratio));
120 const Real lambda_2 = young_modulus/(2.*(1.+poisson_ratio));
121
122 return lambda_1 * kronecker_delta(i,j) * kronecker_delta(k,l) +
123 lambda_2 * (kronecker_delta(i,k) * kronecker_delta(j,l) + kronecker_delta(i,l) * kronecker_delta(j,k));
124 }
LargeDeformationElasticity::kronecker_delta
Real kronecker_delta(unsigned int i, unsigned int j)
Kronecker delta function.
Definition systems_of_equations_ex7.C:102

References kronecker_delta(), and libMesh::Real.

Referenced by compute_stresses(), jacobian(), and residual().

◆ jacobian()

virtual void LargeDeformationElasticity::jacobian ( const NumericVector< Number > &  soln,
SparseMatrix< Number > &  jacobian,
NonlinearImplicitSystem &   
)
inlinevirtual

Evaluate the Jacobian of the nonlinear system.

Definition at line 130 of file systems_of_equations_ex7.C.

133 {
134 const Real young_modulus = es.parameters.get<Real>("young_modulus");
135 const Real poisson_ratio = es.parameters.get<Real>("poisson_ratio");
136
137 const MeshBase & mesh = es.get_mesh();
138 const unsigned int dim = mesh.mesh_dimension();
139
140 NonlinearImplicitSystem & system =
141 es.get_system<NonlinearImplicitSystem>("NonlinearElasticity");
142
143 const unsigned int u_var = system.variable_number ("u");
144
145 const DofMap & dof_map = system.get_dof_map();
146
147 FEType fe_type = dof_map.variable_type(u_var);
148 std::unique_ptr<FEBase> fe (FEBase::build(dim, fe_type));
149 QGauss qrule (dim, fe_type.default_quadrature_order());
150 fe->attach_quadrature_rule (&qrule);
151
152 std::unique_ptr<FEBase> fe_face (FEBase::build(dim, fe_type));
153 QGauss qface (dim-1, fe_type.default_quadrature_order());
154 fe_face->attach_quadrature_rule (&qface);
155
156 const std::vector<Real> & JxW = fe->get_JxW();
157 const std::vector<std::vector<RealGradient>> & dphi = fe->get_dphi();
158
159 DenseMatrix<Number> Ke;
160 DenseSubMatrix<Number> Ke_var[3][3] =
161 {
162 {DenseSubMatrix<Number>(Ke), DenseSubMatrix<Number>(Ke), DenseSubMatrix<Number>(Ke)},
163 {DenseSubMatrix<Number>(Ke), DenseSubMatrix<Number>(Ke), DenseSubMatrix<Number>(Ke)},
164 {DenseSubMatrix<Number>(Ke), DenseSubMatrix<Number>(Ke), DenseSubMatrix<Number>(Ke)}
165 };
166
167 std::vector<dof_id_type> dof_indices;
168 std::vector<std::vector<dof_id_type>> dof_indices_var(3);
169
170 jacobian.zero();
171
172 for (const auto & elem : mesh.active_local_element_ptr_range())
173 {
174 dof_map.dof_indices (elem, dof_indices);
175 for (unsigned int var=0; var<3; var++)
176 dof_map.dof_indices (elem, dof_indices_var[var], var);
177
178 const unsigned int n_dofs = dof_indices.size();
179 const unsigned int n_var_dofs = dof_indices_var[0].size();
180
181 fe->reinit (elem);
182
183 Ke.resize (n_dofs, n_dofs);
184 for (unsigned int var_i=0; var_i<3; var_i++)
185 for (unsigned int var_j=0; var_j<3; var_j++)
186 Ke_var[var_i][var_j].reposition (var_i*n_var_dofs, var_j*n_var_dofs, n_var_dofs, n_var_dofs);
187
188 for (unsigned int qp=0; qp<qrule.n_points(); qp++)
189 {
190 TensorValue<Number> grad_u;
191 for (unsigned int var_i=0; var_i<3; var_i++)
192 {
193 // Row is variable u1, u2, or u3, column is x, y, or z
194 for (unsigned int var_j=0; var_j<3; var_j++)
195 for (unsigned int j=0; j<n_var_dofs; j++)
196 grad_u(var_i,var_j) += dphi[j][qp](var_j)*soln(dof_indices_var[var_i][j]);
197 }
198
199 TensorValue<Number> strain_tensor;
200 for (unsigned int i=0; i<3; i++)
201 for (unsigned int j=0; j<3; j++)
202 {
203 strain_tensor(i,j) += 0.5 * (grad_u(i,j) + grad_u(j,i));
204
205 for (unsigned int k=0; k<3; k++)
206 strain_tensor(i,j) += 0.5 * grad_u(k,i)*grad_u(k,j);
207 }
208
209 // Define the deformation gradient
210 auto F = grad_u;
211 for (unsigned int var=0; var<3; var++)
212 F(var, var) += 1.;
213
214 TensorValue<Number> stress_tensor;
215
216 for (unsigned int i=0; i<3; i++)
217 for (unsigned int j=0; j<3; j++)
218 for (unsigned int k=0; k<3; k++)
219 for (unsigned int l=0; l<3; l++)
220 stress_tensor(i,j) +=
221 elasticity_tensor(young_modulus, poisson_ratio, i, j, k, l) * strain_tensor(k, l);
222
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 {
226 for (unsigned int i=0; i<3; i++)
227 for (unsigned int j=0; j<3; j++)
228 for (unsigned int m=0; m<3; m++)
229 Ke_var[i][i](dof_i,dof_j) += JxW[qp] *
230 (-dphi[dof_j][qp](m) * stress_tensor(m,j) * dphi[dof_i][qp](j));
231
232 for (unsigned int i=0; i<3; i++)
233 for (unsigned int j=0; j<3; j++)
234 for (unsigned int k=0; k<3; k++)
235 for (unsigned int l=0; l<3; l++)
236 {
237 Number FxC_ijkl = 0.;
238 for (unsigned int m=0; m<3; m++)
239 FxC_ijkl += F(i,m) * elasticity_tensor(young_modulus, poisson_ratio, m, j, k, l);
240
241 Ke_var[i][k](dof_i,dof_j) += JxW[qp] *
242 (-0.5 * FxC_ijkl * dphi[dof_j][qp](l) * dphi[dof_i][qp](j));
243
244 Ke_var[i][l](dof_i,dof_j) += JxW[qp] *
245 (-0.5 * FxC_ijkl * dphi[dof_j][qp](k) * dphi[dof_i][qp](j));
246
247 for (unsigned int n=0; n<3; n++)
248 Ke_var[i][n](dof_i,dof_j) += JxW[qp] *
249 (-0.5 * FxC_ijkl * (dphi[dof_j][qp](k) * grad_u(n,l) + dphi[dof_j][qp](l) * grad_u(n,k)) * dphi[dof_i][qp](j));
250 }
251 }
252 }
253
254 dof_map.constrain_element_matrix (Ke, dof_indices);
255 jacobian.add_matrix (Ke, dof_indices);
256 }
257 }
LargeDeformationElasticity::jacobian
virtual void jacobian(const NumericVector< Number > &soln, SparseMatrix< Number > &jacobian, NonlinearImplicitSystem &)
Evaluate the Jacobian of the nonlinear system.
Definition systems_of_equations_ex7.C:130
libMesh::DenseMatrix
Defines a dense matrix for use in Finite Element-type computations.
Definition dense_matrix.h:72
libMesh::DenseMatrix::resize
void resize(const unsigned int new_m, const unsigned int new_n)
Resizes the matrix to the specified size and calls zero().
Definition dense_matrix.h:895
libMesh::DenseSubMatrix
Defines a dense submatrix for use in Finite Element-type computations.
Definition dense_submatrix.h:46
libMesh::DofMap::constrain_element_matrix
void constrain_element_matrix(DenseMatrix< Number > &matrix, std::vector< dof_id_type > &elem_dofs, bool asymmetric_constraint_rows=true) const
Constrains the element matrix.
Definition dof_map.h:2485

References libMesh::FEGenericBase< OutputType >::build(), libMesh::DofMap::constrain_element_matrix(), libMesh::FEType::default_quadrature_order(), dim, libMesh::DofMap::dof_indices(), elasticity_tensor(), es, libMesh::Parameters::get(), libMesh::System::get_dof_map(), libMesh::EquationSystems::get_mesh(), libMesh::EquationSystems::get_system(), jacobian(), mesh, libMesh::MeshBase::mesh_dimension(), libMesh::QBase::n_points(), libMesh::EquationSystems::parameters, libMesh::Real, libMesh::DenseMatrix< T >::resize(), libMesh::System::variable_number(), and libMesh::DofMap::variable_type().

Referenced by jacobian().

◆ kronecker_delta()

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

Kronecker delta function.

Definition at line 102 of file systems_of_equations_ex7.C.

104 {
105 return i == j ? 1. : 0.;
106 }

Referenced by elasticity_tensor().

◆ residual()

virtual void LargeDeformationElasticity::residual ( const NumericVector< Number > &  soln,
NumericVector< Number > &  residual,
NonlinearImplicitSystem &   
)
inlinevirtual

Evaluate the residual of the nonlinear system.

Definition at line 262 of file systems_of_equations_ex7.C.

265 {
266 const Real young_modulus = es.parameters.get<Real>("young_modulus");
267 const Real poisson_ratio = es.parameters.get<Real>("poisson_ratio");
268 const Real forcing_magnitude = es.parameters.get<Real>("forcing_magnitude");
269
270 const MeshBase & mesh = es.get_mesh();
271 const unsigned int dim = mesh.mesh_dimension();
272
273 NonlinearImplicitSystem & system =
274 es.get_system<NonlinearImplicitSystem>("NonlinearElasticity");
275
276 const unsigned int u_var = system.variable_number ("u");
277
278 const DofMap & dof_map = system.get_dof_map();
279
280 FEType fe_type = dof_map.variable_type(u_var);
281 std::unique_ptr<FEBase> fe (FEBase::build(dim, fe_type));
282 QGauss qrule (dim, fe_type.default_quadrature_order());
283 fe->attach_quadrature_rule (&qrule);
284
285 std::unique_ptr<FEBase> fe_face (FEBase::build(dim, fe_type));
286 QGauss qface (dim-1, fe_type.default_quadrature_order());
287 fe_face->attach_quadrature_rule (&qface);
288
289 const std::vector<Real> & JxW = fe->get_JxW();
290 const std::vector<std::vector<Real>> & phi = fe->get_phi();
291 const std::vector<std::vector<RealGradient>> & dphi = fe->get_dphi();
292
293 DenseVector<Number> Re;
294
295 DenseSubVector<Number> Re_var[3] =
296 {DenseSubVector<Number>(Re),
297 DenseSubVector<Number>(Re),
298 DenseSubVector<Number>(Re)};
299
300 std::vector<dof_id_type> dof_indices;
301 std::vector<std::vector<dof_id_type>> dof_indices_var(3);
302
303 residual.zero();
304
305 for (const auto & elem : mesh.active_local_element_ptr_range())
306 {
307 dof_map.dof_indices (elem, dof_indices);
308 for (unsigned int var=0; var<3; var++)
309 dof_map.dof_indices (elem, dof_indices_var[var], var);
310
311 const unsigned int n_dofs = dof_indices.size();
312 const unsigned int n_var_dofs = dof_indices_var[0].size();
313
314 fe->reinit (elem);
315
316 Re.resize (n_dofs);
317 for (unsigned int var=0; var<3; var++)
318 Re_var[var].reposition (var*n_var_dofs, n_var_dofs);
319
320 for (unsigned int qp=0; qp<qrule.n_points(); qp++)
321 {
322 TensorValue<Number> grad_u;
323 for (unsigned int var_i=0; var_i<3; var_i++)
324 {
325 // Row is variable u, v, or w column is x, y, or z
326 for (unsigned int var_j=0; var_j<3; var_j++)
327 for (unsigned int j=0; j<n_var_dofs; j++)
328 grad_u(var_i,var_j) += dphi[j][qp](var_j)*soln(dof_indices_var[var_i][j]);
329 }
330
331 TensorValue<Number> strain_tensor;
332 for (unsigned int i=0; i<3; i++)
333 for (unsigned int j=0; j<3; j++)
334 {
335 strain_tensor(i,j) += 0.5 * (grad_u(i,j) + grad_u(j,i));
336
337 for (unsigned int k=0; k<3; k++)
338 strain_tensor(i,j) += 0.5 * grad_u(k,i)*grad_u(k,j);
339 }
340
341 // Define the deformation gradient
342 auto F = grad_u;
343 for (unsigned int var=0; var<3; var++)
344 F(var, var) += 1.;
345
346 TensorValue<Number> stress_tensor;
347
348 for (unsigned int i=0; i<3; i++)
349 for (unsigned int j=0; j<3; j++)
350 for (unsigned int k=0; k<3; k++)
351 for (unsigned int l=0; l<3; l++)
352 stress_tensor(i,j) +=
353 elasticity_tensor(young_modulus, poisson_ratio, i, j, k, l) * strain_tensor(k,l);
354
355 VectorValue<Number> f_vec(0., 0., -forcing_magnitude);
356
357 for (unsigned int dof_i=0; dof_i<n_var_dofs; dof_i++)
358 for (unsigned int i=0; i<3; i++)
359 {
360 for (unsigned int j=0; j<3; j++)
361 {
362 Number FxStress_ij = 0.;
363 for (unsigned int m=0; m<3; m++)
364 FxStress_ij += F(i,m) * stress_tensor(m,j);
365
366 Re_var[i](dof_i) += JxW[qp] * (-FxStress_ij * dphi[dof_i][qp](j));
367 }
368
369 Re_var[i](dof_i) += JxW[qp] * (f_vec(i) * phi[dof_i][qp]);
370 }
371 }
372
373 dof_map.constrain_element_vector (Re, dof_indices);
374 residual.add_vector (Re, dof_indices);
375 }
376 }
LargeDeformationElasticity::residual
virtual void residual(const NumericVector< Number > &soln, NumericVector< Number > &residual, NonlinearImplicitSystem &)
Evaluate the residual of the nonlinear system.
Definition systems_of_equations_ex7.C:262
libMesh::DenseSubVector
Defines a dense subvector for use in finite element computations.
Definition dense_subvector.h:44
libMesh::DenseVector
Defines a dense vector for use in Finite Element-type computations.
Definition dense_vector.h:59
libMesh::DenseVector::resize
void resize(const unsigned int n)
Resize the vector.
Definition dense_vector.h:396
libMesh::DofMap::constrain_element_vector
void constrain_element_vector(DenseVector< Number > &rhs, std::vector< dof_id_type > &dofs, bool asymmetric_constraint_rows=true) const
Constrains the element vector.
Definition dof_map.h:2494
libMesh::VectorValue
This class defines a vector in LIBMESH_DIM dimensional Real or Complex space.
Definition vector_value.h:48

References libMesh::FEGenericBase< OutputType >::build(), libMesh::DofMap::constrain_element_vector(), libMesh::FEType::default_quadrature_order(), dim, libMesh::DofMap::dof_indices(), elasticity_tensor(), es, libMesh::Parameters::get(), libMesh::System::get_dof_map(), libMesh::EquationSystems::get_mesh(), libMesh::EquationSystems::get_system(), mesh, libMesh::MeshBase::mesh_dimension(), libMesh::QBase::n_points(), libMesh::EquationSystems::parameters, libMesh::Real, residual(), libMesh::DenseVector< T >::resize(), libMesh::System::variable_number(), and libMesh::DofMap::variable_type().

Referenced by residual().

Member Data Documentation

◆ es

EquationSystems& LargeDeformationElasticity::es
private

Definition at line 91 of file systems_of_equations_ex7.C.

Referenced by compute_stresses(), jacobian(), and residual().

The documentation for this class was generated from the following file:
Generated on Fri Aug 21 2026 13:32:35 for libMesh by doxygen 1.9.8