Skip to content

quad_affine

Implementation of Affine Transformation for Quadrilateral Elements.

This module provides functionality for affine transformations of quadrilateral elements in finite element analysis. It implements mapping between reference and physical elements based on the ParMooN project's methodology.

Key functionalities
  • Reference to physical coordinate mapping
  • Jacobian computation
  • First-order derivatives transformation
  • Second-order derivatives transformation

The implementation follows standard finite element mapping techniques with focus on quadrilateral elements. The transformations maintain geometric consistency and numerical accuracy required for FEM computations.

Key classes
  • QuadAffin: Main class implementing affine transformation for quads
Note

This implementation is specifically referenced from ParMooN project's QuadAffine.C file with adaptations for Python and SciREX framework.

References

[1] ParMooN Project: ParMooN/FiniteElement/QuadAffine.C

Authors

Thivin Anandh D (https://thivinanandh.github.io)

Version

27/Dec/2024: Initial version - Thivin Anandh D

QuadAffin

Bases: FETransforamtion2D

Implements affine transformation for quadrilateral elements.

This class provides methods to transform between reference and physical quadrilateral elements using affine mapping. It handles coordinate transformations, Jacobian computations, and derivative mappings.

Attributes:

Name Type Description
co_ordinates

Array of physical element vertex coordinates Shape: (4, 2) for 2D quadrilateral
x0, (x1, x2, x3)

x-coordinates of vertices
y0, (y1, y2, y3)

y-coordinates of vertices
xc0, (xc1, xc2)

x-coordinate transformation coefficients
yc0, (yc1, yc2)

y-coordinate transformation coefficients
detjk (yc1, yc2)

Determinant of the Jacobian
rec_detjk (yc1, yc2)

Reciprocal of Jacobian determinant
Example

coords = np.array([[0,0], [1,0], [1,1], [0,1]]) quad = QuadAffin(coords) ref_point = np.array([0.5, 0.5]) physical_point = quad.get_original_from_ref(*ref_point)

Note

The implementation assumes counterclockwise vertex ordering and non-degenerate quadrilateral elements.

References

[1] ParMooN Project: QuadAffine.C implementation

Source code in scirex/core/sciml/fe/quad_affine.py 
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
class QuadAffin(FETransforamtion2D):
    """
    Implements affine transformation for quadrilateral elements.

    This class provides methods to transform between reference and physical
    quadrilateral elements using affine mapping. It handles coordinate
    transformations, Jacobian computations, and derivative mappings.

    Attributes:
        co_ordinates: Array of physical element vertex coordinates
            Shape: (4, 2) for 2D quadrilateral
        x0, x1, x2, x3: x-coordinates of vertices
        y0, y1, y2, y3: y-coordinates of vertices
        xc0, xc1, xc2: x-coordinate transformation coefficients
        yc0, yc1, yc2: y-coordinate transformation coefficients
        detjk: Determinant of the Jacobian
        rec_detjk: Reciprocal of Jacobian determinant

    Example:
        >>> coords = np.array([[0,0], [1,0], [1,1], [0,1]])
        >>> quad = QuadAffin(coords)
        >>> ref_point = np.array([0.5, 0.5])
        >>> physical_point = quad.get_original_from_ref(*ref_point)

    Note:
        The implementation assumes counterclockwise vertex ordering and
        non-degenerate quadrilateral elements.

    References:
        [1] ParMooN Project: QuadAffine.C implementation
    """

    def __init__(self, co_ordinates: np.ndarray) -> None:
        """
        Constructor for the QuadAffin class.

        Args:
            co_ordinates: Array of physical element vertex coordinates
                Shape: (4, 2) for 2D quadrilateral

        Returns:
            None
        """
        self.co_ordinates = co_ordinates
        self.set_cell()
        self.get_jacobian(
            0, 0
        )  # 0,0 is just a dummy value # this sets the jacobian and the inverse of the jacobian

    def set_cell(self):
        """
        Set the cell coordinates, which will be used to calculate the Jacobian and actual values.

        Args:
            None

        Returns:
            None
        """

        self.x0 = self.co_ordinates[0][0]
        self.x1 = self.co_ordinates[1][0]
        self.x2 = self.co_ordinates[2][0]
        self.x3 = self.co_ordinates[3][0]

        # get the y-coordinates of the cell
        self.y0 = self.co_ordinates[0][1]
        self.y1 = self.co_ordinates[1][1]
        self.y2 = self.co_ordinates[2][1]
        self.y3 = self.co_ordinates[3][1]

        self.xc0 = (self.x1 + self.x3) * 0.5
        self.xc1 = (self.x1 - self.x0) * 0.5
        self.xc2 = (self.x3 - self.x0) * 0.5

        self.yc0 = (self.y1 + self.y3) * 0.5
        self.yc1 = (self.y1 - self.y0) * 0.5
        self.yc2 = (self.y3 - self.y0) * 0.5

    def get_original_from_ref(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        Returns the original coordinates from the reference coordinates.

        Args:
            xi (np.ndarray): The xi coordinate.
            eta (np.ndarray): The eta coordinate.

        Returns:
            np.ndarray: The transformed original coordinates from the reference coordinates.
        """
        x = self.xc0 + self.xc1 * xi + self.xc2 * eta
        y = self.yc0 + self.yc1 * xi + self.yc2 * eta

        return np.array([x, y])

    def get_jacobian(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        Returns the Jacobian of the transformation.

        Args:
            xi (np.ndarray): The xi coordinate.
            eta (np.ndarray): The eta coordinate.

        Returns:
            np.ndarray: The Jacobian of the transformation.
        """
        self.detjk = self.xc1 * self.yc2 - self.xc2 * self.yc1
        self.rec_detjk = 1 / self.detjk

        return abs(self.detjk)

    def get_orig_from_ref_derivative(self, ref_gradx, ref_grady, xi, eta):
        """
        Returns the derivatives of the original coordinates with respect to the reference coordinates.

        Args:
            ref_gradx (np.ndarray): The reference gradient in the x-direction.
            ref_grady (np.ndarray): The reference gradient in the y-direction.
            xi (np.ndarray): The xi coordinate.
            eta (np.ndarray): The eta coordinate.

        Returns:
            tuple: The derivatives of the original coordinates with respect to the reference coordinates.
        """
        gradx_orig = np.zeros(ref_gradx.shape)
        grady_orig = np.zeros(ref_grady.shape)

        for i in range(ref_gradx.shape[0]):
            gradx_orig[i] = (
                self.yc2 * ref_gradx[i] - self.yc1 * ref_grady[i]
            ) * self.rec_detjk
            grady_orig[i] = (
                -self.xc2 * ref_gradx[i] + self.xc1 * ref_grady[i]
            ) * self.rec_detjk

        return gradx_orig, grady_orig

    def get_orig_from_ref_second_derivative(
        self, grad_xx_ref, grad_xy_ref, grad_yy_ref, xi, eta
    ):
        """
        Returns the second derivatives (xx, xy, yy) of the original coordinates with respect to the reference coordinates.

        Args:
            grad_xx_ref (np.ndarray): The reference second derivative in the x-direction.
            grad_xy_ref (np.ndarray): The reference second derivative in the xy-direction.
            grad_yy_ref (np.ndarray): The reference second derivative in the y-direction.
            xi (np.ndarray): The xi coordinate.
            eta (np.ndarray): The eta coordinate.

        Returns:
            tuple: The second derivatives (xx, xy, yy) of the original coordinates with respect to the reference coordinates.
        """
        GeoData = np.zeros((3, 3))
        Eye = np.identity(3)

        # Populate GeoData (assuming xc1, xc2, yc1, yc2 are defined)
        GeoData[0, 0] = self.xc1 * self.xc1
        GeoData[0, 1] = 2 * self.xc1 * self.yc1
        GeoData[0, 2] = self.yc1 * self.yc1
        GeoData[1, 0] = self.xc1 * self.xc2
        GeoData[1, 1] = self.yc1 * self.xc2 + self.xc1 * self.yc2
        GeoData[1, 2] = self.yc1 * self.yc2
        GeoData[2, 0] = self.xc2 * self.xc2
        GeoData[2, 1] = 2 * self.xc2 * self.yc2
        GeoData[2, 2] = self.yc2 * self.yc2

        # solve the linear system
        solution = np.linalg.solve(GeoData, Eye)

        # generate empty arrays for the original second derivatives
        grad_xx_orig = np.zeros(grad_xx_ref.shape)
        grad_xy_orig = np.zeros(grad_xy_ref.shape)
        grad_yy_orig = np.zeros(grad_yy_ref.shape)

        for j in range(grad_xx_ref.shape[0]):
            r20 = grad_xx_ref[j]
            r11 = grad_xy_ref[j]
            r02 = grad_yy_ref[j]

            grad_xx_orig[j] = (
                solution[0, 0] * r20 + solution[0, 1] * r11 + solution[0, 2] * r02
            )
            grad_xy_orig[j] = (
                solution[1, 0] * r20 + solution[1, 1] * r11 + solution[1, 2] * r02
            )
            grad_yy_orig[j] = (
                solution[2, 0] * r20 + solution[2, 1] * r11 + solution[2, 2] * r02
            )

        return grad_xx_orig, grad_xy_orig, grad_yy_orig

__init__(co_ordinates)

Constructor for the QuadAffin class.

Parameters:

Name Type Description Default
co_ordinates ndarray

Array of physical element vertex coordinates Shape: (4, 2) for 2D quadrilateral

 required

Returns:

Type Description
None

None
Source code in scirex/core/sciml/fe/quad_affine.py 
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
def __init__(self, co_ordinates: np.ndarray) -> None:
    """
    Constructor for the QuadAffin class.

    Args:
        co_ordinates: Array of physical element vertex coordinates
            Shape: (4, 2) for 2D quadrilateral

    Returns:
        None
    """
    self.co_ordinates = co_ordinates
    self.set_cell()
    self.get_jacobian(
        0, 0
    )  # 0,0 is just a dummy value # this sets the jacobian and the inverse of the jacobian

get_jacobian(xi, eta)

Returns the Jacobian of the transformation.

Parameters:

Name Type Description Default
xi ndarray

The xi coordinate.

 required
eta ndarray

The eta coordinate.

 required

Returns:

Type Description
ndarray

np.ndarray: The Jacobian of the transformation.
Source code in scirex/core/sciml/fe/quad_affine.py 
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
def get_jacobian(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    Returns the Jacobian of the transformation.

    Args:
        xi (np.ndarray): The xi coordinate.
        eta (np.ndarray): The eta coordinate.

    Returns:
        np.ndarray: The Jacobian of the transformation.
    """
    self.detjk = self.xc1 * self.yc2 - self.xc2 * self.yc1
    self.rec_detjk = 1 / self.detjk

    return abs(self.detjk)

get_orig_from_ref_derivative(ref_gradx, ref_grady, xi, eta)

Returns the derivatives of the original coordinates with respect to the reference coordinates.

Parameters:

Name Type Description Default
ref_gradx ndarray

The reference gradient in the x-direction.

 required
ref_grady ndarray

The reference gradient in the y-direction.

 required
xi ndarray

The xi coordinate.

 required
eta ndarray

The eta coordinate.

 required

Returns:

Name Type Description
tuple

The derivatives of the original coordinates with respect to the reference coordinates.
Source code in scirex/core/sciml/fe/quad_affine.py 
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
def get_orig_from_ref_derivative(self, ref_gradx, ref_grady, xi, eta):
    """
    Returns the derivatives of the original coordinates with respect to the reference coordinates.

    Args:
        ref_gradx (np.ndarray): The reference gradient in the x-direction.
        ref_grady (np.ndarray): The reference gradient in the y-direction.
        xi (np.ndarray): The xi coordinate.
        eta (np.ndarray): The eta coordinate.

    Returns:
        tuple: The derivatives of the original coordinates with respect to the reference coordinates.
    """
    gradx_orig = np.zeros(ref_gradx.shape)
    grady_orig = np.zeros(ref_grady.shape)

    for i in range(ref_gradx.shape[0]):
        gradx_orig[i] = (
            self.yc2 * ref_gradx[i] - self.yc1 * ref_grady[i]
        ) * self.rec_detjk
        grady_orig[i] = (
            -self.xc2 * ref_gradx[i] + self.xc1 * ref_grady[i]
        ) * self.rec_detjk

    return gradx_orig, grady_orig

get_orig_from_ref_second_derivative(grad_xx_ref, grad_xy_ref, grad_yy_ref, xi, eta)

Returns the second derivatives (xx, xy, yy) of the original coordinates with respect to the reference coordinates.

Parameters:

Name Type Description Default
grad_xx_ref ndarray

The reference second derivative in the x-direction.

 required
grad_xy_ref ndarray

The reference second derivative in the xy-direction.

 required
grad_yy_ref ndarray

The reference second derivative in the y-direction.

 required
xi ndarray

The xi coordinate.

 required
eta ndarray

The eta coordinate.

 required

Returns:

Name Type Description
tuple

The second derivatives (xx, xy, yy) of the original coordinates with respect to the reference coordinates.
Source code in scirex/core/sciml/fe/quad_affine.py 
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
def get_orig_from_ref_second_derivative(
    self, grad_xx_ref, grad_xy_ref, grad_yy_ref, xi, eta
):
    """
    Returns the second derivatives (xx, xy, yy) of the original coordinates with respect to the reference coordinates.

    Args:
        grad_xx_ref (np.ndarray): The reference second derivative in the x-direction.
        grad_xy_ref (np.ndarray): The reference second derivative in the xy-direction.
        grad_yy_ref (np.ndarray): The reference second derivative in the y-direction.
        xi (np.ndarray): The xi coordinate.
        eta (np.ndarray): The eta coordinate.

    Returns:
        tuple: The second derivatives (xx, xy, yy) of the original coordinates with respect to the reference coordinates.
    """
    GeoData = np.zeros((3, 3))
    Eye = np.identity(3)

    # Populate GeoData (assuming xc1, xc2, yc1, yc2 are defined)
    GeoData[0, 0] = self.xc1 * self.xc1
    GeoData[0, 1] = 2 * self.xc1 * self.yc1
    GeoData[0, 2] = self.yc1 * self.yc1
    GeoData[1, 0] = self.xc1 * self.xc2
    GeoData[1, 1] = self.yc1 * self.xc2 + self.xc1 * self.yc2
    GeoData[1, 2] = self.yc1 * self.yc2
    GeoData[2, 0] = self.xc2 * self.xc2
    GeoData[2, 1] = 2 * self.xc2 * self.yc2
    GeoData[2, 2] = self.yc2 * self.yc2

    # solve the linear system
    solution = np.linalg.solve(GeoData, Eye)

    # generate empty arrays for the original second derivatives
    grad_xx_orig = np.zeros(grad_xx_ref.shape)
    grad_xy_orig = np.zeros(grad_xy_ref.shape)
    grad_yy_orig = np.zeros(grad_yy_ref.shape)

    for j in range(grad_xx_ref.shape[0]):
        r20 = grad_xx_ref[j]
        r11 = grad_xy_ref[j]
        r02 = grad_yy_ref[j]

        grad_xx_orig[j] = (
            solution[0, 0] * r20 + solution[0, 1] * r11 + solution[0, 2] * r02
        )
        grad_xy_orig[j] = (
            solution[1, 0] * r20 + solution[1, 1] * r11 + solution[1, 2] * r02
        )
        grad_yy_orig[j] = (
            solution[2, 0] * r20 + solution[2, 1] * r11 + solution[2, 2] * r02
        )

    return grad_xx_orig, grad_xy_orig, grad_yy_orig

get_original_from_ref(xi, eta)

Returns the original coordinates from the reference coordinates.

Parameters:

Name Type Description Default
xi ndarray

The xi coordinate.

 required
eta ndarray

The eta coordinate.

 required

Returns:

Type Description
ndarray

np.ndarray: The transformed original coordinates from the reference coordinates.
Source code in scirex/core/sciml/fe/quad_affine.py 
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
def get_original_from_ref(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    Returns the original coordinates from the reference coordinates.

    Args:
        xi (np.ndarray): The xi coordinate.
        eta (np.ndarray): The eta coordinate.

    Returns:
        np.ndarray: The transformed original coordinates from the reference coordinates.
    """
    x = self.xc0 + self.xc1 * xi + self.xc2 * eta
    y = self.yc0 + self.yc1 * xi + self.yc2 * eta

    return np.array([x, y])

set_cell()

Set the cell coordinates, which will be used to calculate the Jacobian and actual values.

Returns:

Type Description

None
Source code in scirex/core/sciml/fe/quad_affine.py 
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
def set_cell(self):
    """
    Set the cell coordinates, which will be used to calculate the Jacobian and actual values.

    Args:
        None

    Returns:
        None
    """

    self.x0 = self.co_ordinates[0][0]
    self.x1 = self.co_ordinates[1][0]
    self.x2 = self.co_ordinates[2][0]
    self.x3 = self.co_ordinates[3][0]

    # get the y-coordinates of the cell
    self.y0 = self.co_ordinates[0][1]
    self.y1 = self.co_ordinates[1][1]
    self.y2 = self.co_ordinates[2][1]
    self.y3 = self.co_ordinates[3][1]

    self.xc0 = (self.x1 + self.x3) * 0.5
    self.xc1 = (self.x1 - self.x0) * 0.5
    self.xc2 = (self.x3 - self.x0) * 0.5

    self.yc0 = (self.y1 + self.y3) * 0.5
    self.yc1 = (self.y1 - self.y0) * 0.5
    self.yc2 = (self.y3 - self.y0) * 0.5