Seismo-VLAB  1.3
An Open-Source Finite Element Software for Meso-Scale Simulations

Background

The PML3DHexa8 class creates a perfectly matched layer (PML) element with eight nodes, which can be used to define absorbing boundary layer to encapsulate the regular domain. This element is based on the symmetric hybrid formulation, which results in mixed finite element formulation (displacement and stress fields) in the PML elements. Each node of the PML element has nine degrees of freedom, three of which for translation in X, Y, ans Z directions and six of which for stress tensor histories.

The key idea in PML elements is the use of complex coordinate stretching. The complex stretching function in direction s-, which can be either x-, y-, or z- directions in the cartesian coordinate system, is defined as:

\[ \epsilon_s(s,\omega) = \alpha_s(s) + \frac{\beta_s(s)}{i\omega} \]

where \(\alpha_s\) and \(\beta_s\) are scaling and attenuation functions, respectively. The form of these functions are:

\[ \alpha_s (s) = 1+\alpha_0\left [ \frac{(s-s_0)n_s}{L_\textrm{pml}}\right ]^m, \quad \beta_s (s) = \beta_0 \left [ \frac{(s-s_0)n_s}{L_\textrm{pml}}\right ]^m \]

where \(m\) is the degree of the polynomial; \(n_s\) is the \(s\)th component of outward normal to the interface between PML region and the regular domain; \(L_\textrm{pml}\) is the thickness of the PML region; and \(s_0\) is the \(s\)th component of the reference point where stretching is defined. \(\alpha_0\) and \(\beta_0\) are defined as:

\[ \alpha_0 = \frac{3b}{2L_\textrm{pml}}\log \frac{1}{R}, \quad \beta_0 = \frac{3V_\textrm{ref}}{2L_\textrm{pml}}\log \frac{1}{R} \]

where \(b\) is a characteristic length and we set it to \(L_\textrm{pml}/10\); \(V_\textrm{ref}\) is a reference velocity, and \(R\) is a user-defined reflection coefficient. The Figure shows PML elements at different positions. The figure below show how to define the reference point (white dot) – to start coordinate stretching – and the normal vector (white arrow) for some element region \(E_i, i = 1,\dots, 17\).

PML3DNormals.png

REFERENCE:

  • Fathi, A., Poursartip, B., & Kallivokas, L. F. (2015). "Time‐domain hybrid formulations for wave simulations in three‐dimensional PML‐truncated heterogeneous media". International Journal for Numerical Methods in Engineering, 101(3), 165-198.

Pre-Analysis

The python Pre-Analysis in the 01-Pre_Process/Method/Attach.py file provides with an interface to populate the Entities dictionary. This file contains several functions to populate specific fields. For example, to create a PML3DHexa8, using json format, use:

  • addElement(tag, name='PML3DHexa8', conn, attributes):

    • tag : The identifier of this element, i.e., tag > -1
    • name : Seismo-VLAB element class name
    • conn : Connectivity array of this element
    • attributes : Specific properties for the created element, for example
      • 'material' is material identifier.
      • 'n' is the polynomial degree.
      • 'L' is the Thickness of the PML region.
      • 'R' is the Reflection coefficient.
      • 'x0' is the PML coordinate of the reference point.
      • 'npml' is the vector of the outward normal.
      • 'np' Number of quadrature-points used for the integration.
      • 'rule' the integration rule: 'Gauss', 'Lobatto'.

    Example

    A PML3DHEXA8 element can be defined using the python interface as follows:
    SVL.addElement(tag=1, name='PML3DHexa8', conn=[1,2,3,4,5,6,7,8], attributes={'material': 1, 'np': 8, 'rule': Gauss, 'x0': [0,0,0], 'npml': [1,0,0], 'L': 10.0})

    Application
    Please refer to the J12-DY_Axial_Load_Long_Rod_PML3D.py file located at 03-Validations/02-Performance/ to see an example on how to define a PML3DHexa8 element.

On the contrary, the 01-Pre_Process/Method/Remove.py file provides with an interface to depopulate the Entities dictionary. For example, to remove an already define Element, use:

  • delElement(tag):
    • tag : The identifier of the element to be removed, i.e., tag > -1

Run-Analysis

The C++ Run-Analysis in the 02-Run_Process/04-Elements/06-Quadrilateral/PML3DHexa8.cpp file provides the class implementation. A PML3DHexa8 element is created using the built-in json parse-structure provided in the Driver.hpp. A PML3DHexa8 is defined inside the "Elements" json field indicating its "Tag" as follows,

  • {
        "Elements": {
            "Tag": {
                "name" : "PML3DHEXA8",
                "conn" : [ ],
                "attributes": {
                    "n": double,
                    "L": double,
                    "R": double,
                    "x0": [ ],
                    "npml": [ ],
                    "np": int,
                    "rule": str,
                    "material": int
                }
            }
        }
    }
    
    Variable Description
    Tag Unique element object identifier.
    conn The element connectivity node array.
    n Polynomial degree.
    L Thickness of the PML region.
    R Reflection coefficient.
    x0 Coordinates \((x_0,y_0)\) of the reference point.
    npml Vector components \((n_x,n_y)\) of the outward normal.
    np Number of quadrature-points used for the integration.
    rule The quadrature integration name rule=GAUSS, LOBATTO.
    materialThe material identifier.

    Attention
    This formulation only allows small deformations.
    Material non-linearity is not tested in this element.
    Currently the number of integration points np= 1,8,27,64,125,216,343 are supported.
    The material is considered to be the same at each quadrature point.
    The domain reduction force method is not compatible with this element.
    For the spatial extent of the PML layer–i.e., \(\textrm{L}_\textrm{pml}\), it is problem specific, but as a rule of thumb, considering one wavelength associated with the dominant frequency should work in general.
    For the PML element size, one can consider about 12 points in one wavelength. Once the mesh size is fixed, satisfying CFL condition to determine \(\Delta_\textrm{t}\) (i.e, \(\Delta_\textrm{t} \leq 0.3\,\Delta_\textrm{h}/V_p\), with \(V_p\) the p-wave velocity) for time marching.
    Example

    Consider the white domain of PML3DHEXA8 elements showed in Figure above, the box has length D in x and y and H in z.
    {
    "Elements": {
    "1": { "name" : "PML3DHEXA8", "conn" : [1,2,3,4,5,6,7,8], "attributes": { "n": 2.0, "L": 5.0, "R": 0.001, "x0": [D,0,0], "npml": [-0.707, 0.0, -0.707], "np": 8, "rule": "GAUSS", "material": 1 } },
    "2": { "name" : "PML2DQUAD8", "conn" : [1,2,3,4,5,6,7,8], "attributes": { "n": 2.0, "L": 5.0, "R": 0.001, "x0": [0,0,0], "npml": [0.0, 0.0, -1.0], "np": 8, "rule": "GAUSS", "material": 1 } },
    "3": { "name" : "PML2DQUAD8", "conn" : [1,2,3,4,5,6,7,8], "attributes": { "n": 2.0, "L": 5.0, "R": 0.001, "x0": [0,D,H/2], "npml": [0.0, 1.0, 0.0], "np": 8, "rule": "GAUSS", "material": 1 } } }
    }