OOF2: The Manual

Name

PhysicalProperty — Intermediate base class for FluxProperty and EqnProperty

Synopses

Extension Properties derived from FluxProperty or EqnProperty must redefine integration_order(), in either C++ or Python. This is the only method in the Property class hierarchy that must be redefined in subclasses.

C++ Synopsis

#include "engine/property.h"
class PhysicalProperty: , public Property {
  virtual int integration_order(const CSubProblem* subproblem,
                                const Element* element) const;

}

C++ Derived Classes

C++ Source Files

  • SRC/engine/property.h: C++ header
  • SRC/engine/property.swg: SWIG source code

Python Synopsis

The Python classes PyFluxProperty and PyEqnProperty are derived from a C++ base class, PyPhysicalPropertyMethods. PyPhysicalPropertyMethods defines only one function, integration_order(), that can be redefined in a Python subclass.[76]

from ooflib.SWIG.engine import pypropertywrapper 
class PyPhysicalPropertyMethods:
  def integration_order(self, subproblem, element)

Python Derived Classes

Source Files

  • SRC/engine/pypropertywrapper.C: C++ code
  • SRC/engine/pypropertywrapper.h: C++ header

Overview

PhysicalProperty is a subclass of Property and a base class of FluxProperty and EqnProperty, the classes that define the physics behavior of materials. Users should not derive new classes from PhysicalProperty directly.

Methods

int integration_order(const CSubProblem *subproblem, const Element *element) const

When OOF2 converts an equation like (2.8) or (2.9) to matrix form, it replaces the Fields by an expansion in the Element's shape functions Ni⁒(𝐱) , multiplies everything by a test function (which is another shape function, Nj⁒(𝐱) ), and then integrates over the entire mesh (i.e, over each element). Divergence terms are integrated by parts, which replaces the test function by its derivative.

The integration over the elements is done by Gaussian quadrature. To choose the appropriate number of Gauss points, it's necessary to know the polynomial degree of the integrand in x and y. Generally, higher order polynomials require more points. Choosing an integration order that is too large will slow down the calculation, but choosing one that is too small could lead to incorrect results.

The method integration_order(), which must be defined in each PhysicalProperty class, returns the polynomial degree of the Property's contribution to the integrand. That is, it depends on the Property before multiplying by the test function and integrating. The degree depends on the details of the Property and on the polynomial degrees of the shape functions and their derivatives, which can be obtained from Element::shapefun_degree() and Element::dshapefun_degree().

For example, a Property that adds a constant to the force term in (2.8) or (2.9) has no dependence on the Fields and therefore has no shape functions in its expansion, so it has order 0. It should define integration_order() like this:

int ForceDensity::integration_order(const CSubProblem*, const Element*) const
{
  return 0;
} 

Sometimes the polynomial order depends on whether the Field is in-plane or out-of-plane. Remember that OOF2 adds an auxiliary field,

ψ⁒(𝐱)β‰‘βˆ‚Ο†β’(𝐱)βˆ‚z
(183)

that contains the z -derivative of Fields that are not in-plane. The gradient of such a Field is

βˆ‚Ο†β’(𝐱)βˆ‚π±=βˆ‚Ο†β’(𝐱)βˆ‚x⁒𝐱^+βˆ‚Ο†β’(𝐱)βˆ‚y⁒𝐲^+ψ⁒(𝐱)⁒𝐳^
(184)

Both of the fields Ο† and ψ are expanded in shape functions

φ⁒(𝐱) =βˆ‘Ξ½Nν⁒(𝐱)⁒φν
ψ⁒(𝐱) =βˆ‘Ξ½Nν⁒(𝐱)⁒ψν
(185)

so the gradient (184) is

βˆ‚Ο†β’(x)βˆ‚π±=βˆ‘Ξ½[βˆ‚Nν⁒(𝐱)βˆ‚x⁒φν⁒𝐱^+βˆ‚Nν⁒(𝐱)βˆ‚y⁒φν⁒𝐲^+Nν⁒(𝐱)⁒ψν⁒𝐱^].
(186)

Therefore if the Field is in-plane and ψ is zero, the polynomial degree of the gradient is given by the derivative of the shapefunction, but if the Field is out-of-plane, the degree is the degree of the shapefunction itself, because the polynomial degree of a function is higher that the degree of its derivative.

The result of all this is that a Property that depends on the gradient of a Field will need to check whether or not the Field is in-plane. For example, here is the integration_order() method from Elasticity:

int Elasticity::integration_order(const CSubProblem *subp, const Element *el)
  const
{
  if(displacement->in_plane(subp)) 1
    return el->dshapefun_degree();
  return el->shapefun_degree();
 } 

1

displacement is a TwoVectorField* stored in the Elasticity class, and initialized in the Elasticity constructor like this:

displacement = dynamic_cast<TwoVectorField*>(Field::getField("Displacement")); 

[Warning] Warning

Do not compute the shapefunction's derivative's degree by subtracting 1 from the shapefunction's degree. For purposes of Gaussian integration, the degree of the shapefunction is sometimes less than its actual polynomial degree. For example, the linear quadrilateral shapefunction (1βˆ’x)⁒(1βˆ’y)/4 for a 2×2 square centered on the origin can be integrated exactly with a single Gauss point at (0,0), although its polynomial degree is 2. For this function, both shapefun_degree() and dshapefun_degree() return 1.



[76] As with PyPropertyMethods, there is no actual Python version of this class. But a class derived from PyFluxProperty or PyEqnProperty will act as if it's derived from PyPhysicalPropertyMethods in some ways.