Background

The PML2DQuad8 class creates a plane strain 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 five degrees of freedom, two of which for translation in X and Y directions and three 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- or y- 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 examples below show how to define the reference point – to start coordinate stretching – and the normal vector for each element \(E_i, i = 1,\dots, 5\).

REFERENCE:

Kucukcoban, S., and Loukas F. Kallivokas. "A symmetric hybrid formulation for transient wave simulations in PML-truncated heterogeneous media." Wave Motion Vol 50. Issue 1 (2013): 57-79.

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 PML2DQuad8, using json format, use:

addElement(tag, name='PML2DQuad8', 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.
  - 'th' is the thickness of a quad element.
  - '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 PML2DQUAD8 element can be defined using the python interface as follows:
SVL.addElement(tag=1, name='PML2DQuad8', conn=[1,2,3,4,5,6,7,8], attributes={'th': 0.01, 'material': 1, 'np': 9, 'rule': Gauss, 'x0': [], 'npml': [], 'L': 10.0})

Application
Please refer to the C05-DY_Lin_2D_ElasticPStrainPML_MUMPS_Quad8.py file located at 03-Validations/02-Performance/ to see an example on how to define a PML2DQuad8 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/PML2DQuad8.cpp file provides the class implementation. A PML2DQuad8 element is created using the built-in json parse-structure provided in the Driver.hpp. A PML2DQuad8 is defined inside the "Elements" json field indicating its "Tag" as follows,

{
    "Elements": {
        "Tag": {
            "name" : "PML2DQUAD8",
            "conn" : [ ],
            "attributes": {
                "th": double,
                "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.
`th`	The thickness of the element.
`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.
`material`	The biaxial 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,4,9,16,25,36,49 are supported, another value assumes np=9.; 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 PML2DQUAD8 elements showed in Figure with D = 5, L = 20, and H = 10.
{
"Elements": {
"1": { "name" : "PML2DQUAD8", "conn" : [1,2,3,4,5,6,7,8], "attributes": { "th": 1.0, "n": 2.0, "L": 5.0, "R": 0.001, "x0": [5,10], "npml": [-1.0, 0.0], "np": 4, "rule": "GAUSS", "material": 1 } },
"2": { "name" : "PML2DQUAD8", "conn" : [1,2,3,4,5,6,7,8], "attributes": { "th": 1.0, "n": 2.0, "L": 5.0, "R": 0.001, "x0": [5,5], "npml": [-0.707, -0.707], "np": 4, "rule": "GAUSS", "material": 1 } },
"3": { "name" : "PML2DQUAD8", "conn" : [1,2,3,4,5,6,7,8], "attributes": { "th": 1.0, "n": 2.0, "L": 5.0, "R": 0.001, "x0": [15,5], "npml": [0.0, -1.0], "np": 4, "rule": "GAUSS", "material": 1 } },
"4": { "name" : "PML2DQUAD8", "conn" : [1,2,3,4,5,6,7,8], "attributes": { "th": 1.0, "n": 2.0, "L": 5.0, "R": 0.001, "x0": [25,5], "npml": [0.707, -0.707], "np": 4, "rule": "GAUSS", "material": 1 } },
"5": { "name" : "PML2DQUAD8", "conn" : [1,2,3,4,5,6,7,8], "attributes": { "th": 1.0, "n": 2.0, "L": 5.0, "R": 0.001, "x0": [25,10], "npml": [1.0, 0.0], "np": 4, "rule": "GAUSS", "material": 1 } } }
}