MFEM v4.8.0
Finite element discretization library
Loading...
Searching...
No Matches
bilininteg_elasticity_pa.cpp
Go to the documentation of this file.
1// Copyright (c) 2010-2025, Lawrence Livermore National Security, LLC. Produced
2// at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3// LICENSE and NOTICE for details. LLNL-CODE-806117.
4//
5// This file is part of the MFEM library. For more information and source code
6// availability visit https://mfem.org.
7//
8// MFEM is free software; you can redistribute it and/or modify it under the
9// terms of the BSD-3 license. We welcome feedback and contributions, see file
10// CONTRIBUTING.md for details.
11
12#include "../bilininteg.hpp"
13#include "../gridfunc.hpp"
14#include "../qfunction.hpp"
16
17namespace mfem
18{
19
20void ElasticityIntegrator::SetUpQuadratureSpaceAndCoefficients(
21 const FiniteElementSpace &fes)
22{
23 if (IntRule == nullptr)
24 {
25 // This is where it's assumed that all elements are the same.
26 const auto &T = *fes.GetMesh()->GetTypicalElementTransformation();
27 int quad_order = 2 * T.OrderGrad(fes.GetTypicalFE());
28 IntRule = &IntRules.Get(T.GetGeometryType(), quad_order);
29 }
30
31 Mesh &mesh = *fespace->GetMesh();
32
33 q_space.reset(new QuadratureSpace(mesh, *IntRule));
34 lambda_quad.reset(new CoefficientVector(lambda, *q_space,
36 mu_quad.reset(new CoefficientVector(mu, *q_space, CoefficientStorage::FULL));
37 q_vec.reset(new QuadratureFunction(*q_space, vdim*vdim));
38}
39
41{
42 MFEM_VERIFY(fes.GetOrdering() == Ordering::byNODES,
43 "Elasticity PA only implemented for byNODES ordering.");
44
45 fespace = &fes;
46 Mesh &mesh = *fespace->GetMesh();
47 MFEM_VERIFY(fespace->GetVDim() == mesh.Dimension(), "");
48 vdim = fespace->GetVDim();
49 ndofs = fespace->GetTypicalFE()->GetDof();
50
51 SetUpQuadratureSpaceAndCoefficients(fes);
52
53 auto ordering = GetEVectorOrdering(*fespace);
54 auto mode = ordering == ElementDofOrdering::NATIVE ? DofToQuad::FULL :
56 maps = &fespace->GetTypicalFE()->GetDofToQuad(*IntRule, mode);
58}
59
61{
62 q_vec->SetVDim(vdim*vdim*vdim*vdim);
63 internal::ElasticityAssembleDiagonalPA(vdim, ndofs, *lambda_quad, *mu_quad,
64 *geom, *maps, *q_vec, diag);
65}
66
68{
69 internal::ElasticityAddMultPA(vdim, ndofs, *fespace, *lambda_quad, *mu_quad,
70 *geom, *maps, x, *q_vec, y);
71}
72
74{
75 AddMultPA(x, y); // Operator is symmetric
76}
77
79{
80 fespace = &fes;
81
82 // Avoid projecting the coefficients more than once. If the coefficients
83 // change, the parent ElasticityIntegrator must be reassembled.
84 if (!parent.q_space)
85 {
86 parent.SetUpQuadratureSpaceAndCoefficients(fes);
87 }
88 else
89 {
90 IntRule = parent.IntRule;
91 }
92
93 auto ordering = GetEVectorOrdering(*fespace);
94 auto mode = ordering == ElementDofOrdering::NATIVE ? DofToQuad::FULL :
96 geom = fes.GetMesh()->GetGeometricFactors(*IntRule,
98 maps = &fespace->GetTypicalFE()->GetDofToQuad(*IntRule, mode);
99}
100
102{
103 internal::ElasticityComponentAddMultPA(
104 parent.vdim, parent.ndofs, *fespace, *parent.lambda_quad, *parent.mu_quad,
105 *geom, *maps, x, *parent.q_vec, y, i_block, j_block);
106}
107
109 Vector &y) const
110{
111 // Each block in the operator is symmetric, so we can just switch the roles
112 // of i_block and j_block
113 internal::ElasticityComponentAddMultPA(
114 parent.vdim, parent.ndofs, *fespace, *parent.lambda_quad, *parent.mu_quad,
115 *geom, *maps, x, *parent.q_vec, y, j_block, i_block);
116}
117
118} // namespace mfem
Header for small strain, isotropic, linear elasticity kernels.
@ FULL
Full multidimensional representation which does not use tensor product structure. The ordering of the...
Definition fe_base.hpp:158
@ LEXICOGRAPHIC_FULL
Full multidimensional representation which does not use tensor product structure. The ordering of the...
Definition fe_base.hpp:170
void AssemblePA(const FiniteElementSpace &fes) override
Method defining partial assembly.
void AddMultTransposePA(const Vector &x, Vector &y) const override
Method for partially assembled transposed action.
void AddMultPA(const Vector &x, Vector &y) const override
Method for partially assembled action.
void AddMultTransposePA(const Vector &x, Vector &y) const override
Method for partially assembled transposed action.
void AddMultPA(const Vector &x, Vector &y) const override
Method for partially assembled action.
void AssemblePA(const FiniteElementSpace &fes) override
Method defining partial assembly.
void AssembleDiagonalPA(Vector &diag) override
Assemble diagonal and add it to Vector diag.
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition fespace.hpp:244
Ordering::Type GetOrdering() const
Return the ordering method.
Definition fespace.hpp:876
Mesh * GetMesh() const
Returns the mesh.
Definition fespace.hpp:679
int GetVDim() const
Returns the vector dimension of the finite element space.
Definition fespace.hpp:841
const FiniteElement * GetTypicalFE() const
Return GetFE(0) if the local mesh is not empty; otherwise return a typical FE based on the Geometry t...
Definition fespace.cpp:3871
virtual const DofToQuad & GetDofToQuad(const IntegrationRule &ir, DofToQuad::Mode mode) const
Return a DofToQuad structure corresponding to the given IntegrationRule using the given DofToQuad::Mo...
Definition fe_base.cpp:365
int GetDof() const
Returns the number of degrees of freedom in the finite element.
Definition fe_base.hpp:334
const IntegrationRule & Get(int GeomType, int Order)
Returns an integration rule for given GeomType and Order.
const IntegrationRule * IntRule
Mesh data type.
Definition mesh.hpp:64
int Dimension() const
Dimension of the reference space used within the elements.
Definition mesh.hpp:1216
const GeometricFactors * GetGeometricFactors(const IntegrationRule &ir, const int flags, MemoryType d_mt=MemoryType::DEFAULT)
Return the mesh geometric factors corresponding to the given integration rule.
Definition mesh.cpp:880
Vector data type.
Definition vector.hpp:82
@ FULL
Store the coefficient as a full QuadratureFunction.
ElementDofOrdering GetEVectorOrdering(const FiniteElementSpace &fes)
Return LEXICOGRAPHIC if mesh contains only one topology and the elements are tensor elements,...
Definition fespace.cpp:4596
@ NATIVE
Native ordering as defined by the FiniteElement.
IntegrationRules IntRules(0, Quadrature1D::GaussLegendre)
A global object with all integration rules (defined in intrules.cpp)
Definition intrules.hpp:492