TET4 element with one quadrature point

[rcontinuum]

Hello,

my aim is to use TET4 elements with one quadrature point. But with the following input fragment moose selects me a quadrature rule with four five quadrature points.
What should I change so that a quadrature rule with only one quadrature point gets selected?

I read a mesh file with one TET4 element, generated by gmsh, as

[GlobalParams]
  displacements = 'disp_x disp_y disp_z'
  large_kinematics = false
[]

[Mesh]
  [msh]
    type = FileMeshGenerator
    file = '1tet4.msh'
    show_info = true
  []
[]

[Variables]
  [disp_x]
    order = FIRST
    family = LAGRANGE
  []
  [disp_y]
    order = FIRST
    family = LAGRANGE
  []
  [disp_z]
    order = FIRST
    family = LAGRANGE
  []
[]

Mesh check

FileMeshGenerator 'msh':  Mesh Subdomains:
FileMeshGenerator 'msh':   Subdomain 1 (Tetrahedron): 1 elems (TET4, 1 active), 4 active nodes
FileMeshGenerator 'msh':    Volume: 0.166667
FileMeshGenerator 'msh':    Bounding box minimum: (x,y,z)=(       0,        0,        0)
FileMeshGenerator 'msh':    Bounding box maximum: (x,y,z)=(       1,        1,        1)
FileMeshGenerator 'msh':    Bounding box delta: (x,y,z)=(       1,        1,        1)
FileMeshGenerator 'msh':   Global mesh volume = 0.166667

But in my src/materials "subroutine" for calculating stresses and tangent operator the variable _qp runs from 0 to 4 3 .
Thus, it's a 4 5-point TET4 quadrature scheme. That works, but I was hoping to be able to select the one point quadrature rule.
Like in libmesh the case FIRST: in line 53 ff. of the file

https://mooseframework.inl.gov/docs/doxygen/libmesh/quadrature__gauss__3D_8C_source.html#l00029

//---------------------------------------------
      // Tetrahedral quadrature rules
    case TET4:
    case TET10:
    case TET14:
      {
        switch(get_order())
          {
            // Taken from pg. 222 of "The finite element method," vol. 1
            // ed. 5 by Zienkiewicz & Taylor
          case CONSTANT:
          case FIRST:
            {
              // Exact for linears
              _points.resize(1);
              _weights.resize(1);
 
 
              _points[0](0) = .25;
              _points[0](1) = .25;
              _points[0](2) = .25;
 
              _weights[0] = Real(1)/6;
 
              return;
            }

If I use [Variables] ... order = CONSTANT I get an error message "Unsupported order: 0"

Any hints?

Edit: four -> five

[GiudGiud]

Hello

You may use a constant monomial for that if you want to have this smaller quadrature all over the domain

[Variables]
  [u]
    family = MONOMIAL
    order = CONSTANT
  []
[]

note that the kernels and boundary conditions are often different for these discontinuous FEM families.

Guillaume

[rcontinuum]

Well, if I do this, then I get the runtime error (as you suspected)

*** ERROR ***
Trying to use nodal boundary condition 'fix_x' on a non-nodal variable 'disp_x'.

Because I want to apply boundary conditions at nodes (typical for Tensor Mechanics) like

[BCs]
  [fix_x]
    type = DirichletBC
    preset = true
    boundary = 'Point_1 Point_3 Point_4'
    variable = disp_x
    value = 0.0
  []
  ...

I mean, the numerical integration of the stiffness matrix of a TET4 element with one quadrature point rule is exact. As the gradients of the shape functions are constants. Sure, for a mass matrix that is not true. But even for the mass matrix it works from a practical point of view.

Do you know any other possibility to select the quadrature rule "manually" for C0 elements like the TET4 element?

Edit: Btw, _qp runs from 0 to 4. So it's actually a 5-point quadrature rule.

[GiudGiud]

are you trying to have a reduced quadrature only on the boundary?

[rcontinuum]

No, I don't have boundary integrals.
Sure, if I would have, e.g., NeumannBCs then one quadrature point for the triangular surface of the element would not be enough.

I only apply DirichletBCs at the nodes (for debugging reasons). Thus I don't have integrals over the boundary of the domain.
For the usual Tensor Mechanics kernel, I use TotalLagrangianStressDivergence but doesn't really matter, we have only one volume integral to compute the stiffness matrix. And that volume integral is computed exact with a one-point quadrature rule for the TET4 element.

[hugary1995]

Do this:

[Executioner]
  [Quadrature]
    order = CONSTANT
  []
[]

[hugary1995]

Ref: https://mooseframework.inl.gov/moose/source/actions/SetupQuadratureAction.html

[rcontinuum]

Thank you!
That works!