Is it possible to get second order accuracy for hexahedral finite elements with fewer than 8 Gauss points without introducing unphysical modes? A single central Gauss point introduces an unphysical shearing mode, and the standard symmetric arrangement of 8 Gauss points is expensive compared to tetrahedral discretizations.

Edit: Someone asked for equations. The equations I'm interested in are nonlinear elasticity, either dynamic or quasistatic. The quasistatic equations are

$$\nabla \cdot P\left(\nabla \phi \right) = 0$$

where $\phi : \Omega \to \mathbf{R}^3$, $\Omega \subset \mathbf{R}^3$, and $P : \mathbf{R}^{3 \times 3} \to \mathbf{R}^{3 \times 3}$ is a hyperelastic first Piola-Kirchoff stress function. An simple example is compressible neo-Hookean, where $$ P(F) = \mu (F - F^{-T}) +\lambda F^{-T} \log \det F$$

  • $\begingroup$ What exactly are you simulating? $\endgroup$ – Dan Jan 3 '12 at 5:49
  • $\begingroup$ Linear elasticity at the moment, but the question is about nonlinear elasticity in general. $\endgroup$ – Geoffrey Irving Jan 3 '12 at 7:41
  • 1
    $\begingroup$ You should probably include the equations you're interested in, since the definition of "unphysical" depends on them. Or at least define precisely the space of functions that are "physical". $\endgroup$ – David Ketcheson Jan 3 '12 at 11:21
  • $\begingroup$ Equations added. $\endgroup$ – Geoffrey Irving Jan 3 '12 at 17:55
  • 1
    $\begingroup$ With dPhi/dx, do you mean the gradient? $\endgroup$ – Wolfgang Bangerth Jan 4 '12 at 1:34

As far as finite element solid mechanics simulations are concerned, you can't use less than 8 quadrature points without using stabilization forces. In case of incompressible material (your case), the best solution for accuracy purpose is to use mixed formulation. You can refer to the book by Simo and Hughes : http://books.google.fr/books/about/Computational_inelasticity.html?hl=fr&id=ftL2AJL8OPYC.


It's relatively obvious that you can't in general get away with fewer quadrature points per cell than there are degrees of freedom. In the case of trilinear elements on a 3d hexahedron, there are 8 degrees of freedom (one per vertex) so the minimum number of quadrature points would be eight as well.

(Now, as for the argument why that is so, here's a sketch: We know that a finite element is only good if it satisfies the "unisolvency" condition, i.e., there needs to be a single solution. The problem with too few quadrature points is that the bilinear form you get using quadrature instead of the exact integral can't distinguish between different shape functions any more. To give an example, think of forming the mass matrix in 1d when using a discontinuous element with basis {1,x} on the reference element [-1,1]. Let's assume that we have a mesh with only one element, which is also the interval [-1,1]. The exact mass matrix is, of course, [[2,0],[0,2/3]], which is invertible as should be. On the other hand, if you approximate the integral using the midpoint rule the mass matrix would be [[2,0],[0,0]], which is not invertible and consequently completely useless. The reason is that a one-point quadrature formula can't distinguish between all linear functions (part of the trial space) that have the same value at the quadrature point; in other words, for the midpoint rule, the shape function 'x' is the same as the function '0' is the same as the function '-x'. In other words, while the trial space has dimension 2 with exact integrals, for the midpoint rule the space has dimension 1, even though there are two degrees of freedom -- that's the definition of a space that's not unisolvent.)

  • $\begingroup$ I think Geoff's question is more subtle. For continuous finite element spaces on tetrahedra in well-shaped domains (e.g. without isolated elements), you can get away with single-point quadratures which is clearly under-integration. The question is whether there it is also possible to under-integrate in some way with hexahedral elements. I don't know the answer, but I'm not sure how big a deal it is since quadrature points do not require extra memory motion. Once you vectorize finite element residual evaluation, it's common for it to be memory bound, so you might be better off using the flops. $\endgroup$ – Jed Brown Jan 4 '12 at 20:48
  • $\begingroup$ Good point about the memory motion. $\endgroup$ – Geoffrey Irving Jan 4 '12 at 21:01
  • 1
    $\begingroup$ To expand on Jed's point: the reason the above "obvious" argument is false is that each quadrature point sees a $3 \times 3$ matrix. For tetrahedra, that covers all motions of the vertices excluding uniform translation, which doesn't affect energy or forces, so one quadrature point is sufficient for first order accuracy. $\endgroup$ – Geoffrey Irving Jan 4 '12 at 21:03
  • 1
    $\begingroup$ Rather inconvenient that comments can't include newlines. $\endgroup$ – Geoffrey Irving Jan 4 '12 at 21:04
  • $\begingroup$ @JedBrown: Good point. The gradient of linear functions on tets are constants and so a single quadrature point is sufficient, following the argument I made for the mass matrix (the stiffness matrix is the mass matrix for the gradients :-). On the other hand, the gradients of trilinear functions on hexahedra are (anisotropic) quadratic functions so one needs certainly more than just one quadrature point per coordinate direction. $\endgroup$ – Wolfgang Bangerth Jan 5 '12 at 13:31

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.