FE2D_Cell

Module: FE2D_Cell.py

This module provides functionality for setting up and managing finite element calculations for individual 2D cells, including basis functions, quadrature rules, and transformations.

Classes:

Name	Description
FE2D_Cell	Main class for managing cell-level FE computations

Dependencies

• basis_function_2d: Base classes for 2D basis functions
• quadratureformulas_quad2d: Quadrature rules for 2D elements
• fe2d_setup_main: Setup utilities for 2D FE calculations
• numpy: Numerical computations

Key Features

• Cell-level finite element value storage
• Basis function evaluation at quadrature points
• Reference to physical domain transformations
• Gradient and derivative computations
• Quadrature rule implementation
• Forcing function integration
• Support for different element types and orders

Authors

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

Version Info

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

Notes

The implementation includes optimization for gradient calculations where grad_x_orig and grad_y_orig store multiplication factors for reference gradients to improve computational efficiency.

Basis2DQNChebyshev2

Bases: BasisFunction2D

A specialized implementation of two-dimensional basis functions using Chebyshev polynomials for Q1 elements.

This class provides a complete implementation for computing basis functions and their derivatives in two dimensions, specifically designed for use in variational physics-informed neural networks (VPINNs) with domain decomposition. The basis functions are constructed using Chebyshev polynomials through Jacobi polynomial representations.

The class inherits from BasisFunction2D and implements all required methods for computing function values and derivatives. The implementation follows the methodology described in hp-VPINNs research by Ehsan Kharazmi et al.

Attributes:

Name	Type	Description
num_shape_functions	int	Total number of shape functions in the 2D element. Must be a perfect square as it represents tensor product of 1D functions.

Methods:

Name	Description
value	Computes values of all basis functions at given points
gradx	Computes x-derivatives of all basis functions
grady	Computes y-derivatives of all basis functions
gradxx	Computes second x-derivatives of all basis functions
gradyy	Computes second y-derivatives of all basis functions
gradxy	Computes mixed xy-derivatives of all basis functions

Implementation Details

• Basis functions are constructed as tensor products of 1D test functions
• Test functions are derived from normalized Jacobi polynomials
• Special cases are handled for first few polynomial degrees in derivatives
• All computations maintain double precision (float64)
• Efficient vectorized operations using numpy arrays

Example

basis = Basis2DQNChebyshev2(num_shape_functions=16)  # Creates 4x4 basis functions
xi = np.linspace(-1, 1, 100)
eta = np.linspace(-1, 1, 100)
values = basis.value(xi, eta)
x_derivatives = basis.gradx(xi, eta)

Notes

• num_shape_functions must be a perfect square
• All coordinate inputs (xi, eta) should be in the range [-1, 1]
• Implementation optimized for vectorized operations on multiple points
• Based on hp-VPINNs methodology: https://github.com/ehsankharazmi/hp-VPINNs/

References

Kharazmi, E., et al. "hp-VPINNs: Variational Physics-Informed Neural Networks With Domain Decomposition"

Source code in scirex/core/sciml/fe/basis_2d_qn_chebyshev_2.py

class Basis2DQNChebyshev2(BasisFunction2D):
    """A specialized implementation of two-dimensional basis functions using Chebyshev polynomials for Q1 elements.

    This class provides a complete implementation for computing basis functions and their derivatives
    in two dimensions, specifically designed for use in variational physics-informed neural networks
    (VPINNs) with domain decomposition. The basis functions are constructed using Chebyshev polynomials
    through Jacobi polynomial representations.

    The class inherits from BasisFunction2D and implements all required methods for computing
    function values and derivatives. The implementation follows the methodology described in
    hp-VPINNs research by Ehsan Kharazmi et al.

    Attributes:
        num_shape_functions (int): Total number of shape functions in the 2D element.
            Must be a perfect square as it represents tensor product of 1D functions.

    Methods:
        value(xi, eta): Computes values of all basis functions at given points
        gradx(xi, eta): Computes x-derivatives of all basis functions
        grady(xi, eta): Computes y-derivatives of all basis functions
        gradxx(xi, eta): Computes second x-derivatives of all basis functions
        gradyy(xi, eta): Computes second y-derivatives of all basis functions
        gradxy(xi, eta): Computes mixed xy-derivatives of all basis functions

    Implementation Details:
        - Basis functions are constructed as tensor products of 1D test functions
        - Test functions are derived from normalized Jacobi polynomials
        - Special cases are handled for first few polynomial degrees in derivatives
        - All computations maintain double precision (float64)
        - Efficient vectorized operations using numpy arrays

    Example:
        ```python
        basis = Basis2DQNChebyshev2(num_shape_functions=16)  # Creates 4x4 basis functions
        xi = np.linspace(-1, 1, 100)
        eta = np.linspace(-1, 1, 100)
        values = basis.value(xi, eta)
        x_derivatives = basis.gradx(xi, eta)
        ```

    Notes:
        - num_shape_functions must be a perfect square
        - All coordinate inputs (xi, eta) should be in the range [-1, 1]
        - Implementation optimized for vectorized operations on multiple points
        - Based on hp-VPINNs methodology: https://github.com/ehsankharazmi/hp-VPINNs/

    References:
        Kharazmi, E., et al. "hp-VPINNs: Variational Physics-Informed Neural Networks
        With Domain Decomposition"
    """

    def __init__(self, num_shape_functions: int):
        super().__init__(num_shape_functions)

    def jacobi_wrapper(self, n: int, a: int, b: int, x: np.ndarray) -> np.ndarray:
        """Evaluates Jacobi polynomial at specified points.

        Computes values of nth degree Jacobi polynomial with parameters (a,b)
        at given points x.

        Args:
            n: Degree of Jacobi polynomial. Must be non-negative integer.
            a: First parameter of Jacobi polynomial
            b: Second parameter of Jacobi polynomial
            x: Points at which to evaluate polynomial
                Shape: (n_points,)

        Returns:
            np.ndarray: Values of Jacobi polynomial at input points
                Shape: Same as input x

        Notes:
            Wrapper around scipy.special.jacobi that ensures float64 precision
            and proper array handling.
        """
        x = np.array(x, dtype=np.float64)
        return jacobi(n, a, b)(x)

    ## Helper Function
    def test_fcnx(self, n_test: int, x: np.ndarray) -> np.ndarray:
        """Computes x-component test functions.

        Evaluates the x-direction test functions constructed as differences
        of normalized Jacobi polynomials.

        Args:
            n_test: Number of test functions to compute
            x: Points at which to evaluate functions
                Shape: (n_points,)

        Returns:
            np.ndarray: Values of test functions at input points
                Shape: (n_test, n_points)

        Notes:
            Test functions are constructed as differences of normalized Jacobi
            polynomials following hp-VPINNs methodology.
        """
        test_total = []
        for n in range(1, n_test + 1):
            test = self.jacobi_wrapper(n + 1, -1 / 2, -1 / 2, x) / self.jacobi_wrapper(
                n + 1, -1 / 2, -1 / 2, 1
            ) - self.jacobi_wrapper(n - 1, -1 / 2, -1 / 2, x) / self.jacobi_wrapper(
                n - 1, -1 / 2, -1 / 2, 1
            )
            test_total.append(test)
        return np.asarray(test_total, np.float64)

    def test_fcny(self, n_test: int, y: np.ndarray) -> np.ndarray:
        """Computes y-component test functions.

        Evaluates the y-direction test functions constructed as differences
        of normalized Jacobi polynomials.

        Args:
            n_test: Number of test functions to compute
            y: Points at which to evaluate functions
                Shape: (n_points,)

        Returns:
            np.ndarray: Values of test functions at input points
                Shape: (n_test, n_points)

        Notes:
            Test functions are constructed as differences of normalized Jacobi
            polynomials following hp-VPINNs methodology.
        """
        test_total = []
        for n in range(1, n_test + 1):
            test = self.jacobi_wrapper(n + 1, -1 / 2, -1 / 2, y) / self.jacobi_wrapper(
                n + 1, -1 / 2, -1 / 2, 1
            ) - self.jacobi_wrapper(n - 1, -1 / 2, -1 / 2, y) / self.jacobi_wrapper(
                n - 1, -1 / 2, -1 / 2, 1
            )
            test_total.append(test)
        return np.asarray(test_total, np.float64)

    def dtest_fcn(self, n_test: int, x: np.ndarray) -> np.ndarray:
        """Computes first and second derivatives of test functions.

        Calculates derivatives of test functions constructed from Jacobi
        polynomials, handling special cases for n=1,2 separately.

        Args:
            n_test: Number of test functions
            x: Points at which to evaluate derivatives
                Shape: (n_points,)

        Returns:
            tuple(np.ndarray, np.ndarray): First and second derivatives
                First element: First derivatives, shape (n_test, n_points)
                Second element: Second derivatives, shape (n_test, n_points)

        Notes:
            Special cases for n=1,2 ensure proper derivative calculations
            following hp-VPINNs methodology.
        """
        d1test_total = []
        d2test_total = []
        for n in range(1, n_test + 1):
            if n == 1:
                d1test = (
                    ((n + 1) / 2)
                    * self.jacobi_wrapper(n, 1 / 2, 1 / 2, x)
                    / self.jacobi_wrapper(n + 1, -1 / 2, -1 / 2, 1)
                )
                d2test = (
                    ((n + 2) * (n + 1) / (2 * 2))
                    * self.jacobi_wrapper(n - 1, 3 / 2, 3 / 2, x)
                    / self.jacobi_wrapper(n + 1, -1 / 2, -1 / 2, 1)
                )
                d1test_total.append(d1test)
                d2test_total.append(d2test)
            elif n == 2:
                d1test = ((n + 1) / 2) * self.jacobi_wrapper(
                    n, 1 / 2, 1 / 2, x
                ) / self.jacobi_wrapper(n + 1, -1 / 2, -1 / 2, 1) - (
                    (n - 1) / 2
                ) * self.jacobi_wrapper(
                    n - 2, 1 / 2, 1 / 2, x
                ) / self.jacobi_wrapper(
                    n - 1, -1 / 2, -1 / 2, 1
                )
                d2test = (
                    ((n + 2) * (n + 1) / (2 * 2))
                    * self.jacobi_wrapper(n - 1, 3 / 2, 3 / 2, x)
                    / self.jacobi_wrapper(n + 1, -1 / 2, -1 / 2, 1)
                )
                d1test_total.append(d1test)
                d2test_total.append(d2test)
            else:
                d1test = ((n + 1) / 2) * self.jacobi_wrapper(
                    n, 1 / 2, 1 / 2, x
                ) / self.jacobi_wrapper(n + 1, -1 / 2, -1 / 2, 1) - (
                    (n - 1) / 2
                ) * self.jacobi_wrapper(
                    n - 2, 1 / 2, 1 / 2, x
                ) / self.jacobi_wrapper(
                    n - 1, -1 / 2, -1 / 2, 1
                )
                d2test = ((n + 2) * (n + 1) / (2 * 2)) * self.jacobi_wrapper(
                    n - 1, 3 / 2, 3 / 2, x
                ) / self.jacobi_wrapper(n + 1, -1 / 2, -1 / 2, 1) - (
                    (n) * (n - 1) / (2 * 2)
                ) * self.jacobi_wrapper(
                    n - 3, 3 / 2, 3 / 2, x
                ) / self.jacobi_wrapper(
                    n - 1, -1 / 2, -1 / 2, 1
                )
                d1test_total.append(d1test)
                d2test_total.append(d2test)
        return np.asarray(d1test_total), np.asarray(d2test_total)

    def value(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """Evaluates basis functions at given coordinates.

        Computes values of all basis functions at specified (xi,eta) points
        using tensor product of 1D test functions.

        Args:
            xi: x-coordinates at which to evaluate functions
                Shape: (n_points,)
            eta: y-coordinates at which to evaluate functions
                Shape: (n_points,)

        Returns:
            np.ndarray: Values of all basis functions
                Shape: (num_shape_functions, n_points)

        Notes:
            Basis functions are constructed as products of 1D test functions
            in x and y directions.
        """
        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
        test_x = self.test_fcnx(num_shape_func_in_1d, xi)
        test_y = self.test_fcny(num_shape_func_in_1d, eta)
        values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

        for i in range(num_shape_func_in_1d):
            values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
                test_x[i, :] * test_y
            )

        return values

    def gradx(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """Computes x-derivatives of basis functions.

        Evaluates partial derivatives with respect to x of all basis
        functions at given coordinates.

        Args:
            xi: x-coordinates at which to evaluate derivatives
                Shape: (n_points,)
            eta: y-coordinates at which to evaluate derivatives
                Shape: (n_points,)

        Returns:
            np.ndarray: Values of x-derivatives
                Shape: (num_shape_functions, n_points)

        Notes:
            Uses product rule with x-derivatives of test functions in
            x-direction and values in y-direction.
        """
        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
        grad_test_x = self.dtest_fcn(num_shape_func_in_1d, xi)[0]
        test_y = self.test_fcny(num_shape_func_in_1d, eta)
        values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

        for i in range(num_shape_func_in_1d):
            values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
                grad_test_x[i, :] * test_y
            )

        return values

    def grady(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """Computes y-derivatives of basis functions.

        Evaluates partial derivatives with respect to y of all basis
        functions at given coordinates.

        Args:
            xi: x-coordinates at which to evaluate derivatives
                Shape: (n_points,)
            eta: y-coordinates at which to evaluate derivatives
                Shape: (n_points,)

        Returns:
            np.ndarray: Values of y-derivatives
                Shape: (num_shape_functions, n_points)

        Notes:
            Uses product rule with values in x-direction and y-derivatives
            of test functions in y-direction.
        """
        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
        test_x = self.test_fcnx(num_shape_func_in_1d, xi)
        grad_test_y = self.dtest_fcn(num_shape_func_in_1d, eta)[0]
        values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

        for i in range(num_shape_func_in_1d):
            values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
                test_x[i, :] * grad_test_y
            )

        return values

    def gradxx(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """Computes second x-derivatives of basis functions.

        Evaluates second partial derivatives with respect to x of all basis
        functions at given coordinates.

        Args:
            xi: x-coordinates at which to evaluate derivatives
                Shape: (n_points,)
            eta: y-coordinates at which to evaluate derivatives
                Shape: (n_points,)

        Returns:
            np.ndarray: Values of second x-derivatives
                Shape: (num_shape_functions, n_points)

        Notes:
            Uses product rule with second x-derivatives of test functions in
            x-direction and values in y-direction.
        """
        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
        grad_grad_x = self.dtest_fcn(num_shape_func_in_1d, xi)[1]
        test_y = self.test_fcny(num_shape_func_in_1d, eta)
        values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

        for i in range(num_shape_func_in_1d):
            values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
                grad_grad_x[i, :] * test_y
            )

        return values

    def gradxy(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """Computes second x-derivatives of basis functions.

        Evaluates second partial derivatives with respect to x of all basis
        functions at given coordinates.

        Args:
            xi: x-coordinates at which to evaluate derivatives
                Shape: (n_points,)
            eta: y-coordinates at which to evaluate derivatives
                Shape: (n_points,)

        Returns:
            np.ndarray: Values of second x-derivatives
                Shape: (num_shape_functions, n_points)

        Notes:
            Uses product rule with second x-derivatives of test functions in
            x-direction and y derivative values in y-direction.
        """
        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
        grad_test_x = self.dtest_fcn(num_shape_func_in_1d, xi)[0]
        grad_test_y = self.dtest_fcn(num_shape_func_in_1d, eta)[0]
        values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

        for i in range(num_shape_func_in_1d):
            values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
                grad_test_x[i, :] * grad_test_y
            )

        return values

    def gradyy(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """Computes second x-derivatives of basis functions.

        Evaluates second partial derivatives with respect to x of all basis
        functions at given coordinates.

        Args:
            xi: x-coordinates at which to evaluate derivatives
                Shape: (n_points,)
            eta: y-coordinates at which to evaluate derivatives
                Shape: (n_points,)

        Returns:
            np.ndarray: Values of second x-derivatives
                Shape: (num_shape_functions, n_points)

        Notes:
            Uses product rule with second y-derivatives of test functions in
            x-direction and values in y-direction.
        """
        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
        test_x = self.test_fcnx(num_shape_func_in_1d, xi)
        grad_grad_y = self.dtest_fcn(num_shape_func_in_1d, eta)[1]
        values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

        for i in range(num_shape_func_in_1d):
            values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
                test_x[i, :] * grad_grad_y
            )

        return values

dtest_fcn(n_test, x)

Computes first and second derivatives of test functions.

Calculates derivatives of test functions constructed from Jacobi polynomials, handling special cases for n=1,2 separately.

Parameters:

Name	Type	Description	Default
n_test	int	Number of test functions	required
x	ndarray	Points at which to evaluate derivatives Shape: (n_points,)	required

Returns:

Name	Type	Description
tuple	(ndarray, ndarray)	First and second derivatives First element: First derivatives, shape (n_test, n_points) Second element: Second derivatives, shape (n_test, n_points)

Notes

Special cases for n=1,2 ensure proper derivative calculations following hp-VPINNs methodology.

Source code in scirex/core/sciml/fe/basis_2d_qn_chebyshev_2.py

def dtest_fcn(self, n_test: int, x: np.ndarray) -> np.ndarray:
    """Computes first and second derivatives of test functions.

    Calculates derivatives of test functions constructed from Jacobi
    polynomials, handling special cases for n=1,2 separately.

    Args:
        n_test: Number of test functions
        x: Points at which to evaluate derivatives
            Shape: (n_points,)

    Returns:
        tuple(np.ndarray, np.ndarray): First and second derivatives
            First element: First derivatives, shape (n_test, n_points)
            Second element: Second derivatives, shape (n_test, n_points)

    Notes:
        Special cases for n=1,2 ensure proper derivative calculations
        following hp-VPINNs methodology.
    """
    d1test_total = []
    d2test_total = []
    for n in range(1, n_test + 1):
        if n == 1:
            d1test = (
                ((n + 1) / 2)
                * self.jacobi_wrapper(n, 1 / 2, 1 / 2, x)
                / self.jacobi_wrapper(n + 1, -1 / 2, -1 / 2, 1)
            )
            d2test = (
                ((n + 2) * (n + 1) / (2 * 2))
                * self.jacobi_wrapper(n - 1, 3 / 2, 3 / 2, x)
                / self.jacobi_wrapper(n + 1, -1 / 2, -1 / 2, 1)
            )
            d1test_total.append(d1test)
            d2test_total.append(d2test)
        elif n == 2:
            d1test = ((n + 1) / 2) * self.jacobi_wrapper(
                n, 1 / 2, 1 / 2, x
            ) / self.jacobi_wrapper(n + 1, -1 / 2, -1 / 2, 1) - (
                (n - 1) / 2
            ) * self.jacobi_wrapper(
                n - 2, 1 / 2, 1 / 2, x
            ) / self.jacobi_wrapper(
                n - 1, -1 / 2, -1 / 2, 1
            )
            d2test = (
                ((n + 2) * (n + 1) / (2 * 2))
                * self.jacobi_wrapper(n - 1, 3 / 2, 3 / 2, x)
                / self.jacobi_wrapper(n + 1, -1 / 2, -1 / 2, 1)
            )
            d1test_total.append(d1test)
            d2test_total.append(d2test)
        else:
            d1test = ((n + 1) / 2) * self.jacobi_wrapper(
                n, 1 / 2, 1 / 2, x
            ) / self.jacobi_wrapper(n + 1, -1 / 2, -1 / 2, 1) - (
                (n - 1) / 2
            ) * self.jacobi_wrapper(
                n - 2, 1 / 2, 1 / 2, x
            ) / self.jacobi_wrapper(
                n - 1, -1 / 2, -1 / 2, 1
            )
            d2test = ((n + 2) * (n + 1) / (2 * 2)) * self.jacobi_wrapper(
                n - 1, 3 / 2, 3 / 2, x
            ) / self.jacobi_wrapper(n + 1, -1 / 2, -1 / 2, 1) - (
                (n) * (n - 1) / (2 * 2)
            ) * self.jacobi_wrapper(
                n - 3, 3 / 2, 3 / 2, x
            ) / self.jacobi_wrapper(
                n - 1, -1 / 2, -1 / 2, 1
            )
            d1test_total.append(d1test)
            d2test_total.append(d2test)
    return np.asarray(d1test_total), np.asarray(d2test_total)

gradx(xi, eta)

Computes x-derivatives of basis functions.

Evaluates partial derivatives with respect to x of all basis functions at given coordinates.

Parameters:

Name	Type	Description	Default
xi	ndarray	x-coordinates at which to evaluate derivatives Shape: (n_points,)	required
eta	ndarray	y-coordinates at which to evaluate derivatives Shape: (n_points,)	required

Returns:

Type	Description
ndarray	np.ndarray: Values of x-derivatives Shape: (num_shape_functions, n_points)

Notes

Uses product rule with x-derivatives of test functions in x-direction and values in y-direction.

Source code in scirex/core/sciml/fe/basis_2d_qn_chebyshev_2.py

def gradx(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """Computes x-derivatives of basis functions.

    Evaluates partial derivatives with respect to x of all basis
    functions at given coordinates.

    Args:
        xi: x-coordinates at which to evaluate derivatives
            Shape: (n_points,)
        eta: y-coordinates at which to evaluate derivatives
            Shape: (n_points,)

    Returns:
        np.ndarray: Values of x-derivatives
            Shape: (num_shape_functions, n_points)

    Notes:
        Uses product rule with x-derivatives of test functions in
        x-direction and values in y-direction.
    """
    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
    grad_test_x = self.dtest_fcn(num_shape_func_in_1d, xi)[0]
    test_y = self.test_fcny(num_shape_func_in_1d, eta)
    values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

    for i in range(num_shape_func_in_1d):
        values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
            grad_test_x[i, :] * test_y
        )

    return values

gradxx(xi, eta)

Computes second x-derivatives of basis functions.

Evaluates second partial derivatives with respect to x of all basis functions at given coordinates.

Parameters:

Name	Type	Description	Default
xi	ndarray	x-coordinates at which to evaluate derivatives Shape: (n_points,)	required
eta	ndarray	y-coordinates at which to evaluate derivatives Shape: (n_points,)	required

Returns:

Type	Description
ndarray	np.ndarray: Values of second x-derivatives Shape: (num_shape_functions, n_points)

Notes

Uses product rule with second x-derivatives of test functions in x-direction and values in y-direction.

Source code in scirex/core/sciml/fe/basis_2d_qn_chebyshev_2.py

def gradxx(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """Computes second x-derivatives of basis functions.

    Evaluates second partial derivatives with respect to x of all basis
    functions at given coordinates.

    Args:
        xi: x-coordinates at which to evaluate derivatives
            Shape: (n_points,)
        eta: y-coordinates at which to evaluate derivatives
            Shape: (n_points,)

    Returns:
        np.ndarray: Values of second x-derivatives
            Shape: (num_shape_functions, n_points)

    Notes:
        Uses product rule with second x-derivatives of test functions in
        x-direction and values in y-direction.
    """
    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
    grad_grad_x = self.dtest_fcn(num_shape_func_in_1d, xi)[1]
    test_y = self.test_fcny(num_shape_func_in_1d, eta)
    values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

    for i in range(num_shape_func_in_1d):
        values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
            grad_grad_x[i, :] * test_y
        )

    return values

gradxy(xi, eta)

Computes second x-derivatives of basis functions.

Evaluates second partial derivatives with respect to x of all basis functions at given coordinates.

Parameters:

Name	Type	Description	Default
xi	ndarray	x-coordinates at which to evaluate derivatives Shape: (n_points,)	required
eta	ndarray	y-coordinates at which to evaluate derivatives Shape: (n_points,)	required

Returns:

Type	Description
ndarray	np.ndarray: Values of second x-derivatives Shape: (num_shape_functions, n_points)

Notes

Uses product rule with second x-derivatives of test functions in x-direction and y derivative values in y-direction.

Source code in scirex/core/sciml/fe/basis_2d_qn_chebyshev_2.py

def gradxy(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """Computes second x-derivatives of basis functions.

    Evaluates second partial derivatives with respect to x of all basis
    functions at given coordinates.

    Args:
        xi: x-coordinates at which to evaluate derivatives
            Shape: (n_points,)
        eta: y-coordinates at which to evaluate derivatives
            Shape: (n_points,)

    Returns:
        np.ndarray: Values of second x-derivatives
            Shape: (num_shape_functions, n_points)

    Notes:
        Uses product rule with second x-derivatives of test functions in
        x-direction and y derivative values in y-direction.
    """
    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
    grad_test_x = self.dtest_fcn(num_shape_func_in_1d, xi)[0]
    grad_test_y = self.dtest_fcn(num_shape_func_in_1d, eta)[0]
    values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

    for i in range(num_shape_func_in_1d):
        values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
            grad_test_x[i, :] * grad_test_y
        )

    return values

grady(xi, eta)

Computes y-derivatives of basis functions.

Evaluates partial derivatives with respect to y of all basis functions at given coordinates.

Parameters:

Name	Type	Description	Default
xi	ndarray	x-coordinates at which to evaluate derivatives Shape: (n_points,)	required
eta	ndarray	y-coordinates at which to evaluate derivatives Shape: (n_points,)	required

Returns:

Type	Description
ndarray	np.ndarray: Values of y-derivatives Shape: (num_shape_functions, n_points)

Notes

Uses product rule with values in x-direction and y-derivatives of test functions in y-direction.

Source code in scirex/core/sciml/fe/basis_2d_qn_chebyshev_2.py

def grady(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """Computes y-derivatives of basis functions.

    Evaluates partial derivatives with respect to y of all basis
    functions at given coordinates.

    Args:
        xi: x-coordinates at which to evaluate derivatives
            Shape: (n_points,)
        eta: y-coordinates at which to evaluate derivatives
            Shape: (n_points,)

    Returns:
        np.ndarray: Values of y-derivatives
            Shape: (num_shape_functions, n_points)

    Notes:
        Uses product rule with values in x-direction and y-derivatives
        of test functions in y-direction.
    """
    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
    test_x = self.test_fcnx(num_shape_func_in_1d, xi)
    grad_test_y = self.dtest_fcn(num_shape_func_in_1d, eta)[0]
    values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

    for i in range(num_shape_func_in_1d):
        values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
            test_x[i, :] * grad_test_y
        )

    return values

gradyy(xi, eta)

Computes second x-derivatives of basis functions.

Evaluates second partial derivatives with respect to x of all basis functions at given coordinates.

Parameters:

Name	Type	Description	Default
xi	ndarray	x-coordinates at which to evaluate derivatives Shape: (n_points,)	required
eta	ndarray	y-coordinates at which to evaluate derivatives Shape: (n_points,)	required

Returns:

Type	Description
ndarray	np.ndarray: Values of second x-derivatives Shape: (num_shape_functions, n_points)

Notes

Uses product rule with second y-derivatives of test functions in x-direction and values in y-direction.

Source code in scirex/core/sciml/fe/basis_2d_qn_chebyshev_2.py

def gradyy(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """Computes second x-derivatives of basis functions.

    Evaluates second partial derivatives with respect to x of all basis
    functions at given coordinates.

    Args:
        xi: x-coordinates at which to evaluate derivatives
            Shape: (n_points,)
        eta: y-coordinates at which to evaluate derivatives
            Shape: (n_points,)

    Returns:
        np.ndarray: Values of second x-derivatives
            Shape: (num_shape_functions, n_points)

    Notes:
        Uses product rule with second y-derivatives of test functions in
        x-direction and values in y-direction.
    """
    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
    test_x = self.test_fcnx(num_shape_func_in_1d, xi)
    grad_grad_y = self.dtest_fcn(num_shape_func_in_1d, eta)[1]
    values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

    for i in range(num_shape_func_in_1d):
        values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
            test_x[i, :] * grad_grad_y
        )

    return values

jacobi_wrapper(n, a, b, x)

Evaluates Jacobi polynomial at specified points.

Computes values of nth degree Jacobi polynomial with parameters (a,b) at given points x.

Parameters:

Name	Type	Description	Default
n	int	Degree of Jacobi polynomial. Must be non-negative integer.	required
a	int	First parameter of Jacobi polynomial	required
b	int	Second parameter of Jacobi polynomial	required
x	ndarray	Points at which to evaluate polynomial Shape: (n_points,)	required

Returns:

Type	Description
ndarray	np.ndarray: Values of Jacobi polynomial at input points Shape: Same as input x

Notes

Wrapper around scipy.special.jacobi that ensures float64 precision and proper array handling.

Source code in scirex/core/sciml/fe/basis_2d_qn_chebyshev_2.py

def jacobi_wrapper(self, n: int, a: int, b: int, x: np.ndarray) -> np.ndarray:
    """Evaluates Jacobi polynomial at specified points.

    Computes values of nth degree Jacobi polynomial with parameters (a,b)
    at given points x.

    Args:
        n: Degree of Jacobi polynomial. Must be non-negative integer.
        a: First parameter of Jacobi polynomial
        b: Second parameter of Jacobi polynomial
        x: Points at which to evaluate polynomial
            Shape: (n_points,)

    Returns:
        np.ndarray: Values of Jacobi polynomial at input points
            Shape: Same as input x

    Notes:
        Wrapper around scipy.special.jacobi that ensures float64 precision
        and proper array handling.
    """
    x = np.array(x, dtype=np.float64)
    return jacobi(n, a, b)(x)

test_fcnx(n_test, x)

Computes x-component test functions.

Evaluates the x-direction test functions constructed as differences of normalized Jacobi polynomials.

Parameters:

Name	Type	Description	Default
n_test	int	Number of test functions to compute	required
x	ndarray	Points at which to evaluate functions Shape: (n_points,)	required

Returns:

Type	Description
ndarray	np.ndarray: Values of test functions at input points Shape: (n_test, n_points)

Notes

Test functions are constructed as differences of normalized Jacobi polynomials following hp-VPINNs methodology.

Source code in scirex/core/sciml/fe/basis_2d_qn_chebyshev_2.py

def test_fcnx(self, n_test: int, x: np.ndarray) -> np.ndarray:
    """Computes x-component test functions.

    Evaluates the x-direction test functions constructed as differences
    of normalized Jacobi polynomials.

    Args:
        n_test: Number of test functions to compute
        x: Points at which to evaluate functions
            Shape: (n_points,)

    Returns:
        np.ndarray: Values of test functions at input points
            Shape: (n_test, n_points)

    Notes:
        Test functions are constructed as differences of normalized Jacobi
        polynomials following hp-VPINNs methodology.
    """
    test_total = []
    for n in range(1, n_test + 1):
        test = self.jacobi_wrapper(n + 1, -1 / 2, -1 / 2, x) / self.jacobi_wrapper(
            n + 1, -1 / 2, -1 / 2, 1
        ) - self.jacobi_wrapper(n - 1, -1 / 2, -1 / 2, x) / self.jacobi_wrapper(
            n - 1, -1 / 2, -1 / 2, 1
        )
        test_total.append(test)
    return np.asarray(test_total, np.float64)

test_fcny(n_test, y)

Computes y-component test functions.

Evaluates the y-direction test functions constructed as differences of normalized Jacobi polynomials.

Parameters:

Name	Type	Description	Default
n_test	int	Number of test functions to compute	required
y	ndarray	Points at which to evaluate functions Shape: (n_points,)	required

Returns:

Type	Description
ndarray	np.ndarray: Values of test functions at input points Shape: (n_test, n_points)

Notes

Test functions are constructed as differences of normalized Jacobi polynomials following hp-VPINNs methodology.

Source code in scirex/core/sciml/fe/basis_2d_qn_chebyshev_2.py

def test_fcny(self, n_test: int, y: np.ndarray) -> np.ndarray:
    """Computes y-component test functions.

    Evaluates the y-direction test functions constructed as differences
    of normalized Jacobi polynomials.

    Args:
        n_test: Number of test functions to compute
        y: Points at which to evaluate functions
            Shape: (n_points,)

    Returns:
        np.ndarray: Values of test functions at input points
            Shape: (n_test, n_points)

    Notes:
        Test functions are constructed as differences of normalized Jacobi
        polynomials following hp-VPINNs methodology.
    """
    test_total = []
    for n in range(1, n_test + 1):
        test = self.jacobi_wrapper(n + 1, -1 / 2, -1 / 2, y) / self.jacobi_wrapper(
            n + 1, -1 / 2, -1 / 2, 1
        ) - self.jacobi_wrapper(n - 1, -1 / 2, -1 / 2, y) / self.jacobi_wrapper(
            n - 1, -1 / 2, -1 / 2, 1
        )
        test_total.append(test)
    return np.asarray(test_total, np.float64)

value(xi, eta)

Evaluates basis functions at given coordinates.

Computes values of all basis functions at specified (xi,eta) points using tensor product of 1D test functions.

Parameters:

Name	Type	Description	Default
xi	ndarray	x-coordinates at which to evaluate functions Shape: (n_points,)	required
eta	ndarray	y-coordinates at which to evaluate functions Shape: (n_points,)	required

Returns:

Type	Description
ndarray	np.ndarray: Values of all basis functions Shape: (num_shape_functions, n_points)

Notes

Basis functions are constructed as products of 1D test functions in x and y directions.

Source code in scirex/core/sciml/fe/basis_2d_qn_chebyshev_2.py

def value(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """Evaluates basis functions at given coordinates.

    Computes values of all basis functions at specified (xi,eta) points
    using tensor product of 1D test functions.

    Args:
        xi: x-coordinates at which to evaluate functions
            Shape: (n_points,)
        eta: y-coordinates at which to evaluate functions
            Shape: (n_points,)

    Returns:
        np.ndarray: Values of all basis functions
            Shape: (num_shape_functions, n_points)

    Notes:
        Basis functions are constructed as products of 1D test functions
        in x and y directions.
    """
    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
    test_x = self.test_fcnx(num_shape_func_in_1d, xi)
    test_y = self.test_fcny(num_shape_func_in_1d, eta)
    values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

    for i in range(num_shape_func_in_1d):
        values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
            test_x[i, :] * test_y
        )

    return values

Basis2DQNJacobi

Bases: BasisFunction2D

A specialized implementation of two-dimensional basis functions using Jacobi polynomials for QN elements.

This class provides a complete implementation for computing basis functions and their derivatives in two dimensions, specifically designed for use in variational physics-informed neural networks (VPINNs) with domain decomposition. The basis functions are constructed using Jacobi polynomials with efficient derivative computations.

The class inherits from BasisFunction2D and implements all required methods for computing function values and first/second order derivatives. The implementation follows the methodology described in hp-VPINNs research.

Attributes:

Name	Type	Description
num_shape_functions	int	Total number of shape functions in the 2D element. Must be a perfect square as it represents tensor product of 1D functions.

Methods:

Name	Description
jacobi_wrapper	Evaluates Jacobi polynomial at given points
djacobi	Computes kth derivative of Jacobi polynomial
test_fcnx	Computes x-component test functions
test_fcny	Computes y-component test functions
dtest_fcn	Computes first derivatives of test functions
ddtest_fcn	Computes second derivatives of test functions
value	Computes values of all basis functions
gradx	Computes x-derivatives of all basis functions
grady	Computes y-derivatives of all basis functions
gradxx	Computes second x-derivatives of all basis functions
gradyy	Computes second y-derivatives of all basis functions
gradxy	Computes mixed xy-derivatives of all basis functions

Implementation Details

• Basis functions are constructed as tensor products of 1D test functions
• Test functions are derived from Jacobi polynomials with parameters (0,0)
• All computations maintain double precision (float64)
• Efficient vectorized operations using numpy arrays

Example

basis = Basis2DQNJacobi(num_shape_functions=16)  # Creates 4x4 basis functions
xi = np.linspace(-1, 1, 100)
eta = np.linspace(-1, 1, 100)
values = basis.value(xi, eta)
x_derivatives = basis.gradx(xi, eta)

Source code in scirex/core/sciml/fe/basis_2d_qn_jacobi.py

class Basis2DQNJacobi(BasisFunction2D):
    """
    A specialized implementation of two-dimensional basis functions using Jacobi polynomials for QN elements.

    This class provides a complete implementation for computing basis functions and their derivatives
    in two dimensions, specifically designed for use in variational physics-informed neural networks
    (VPINNs) with domain decomposition. The basis functions are constructed using Jacobi polynomials
    with efficient derivative computations.

    The class inherits from BasisFunction2D and implements all required methods for computing
    function values and first/second order derivatives. The implementation follows the methodology
    described in hp-VPINNs research.

    Attributes:
        num_shape_functions (int): Total number of shape functions in the 2D element.
            Must be a perfect square as it represents tensor product of 1D functions.

    Methods:
        jacobi_wrapper(n, a, b, x): Evaluates Jacobi polynomial at given points
        djacobi(n, a, b, x, k): Computes kth derivative of Jacobi polynomial
        test_fcnx(n_test, x): Computes x-component test functions
        test_fcny(n_test, y): Computes y-component test functions
        dtest_fcn(n_test, x): Computes first derivatives of test functions
        ddtest_fcn(n_test, x): Computes second derivatives of test functions
        value(xi, eta): Computes values of all basis functions
        gradx(xi, eta): Computes x-derivatives of all basis functions
        grady(xi, eta): Computes y-derivatives of all basis functions
        gradxx(xi, eta): Computes second x-derivatives of all basis functions
        gradyy(xi, eta): Computes second y-derivatives of all basis functions
        gradxy(xi, eta): Computes mixed xy-derivatives of all basis functions

    Implementation Details:
        - Basis functions are constructed as tensor products of 1D test functions
        - Test functions are derived from Jacobi polynomials with parameters (0,0)
        - All computations maintain double precision (float64)
        - Efficient vectorized operations using numpy arrays

    Example:
        ```python
        basis = Basis2DQNJacobi(num_shape_functions=16)  # Creates 4x4 basis functions
        xi = np.linspace(-1, 1, 100)
        eta = np.linspace(-1, 1, 100)
        values = basis.value(xi, eta)
        x_derivatives = basis.gradx(xi, eta)
        ```
    """

    def __init__(self, num_shape_functions: int):
        super().__init__(num_shape_functions)

    def jacobi_wrapper(self, n: int, a: int, b: int, x: np.ndarray) -> np.ndarray:
        """
        Evaluate the Jacobi polynomial of degree `n` with parameters `a` and `b` at the given points `x`.

        Args:
            n (int): Degree of the Jacobi polynomial.
            a (float): First parameter of the Jacobi polynomial.
            b (float): Second parameter of the Jacobi polynomial.
            x (np.ndarray): Points at which to evaluate the Jacobi polynomial.

        Returns:
            np.ndarray: Values of the Jacobi polynomial at the given points `x`.
        """

        x = np.array(x, dtype=np.float64)
        return jacobi(n, a, b)(x)

    # Derivative of the Jacobi polynomials
    def djacobi(self, n: int, a: int, b: int, x: np.ndarray, k: int) -> np.ndarray:
        """
        Evaluate the k-th derivative of the Jacobi polynomial of degree n with parameters a and b at the given points x.

        Args:
            n (int): Degree of the Jacobi polynomial.
            a (float): First parameter of the Jacobi polynomial.
            b (float): Second parameter of the Jacobi polynomial.
            x (np.ndarray): Points at which to evaluate the Jacobi polynomial.
            k (int): Order of the derivative.

        Returns:
            np.ndarray: Values of the k-th derivative of the Jacobi polynomial at the given points x.

        Raises:
            ValueError: If the derivative order is not 1 or 2
        """
        x = np.array(x, dtype=np.float64)
        if k == 1:
            return jacobi(n, a, b).deriv()(x)
        if k == 2:
            return jacobi(n, a, b).deriv(2)(x)
        else:
            print(f"Invalid derivative order {k} in {__name__}.")
            raise ValueError("Derivative order should be 1 or 2.")

    ## Helper Function
    def test_fcnx(self, n_test: int, x: np.ndarray) -> np.ndarray:
        """
        Compute the x-component of the test functions for a given number of test functions and x-coordinates.

        Args:
            n_test (int): Number of test functions.
            x (np.ndarray): x-coordinates at which to evaluate the test functions.

        Returns:
            np.ndarray: Values of the x-component of the test functions.

        """
        test_total = []
        for n in range(1, n_test + 1):
            test = self.jacobi_wrapper(n - 1, 0, 0, x)
            test_total.append(test)
        return np.asarray(test_total, np.float64)

    def test_fcny(self, n_test: int, y: np.ndarray) -> np.ndarray:
        """
        Compute the y-component of the test functions for a given number of test functions and y-coordinates.

        Args:
            n_test (int): Number of test functions.
            y (np.ndarray): y-coordinates at which to evaluate the test functions.

        Returns:
            np.ndarray: Values of the y-component of the test functions.

        """
        test_total = []
        for n in range(1, n_test + 1):
            test = self.jacobi_wrapper(n - 1, 0, 0, y)
            test_total.append(test)
        return np.asarray(test_total, np.float64)

    def dtest_fcn(self, n_test: int, x: np.ndarray) -> np.ndarray:
        """
        Compute the x-derivatives of the test functions for a given number of test functions and x-coordinates.

        Args:
            n_test (int): Number of test functions.
            x (np.ndarray): x-coordinates at which to evaluate the test functions.

        Returns:
            np.ndarray: Values of the x-derivatives of the test functions.
        """
        d1test_total = []
        for n in range(1, n_test + 1):
            d1test = self.djacobi(n - 1, 0, 0, x, 1)
            d1test_total.append(d1test)
        return np.asarray(d1test_total)

    def ddtest_fcn(self, n_test: int, x: np.ndarray) -> np.ndarray:
        """
        Compute the x-derivatives of the test functions for a given number of test functions and x-coordinates.

        Args:
            n_test (int): Number of test functions.
            x (np.ndarray): x-coordinates at which to evaluate the test functions.

        Returns:
            np.ndarray: Values of the x-derivatives of the test functions.
        """
        d1test_total = []
        for n in range(1, n_test + 1):
            d1test = self.djacobi(n - 1, 0, 0, x, 2)
            d1test_total.append(d1test)
        return np.asarray(d1test_total)

    def value(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        This method returns the values of the basis functions at the given (xi, eta) coordinates.

        Args:
            xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
            eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

        Returns:
            np.ndarray: Values of the basis functions.
        """
        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
        test_x = self.test_fcnx(num_shape_func_in_1d, xi)
        test_y = self.test_fcny(num_shape_func_in_1d, eta)
        values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

        for i in range(num_shape_func_in_1d):
            values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
                test_x[i, :] * test_y
            )

        return values

    def gradx(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        This method returns the x-derivatives of the basis functions at the given (xi, eta) coordinates.

        Args:
            xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
            eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

        Returns:
            np.ndarray: Values of the x-derivatives of the basis functions.
        """
        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
        grad_test_x = self.dtest_fcn(num_shape_func_in_1d, xi)
        test_y = self.test_fcny(num_shape_func_in_1d, eta)
        values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

        for i in range(num_shape_func_in_1d):
            values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
                grad_test_x[i, :] * test_y
            )

        return values

    def grady(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        This method returns the y-derivatives of the basis functions at the given (xi, eta) coordinates.

        Args:
            xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
            eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

        Returns:
            np.ndarray: Values of the y-derivatives of the basis functions.
        """
        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
        test_x = self.test_fcnx(num_shape_func_in_1d, xi)
        grad_test_y = self.dtest_fcn(num_shape_func_in_1d, eta)
        values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

        for i in range(num_shape_func_in_1d):
            values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
                test_x[i, :] * grad_test_y
            )

        return values

    def gradxx(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        This method returns the xx-derivatives of the basis functions at the given (xi, eta) coordinates.

        Args:
            xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
            eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

        Returns:
            np.ndarray: Values of the xx-derivatives of the basis functions.
        """
        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
        grad_grad_x = self.ddtest_fcn(num_shape_func_in_1d, xi)
        test_y = self.test_fcny(num_shape_func_in_1d, eta)
        values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

        for i in range(num_shape_func_in_1d):
            values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
                grad_grad_x[i, :] * test_y
            )

        return values

    def gradxy(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        This method returns the xy-derivatives of the basis functions at the given (xi, eta) coordinates.

        Args:
            xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
            eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

        Returns:
            np.ndarray: Values of the xy-derivatives of the basis functions.
        """
        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
        grad_test_x = self.dtest_fcn(num_shape_func_in_1d, xi)
        grad_test_y = self.dtest_fcn(num_shape_func_in_1d, eta)
        values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

        for i in range(num_shape_func_in_1d):
            values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
                grad_test_x[i, :] * grad_test_y
            )

        return values

    def gradyy(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        This method returns the yy-derivatives of the basis functions at the given (xi, eta) coordinates.

        Args:
            xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
            eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

        Returns:
            np.ndarray: Values of the yy-derivatives of the basis functions.
        """
        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
        test_x = self.test_fcnx(num_shape_func_in_1d, xi)
        grad_grad_y = self.ddtest_fcn(num_shape_func_in_1d, eta)
        values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

        for i in range(num_shape_func_in_1d):
            values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
                test_x[i, :] * grad_grad_y
            )

        return values

ddtest_fcn(n_test, x)

Compute the x-derivatives of the test functions for a given number of test functions and x-coordinates.

Parameters:

Name	Type	Description	Default
n_test	int	Number of test functions.	required
x	ndarray	x-coordinates at which to evaluate the test functions.	required

Returns:

Type	Description
ndarray	np.ndarray: Values of the x-derivatives of the test functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_jacobi.py

def ddtest_fcn(self, n_test: int, x: np.ndarray) -> np.ndarray:
    """
    Compute the x-derivatives of the test functions for a given number of test functions and x-coordinates.

    Args:
        n_test (int): Number of test functions.
        x (np.ndarray): x-coordinates at which to evaluate the test functions.

    Returns:
        np.ndarray: Values of the x-derivatives of the test functions.
    """
    d1test_total = []
    for n in range(1, n_test + 1):
        d1test = self.djacobi(n - 1, 0, 0, x, 2)
        d1test_total.append(d1test)
    return np.asarray(d1test_total)

djacobi(n, a, b, x, k)

Evaluate the k-th derivative of the Jacobi polynomial of degree n with parameters a and b at the given points x.

Parameters:

Name	Type	Description	Default
n	int	Degree of the Jacobi polynomial.	required
a	float	First parameter of the Jacobi polynomial.	required
b	float	Second parameter of the Jacobi polynomial.	required
x	ndarray	Points at which to evaluate the Jacobi polynomial.	required
k	int	Order of the derivative.	required

Returns:

Type	Description
ndarray	np.ndarray: Values of the k-th derivative of the Jacobi polynomial at the given points x.

Raises:

Type	Description
ValueError	If the derivative order is not 1 or 2

Source code in scirex/core/sciml/fe/basis_2d_qn_jacobi.py

def djacobi(self, n: int, a: int, b: int, x: np.ndarray, k: int) -> np.ndarray:
    """
    Evaluate the k-th derivative of the Jacobi polynomial of degree n with parameters a and b at the given points x.

    Args:
        n (int): Degree of the Jacobi polynomial.
        a (float): First parameter of the Jacobi polynomial.
        b (float): Second parameter of the Jacobi polynomial.
        x (np.ndarray): Points at which to evaluate the Jacobi polynomial.
        k (int): Order of the derivative.

    Returns:
        np.ndarray: Values of the k-th derivative of the Jacobi polynomial at the given points x.

    Raises:
        ValueError: If the derivative order is not 1 or 2
    """
    x = np.array(x, dtype=np.float64)
    if k == 1:
        return jacobi(n, a, b).deriv()(x)
    if k == 2:
        return jacobi(n, a, b).deriv(2)(x)
    else:
        print(f"Invalid derivative order {k} in {__name__}.")
        raise ValueError("Derivative order should be 1 or 2.")

dtest_fcn(n_test, x)

Compute the x-derivatives of the test functions for a given number of test functions and x-coordinates.

Parameters:

Name	Type	Description	Default
n_test	int	Number of test functions.	required
x	ndarray	x-coordinates at which to evaluate the test functions.	required

Returns:

Type	Description
ndarray	np.ndarray: Values of the x-derivatives of the test functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_jacobi.py

def dtest_fcn(self, n_test: int, x: np.ndarray) -> np.ndarray:
    """
    Compute the x-derivatives of the test functions for a given number of test functions and x-coordinates.

    Args:
        n_test (int): Number of test functions.
        x (np.ndarray): x-coordinates at which to evaluate the test functions.

    Returns:
        np.ndarray: Values of the x-derivatives of the test functions.
    """
    d1test_total = []
    for n in range(1, n_test + 1):
        d1test = self.djacobi(n - 1, 0, 0, x, 1)
        d1test_total.append(d1test)
    return np.asarray(d1test_total)

gradx(xi, eta)

This method returns the x-derivatives of the basis functions at the given (xi, eta) coordinates.

Parameters:

Name	Type	Description	Default
xi	ndarray	x-coordinates at which to evaluate the basis functions.	required
eta	ndarray	y-coordinates at which to evaluate the basis functions.	required

Returns:

Type	Description
ndarray	np.ndarray: Values of the x-derivatives of the basis functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_jacobi.py

def gradx(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    This method returns the x-derivatives of the basis functions at the given (xi, eta) coordinates.

    Args:
        xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
        eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

    Returns:
        np.ndarray: Values of the x-derivatives of the basis functions.
    """
    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
    grad_test_x = self.dtest_fcn(num_shape_func_in_1d, xi)
    test_y = self.test_fcny(num_shape_func_in_1d, eta)
    values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

    for i in range(num_shape_func_in_1d):
        values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
            grad_test_x[i, :] * test_y
        )

    return values

gradxx(xi, eta)

This method returns the xx-derivatives of the basis functions at the given (xi, eta) coordinates.

Parameters:

Name	Type	Description	Default
xi	ndarray	x-coordinates at which to evaluate the basis functions.	required
eta	ndarray	y-coordinates at which to evaluate the basis functions.	required

Returns:

Type	Description
ndarray	np.ndarray: Values of the xx-derivatives of the basis functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_jacobi.py

def gradxx(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    This method returns the xx-derivatives of the basis functions at the given (xi, eta) coordinates.

    Args:
        xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
        eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

    Returns:
        np.ndarray: Values of the xx-derivatives of the basis functions.
    """
    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
    grad_grad_x = self.ddtest_fcn(num_shape_func_in_1d, xi)
    test_y = self.test_fcny(num_shape_func_in_1d, eta)
    values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

    for i in range(num_shape_func_in_1d):
        values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
            grad_grad_x[i, :] * test_y
        )

    return values

gradxy(xi, eta)

This method returns the xy-derivatives of the basis functions at the given (xi, eta) coordinates.

Parameters:

Name	Type	Description	Default
xi	ndarray	x-coordinates at which to evaluate the basis functions.	required
eta	ndarray	y-coordinates at which to evaluate the basis functions.	required

Returns:

Type	Description
ndarray	np.ndarray: Values of the xy-derivatives of the basis functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_jacobi.py

def gradxy(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    This method returns the xy-derivatives of the basis functions at the given (xi, eta) coordinates.

    Args:
        xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
        eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

    Returns:
        np.ndarray: Values of the xy-derivatives of the basis functions.
    """
    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
    grad_test_x = self.dtest_fcn(num_shape_func_in_1d, xi)
    grad_test_y = self.dtest_fcn(num_shape_func_in_1d, eta)
    values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

    for i in range(num_shape_func_in_1d):
        values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
            grad_test_x[i, :] * grad_test_y
        )

    return values

grady(xi, eta)

This method returns the y-derivatives of the basis functions at the given (xi, eta) coordinates.

Parameters:

Name	Type	Description	Default
xi	ndarray	x-coordinates at which to evaluate the basis functions.	required
eta	ndarray	y-coordinates at which to evaluate the basis functions.	required

Returns:

Type	Description
ndarray	np.ndarray: Values of the y-derivatives of the basis functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_jacobi.py

def grady(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    This method returns the y-derivatives of the basis functions at the given (xi, eta) coordinates.

    Args:
        xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
        eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

    Returns:
        np.ndarray: Values of the y-derivatives of the basis functions.
    """
    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
    test_x = self.test_fcnx(num_shape_func_in_1d, xi)
    grad_test_y = self.dtest_fcn(num_shape_func_in_1d, eta)
    values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

    for i in range(num_shape_func_in_1d):
        values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
            test_x[i, :] * grad_test_y
        )

    return values

gradyy(xi, eta)

This method returns the yy-derivatives of the basis functions at the given (xi, eta) coordinates.

Parameters:

Name	Type	Description	Default
xi	ndarray	x-coordinates at which to evaluate the basis functions.	required
eta	ndarray	y-coordinates at which to evaluate the basis functions.	required

Returns:

Type	Description
ndarray	np.ndarray: Values of the yy-derivatives of the basis functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_jacobi.py

def gradyy(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    This method returns the yy-derivatives of the basis functions at the given (xi, eta) coordinates.

    Args:
        xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
        eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

    Returns:
        np.ndarray: Values of the yy-derivatives of the basis functions.
    """
    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
    test_x = self.test_fcnx(num_shape_func_in_1d, xi)
    grad_grad_y = self.ddtest_fcn(num_shape_func_in_1d, eta)
    values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

    for i in range(num_shape_func_in_1d):
        values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
            test_x[i, :] * grad_grad_y
        )

    return values

jacobi_wrapper(n, a, b, x)

Evaluate the Jacobi polynomial of degree n with parameters a and b at the given points x.

Parameters:

Name	Type	Description	Default
n	int	Degree of the Jacobi polynomial.	required
a	float	First parameter of the Jacobi polynomial.	required
b	float	Second parameter of the Jacobi polynomial.	required
x	ndarray	Points at which to evaluate the Jacobi polynomial.	required

Returns:

Type	Description
ndarray	np.ndarray: Values of the Jacobi polynomial at the given points x.

Source code in scirex/core/sciml/fe/basis_2d_qn_jacobi.py

def jacobi_wrapper(self, n: int, a: int, b: int, x: np.ndarray) -> np.ndarray:
    """
    Evaluate the Jacobi polynomial of degree `n` with parameters `a` and `b` at the given points `x`.

    Args:
        n (int): Degree of the Jacobi polynomial.
        a (float): First parameter of the Jacobi polynomial.
        b (float): Second parameter of the Jacobi polynomial.
        x (np.ndarray): Points at which to evaluate the Jacobi polynomial.

    Returns:
        np.ndarray: Values of the Jacobi polynomial at the given points `x`.
    """

    x = np.array(x, dtype=np.float64)
    return jacobi(n, a, b)(x)

test_fcnx(n_test, x)

Compute the x-component of the test functions for a given number of test functions and x-coordinates.

Parameters:

Name	Type	Description	Default
n_test	int	Number of test functions.	required
x	ndarray	x-coordinates at which to evaluate the test functions.	required

Returns:

Type	Description
ndarray	np.ndarray: Values of the x-component of the test functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_jacobi.py

def test_fcnx(self, n_test: int, x: np.ndarray) -> np.ndarray:
    """
    Compute the x-component of the test functions for a given number of test functions and x-coordinates.

    Args:
        n_test (int): Number of test functions.
        x (np.ndarray): x-coordinates at which to evaluate the test functions.

    Returns:
        np.ndarray: Values of the x-component of the test functions.

    """
    test_total = []
    for n in range(1, n_test + 1):
        test = self.jacobi_wrapper(n - 1, 0, 0, x)
        test_total.append(test)
    return np.asarray(test_total, np.float64)

test_fcny(n_test, y)

Compute the y-component of the test functions for a given number of test functions and y-coordinates.

Parameters:

Name	Type	Description	Default
n_test	int	Number of test functions.	required
y	ndarray	y-coordinates at which to evaluate the test functions.	required

Returns:

Type	Description
ndarray	np.ndarray: Values of the y-component of the test functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_jacobi.py

def test_fcny(self, n_test: int, y: np.ndarray) -> np.ndarray:
    """
    Compute the y-component of the test functions for a given number of test functions and y-coordinates.

    Args:
        n_test (int): Number of test functions.
        y (np.ndarray): y-coordinates at which to evaluate the test functions.

    Returns:
        np.ndarray: Values of the y-component of the test functions.

    """
    test_total = []
    for n in range(1, n_test + 1):
        test = self.jacobi_wrapper(n - 1, 0, 0, y)
        test_total.append(test)
    return np.asarray(test_total, np.float64)

value(xi, eta)

This method returns the values of the basis functions at the given (xi, eta) coordinates.

Parameters:

Name	Type	Description	Default
xi	ndarray	x-coordinates at which to evaluate the basis functions.	required
eta	ndarray	y-coordinates at which to evaluate the basis functions.	required

Returns:

Type	Description
ndarray	np.ndarray: Values of the basis functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_jacobi.py

def value(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    This method returns the values of the basis functions at the given (xi, eta) coordinates.

    Args:
        xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
        eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

    Returns:
        np.ndarray: Values of the basis functions.
    """
    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
    test_x = self.test_fcnx(num_shape_func_in_1d, xi)
    test_y = self.test_fcny(num_shape_func_in_1d, eta)
    values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

    for i in range(num_shape_func_in_1d):
        values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
            test_x[i, :] * test_y
        )

    return values

Basis2DQNLegendre

Bases: BasisFunction2D

A specialized implementation of two-dimensional basis functions using Legendre polynomials for Q1 elements.

This class provides a complete implementation for computing basis functions and their derivatives in two dimensions, specifically designed for use in variational physics-informed neural networks (VPINNs) with domain decomposition. The basis functions are constructed using Legendre polynomials implemented through Jacobi polynomial representations with parameters (0,0).

The class inherits from BasisFunction2D and implements all required methods for computing function values and derivatives up to second order.

Attributes:

Name	Type	Description
num_shape_functions	int	Total number of shape functions in the 2D element. Must be a perfect square as it represents tensor product of 1D functions.

Methods:

Name	Description
jacobi_wrapper	Evaluates Jacobi polynomial at given points
test_fcnx	Computes x-component test functions
test_fcny	Computes y-component test functions
dtest_fcn	Computes first and second derivatives of test functions
value	Computes values of all basis functions
gradx	Computes x-derivatives of all basis functions
grady	Computes y-derivatives of all basis functions
gradxx	Computes second x-derivatives of all basis functions
gradyy	Computes second y-derivatives of all basis functions
gradxy	Computes mixed xy-derivatives of all basis functions

Implementation Details

• Basis functions are constructed as tensor products of 1D test functions
• Test functions use Legendre polynomials via Jacobi polynomials with (0,0) parameters
• Special cases handled for n=1,2 in derivative calculations
• All computations maintain double precision (float64)
• Efficient vectorized operations using numpy arrays

Example

basis = Basis2DQNLegendre(num_shape_functions=16)  # Creates 4x4 basis functions
xi = np.linspace(-1, 1, 100)
eta = np.linspace(-1, 1, 100)
values = basis.value(xi, eta)
x_derivatives = basis.gradx(xi, eta)

Source code in scirex/core/sciml/fe/basis_2d_qn_legendre.py

class Basis2DQNLegendre(BasisFunction2D):
    """
    A specialized implementation of two-dimensional basis functions using Legendre polynomials for Q1 elements.

    This class provides a complete implementation for computing basis functions and their derivatives
    in two dimensions, specifically designed for use in variational physics-informed neural networks
    (VPINNs) with domain decomposition. The basis functions are constructed using Legendre polynomials
    implemented through Jacobi polynomial representations with parameters (0,0).

    The class inherits from BasisFunction2D and implements all required methods for computing
    function values and derivatives up to second order.

    Attributes:
        num_shape_functions (int): Total number of shape functions in the 2D element.
            Must be a perfect square as it represents tensor product of 1D functions.

    Methods:
        jacobi_wrapper(n, a, b, x): Evaluates Jacobi polynomial at given points
        test_fcnx(n_test, x): Computes x-component test functions
        test_fcny(n_test, y): Computes y-component test functions
        dtest_fcn(n_test, x): Computes first and second derivatives of test functions
        value(xi, eta): Computes values of all basis functions
        gradx(xi, eta): Computes x-derivatives of all basis functions
        grady(xi, eta): Computes y-derivatives of all basis functions
        gradxx(xi, eta): Computes second x-derivatives of all basis functions
        gradyy(xi, eta): Computes second y-derivatives of all basis functions
        gradxy(xi, eta): Computes mixed xy-derivatives of all basis functions

    Implementation Details:
        - Basis functions are constructed as tensor products of 1D test functions
        - Test functions use Legendre polynomials via Jacobi polynomials with (0,0) parameters
        - Special cases handled for n=1,2 in derivative calculations
        - All computations maintain double precision (float64)
        - Efficient vectorized operations using numpy arrays

    Example:
        ```python
        basis = Basis2DQNLegendre(num_shape_functions=16)  # Creates 4x4 basis functions
        xi = np.linspace(-1, 1, 100)
        eta = np.linspace(-1, 1, 100)
        values = basis.value(xi, eta)
        x_derivatives = basis.gradx(xi, eta)
        ```
    """

    def __init__(self, num_shape_functions: int):
        super().__init__(num_shape_functions)

    def jacobi_wrapper(self, n: int, a: int, b: int, x: np.ndarray) -> np.ndarray:
        """
        Evaluate the Jacobi polynomial of degree n with parameters a and b at the given points x.

        Args:
            n (int): Degree of the Jacobi polynomial.
            a (int): First parameter of the Jacobi polynomial.
            b (int): Second parameter of the Jacobi polynomial.
            x (np.ndarray): Points at which to evaluate the Jacobi polynomial.

        Returns:
            np.ndarray: Values of the Jacobi polynomial at the given points.
        """
        x = np.array(x, dtype=np.float64)
        return jacobi(n, a, b)(x)

    ## Helper Function
    def test_fcnx(self, n_test: int, x: np.ndarray) -> np.ndarray:
        """
        Compute the x-component of the test functions for a given number of test functions and x-coordinates.

        Args:
            n_test (int): Number of test functions.
            x (np.ndarray): x-coordinates at which to evaluate the test functions.

        Returns:
            np.ndarray: Values of the x-component of the test functions.
        """
        test_total = []
        for n in range(1, n_test + 1):
            test = self.jacobi_wrapper(n + 1, 0, 0, x) - self.jacobi_wrapper(
                n - 1, 0, 0, x
            )
            test_total.append(test)
        return np.asarray(test_total, np.float64)

    def test_fcny(self, n_test: int, y: np.ndarray) -> np.ndarray:
        """
        Compute the y-component of the test functions for a given number of test functions and y-coordinates.

        Args:
            n_test (int): Number of test functions.
            y (np.ndarray): y-coordinates at which to evaluate the test functions.

        Returns:
            np.ndarray: Values of the y-component of the test functions.
        """
        test_total = []
        for n in range(1, n_test + 1):
            test = self.jacobi_wrapper(n + 1, 0, 0, y) - self.jacobi_wrapper(
                n - 1, 0, 0, y
            )
            test_total.append(test)
        return np.asarray(test_total, np.float64)

    def dtest_fcn(self, n_test: int, x: np.ndarray) -> np.ndarray:
        """
        Compute the x-derivatives of the test functions for a given number of test functions and x-coordinates.

        Args:
            n_test (int): Number of test functions.
            x (np.ndarray): x-coordinates at which to evaluate the test functions.

        Returns:
            np.ndarray: Values of the x-derivatives of the test functions.
        """
        d1test_total = []
        d2test_total = []
        for n in range(1, n_test + 1):
            if n == 1:
                d1test = ((n + 2) / 2) * self.jacobi_wrapper(n, 1, 1, x)
                d2test = ((n + 2) * (n + 3) / (2 * 2)) * self.jacobi_wrapper(
                    n - 1, 2, 2, x
                )
                d1test_total.append(d1test)
                d2test_total.append(d2test)
            elif n == 2:
                d1test = ((n + 2) / 2) * self.jacobi_wrapper(n, 1, 1, x) - (
                    (n) / 2
                ) * self.jacobi_wrapper(n - 2, 1, 1, x)
                d2test = ((n + 2) * (n + 3) / (2 * 2)) * self.jacobi_wrapper(
                    n - 1, 2, 2, x
                )
                d1test_total.append(d1test)
                d2test_total.append(d2test)
            else:
                d1test = ((n + 2) / 2) * self.jacobi_wrapper(n, 1, 1, x) - (
                    (n) / 2
                ) * self.jacobi_wrapper(n - 2, 1, 1, x)
                d2test = ((n + 2) * (n + 3) / (2 * 2)) * self.jacobi_wrapper(
                    n - 1, 2, 2, x
                ) - ((n) * (n + 1) / (2 * 2)) * self.jacobi_wrapper(n - 3, 2, 2, x)
                d1test_total.append(d1test)
                d2test_total.append(d2test)
        return np.asarray(d1test_total), np.asarray(d2test_total)

    def value(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        This method returns the values of the basis functions at the given (xi, eta) coordinates.

        Args:
            xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
            eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

        Returns:
            np.ndarray: Values of the basis functions.
        """
        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
        test_x = self.test_fcnx(num_shape_func_in_1d, xi)
        test_y = self.test_fcny(num_shape_func_in_1d, eta)
        values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

        for i in range(num_shape_func_in_1d):
            values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
                test_x[i, :] * test_y
            )

        return values

    def gradx(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        This method returns the x-derivatives of the basis functions at the given (xi, eta) coordinates.

        Args:
            xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
            eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

        Returns:
            np.ndarray: Values of the x-derivatives of the basis functions.
        """
        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
        grad_test_x = self.dtest_fcn(num_shape_func_in_1d, xi)[0]
        test_y = self.test_fcny(num_shape_func_in_1d, eta)
        values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

        for i in range(num_shape_func_in_1d):
            values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
                grad_test_x[i, :] * test_y
            )

        return values

    def grady(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        This method returns the y-derivatives of the basis functions at the given (xi, eta) coordinates.

        Args:
            xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
            eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

        Returns:
            np.ndarray: Values of the y-derivatives of the basis functions.
        """
        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
        test_x = self.test_fcnx(num_shape_func_in_1d, xi)
        grad_test_y = self.dtest_fcn(num_shape_func_in_1d, eta)[0]
        values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

        for i in range(num_shape_func_in_1d):
            values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
                test_x[i, :] * grad_test_y
            )

        return values

    def gradxx(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        This method returns the xx-derivatives of the basis functions at the given (xi, eta) coordinates.

        Args:
            xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
            eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

        Returns:
            np.ndarray: Values of the xx-derivatives of the basis functions.
        """
        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
        grad_grad_x = self.dtest_fcn(num_shape_func_in_1d, xi)[1]
        test_y = self.test_fcny(num_shape_func_in_1d, eta)
        values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

        for i in range(num_shape_func_in_1d):
            values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
                grad_grad_x[i, :] * test_y
            )

        return values

    def gradxy(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        This method returns the xy-derivatives of the basis functions at the given (xi, eta) coordinates.

        Args:
            xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
            eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

        Returns:
            np.ndarray: Values of the xy-derivatives of the basis functions.
        """
        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
        grad_test_x = self.dtest_fcn(num_shape_func_in_1d, xi)[0]
        grad_test_y = self.dtest_fcn(num_shape_func_in_1d, eta)[0]
        values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

        for i in range(num_shape_func_in_1d):
            values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
                grad_test_x[i, :] * grad_test_y
            )

        return values

    def gradyy(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        This method returns the yy-derivatives of the basis functions at the given (xi, eta) coordinates.

        Args:
            xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
            eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

        Returns:
            np.ndarray: Values of the yy-derivatives of the basis functions.
        """
        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
        test_x = self.test_fcnx(num_shape_func_in_1d, xi)
        grad_grad_y = self.dtest_fcn(num_shape_func_in_1d, eta)[1]
        values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

        for i in range(num_shape_func_in_1d):
            values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
                test_x[i, :] * grad_grad_y
            )

        return values

dtest_fcn(n_test, x)

Compute the x-derivatives of the test functions for a given number of test functions and x-coordinates.

Parameters:

Name	Type	Description	Default
n_test	int	Number of test functions.	required
x	ndarray	x-coordinates at which to evaluate the test functions.	required

Returns:

Type	Description
ndarray	np.ndarray: Values of the x-derivatives of the test functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_legendre.py

def dtest_fcn(self, n_test: int, x: np.ndarray) -> np.ndarray:
    """
    Compute the x-derivatives of the test functions for a given number of test functions and x-coordinates.

    Args:
        n_test (int): Number of test functions.
        x (np.ndarray): x-coordinates at which to evaluate the test functions.

    Returns:
        np.ndarray: Values of the x-derivatives of the test functions.
    """
    d1test_total = []
    d2test_total = []
    for n in range(1, n_test + 1):
        if n == 1:
            d1test = ((n + 2) / 2) * self.jacobi_wrapper(n, 1, 1, x)
            d2test = ((n + 2) * (n + 3) / (2 * 2)) * self.jacobi_wrapper(
                n - 1, 2, 2, x
            )
            d1test_total.append(d1test)
            d2test_total.append(d2test)
        elif n == 2:
            d1test = ((n + 2) / 2) * self.jacobi_wrapper(n, 1, 1, x) - (
                (n) / 2
            ) * self.jacobi_wrapper(n - 2, 1, 1, x)
            d2test = ((n + 2) * (n + 3) / (2 * 2)) * self.jacobi_wrapper(
                n - 1, 2, 2, x
            )
            d1test_total.append(d1test)
            d2test_total.append(d2test)
        else:
            d1test = ((n + 2) / 2) * self.jacobi_wrapper(n, 1, 1, x) - (
                (n) / 2
            ) * self.jacobi_wrapper(n - 2, 1, 1, x)
            d2test = ((n + 2) * (n + 3) / (2 * 2)) * self.jacobi_wrapper(
                n - 1, 2, 2, x
            ) - ((n) * (n + 1) / (2 * 2)) * self.jacobi_wrapper(n - 3, 2, 2, x)
            d1test_total.append(d1test)
            d2test_total.append(d2test)
    return np.asarray(d1test_total), np.asarray(d2test_total)

gradx(xi, eta)

This method returns the x-derivatives of the basis functions at the given (xi, eta) coordinates.

Parameters:

Name	Type	Description	Default
xi	ndarray	x-coordinates at which to evaluate the basis functions.	required
eta	ndarray	y-coordinates at which to evaluate the basis functions.	required

Returns:

Type	Description
ndarray	np.ndarray: Values of the x-derivatives of the basis functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_legendre.py

def gradx(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    This method returns the x-derivatives of the basis functions at the given (xi, eta) coordinates.

    Args:
        xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
        eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

    Returns:
        np.ndarray: Values of the x-derivatives of the basis functions.
    """
    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
    grad_test_x = self.dtest_fcn(num_shape_func_in_1d, xi)[0]
    test_y = self.test_fcny(num_shape_func_in_1d, eta)
    values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

    for i in range(num_shape_func_in_1d):
        values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
            grad_test_x[i, :] * test_y
        )

    return values

gradxx(xi, eta)

This method returns the xx-derivatives of the basis functions at the given (xi, eta) coordinates.

Parameters:

Name	Type	Description	Default
xi	ndarray	x-coordinates at which to evaluate the basis functions.	required
eta	ndarray	y-coordinates at which to evaluate the basis functions.	required

Returns:

Type	Description
ndarray	np.ndarray: Values of the xx-derivatives of the basis functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_legendre.py

def gradxx(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    This method returns the xx-derivatives of the basis functions at the given (xi, eta) coordinates.

    Args:
        xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
        eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

    Returns:
        np.ndarray: Values of the xx-derivatives of the basis functions.
    """
    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
    grad_grad_x = self.dtest_fcn(num_shape_func_in_1d, xi)[1]
    test_y = self.test_fcny(num_shape_func_in_1d, eta)
    values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

    for i in range(num_shape_func_in_1d):
        values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
            grad_grad_x[i, :] * test_y
        )

    return values

gradxy(xi, eta)

This method returns the xy-derivatives of the basis functions at the given (xi, eta) coordinates.

Parameters:

Name	Type	Description	Default
xi	ndarray	x-coordinates at which to evaluate the basis functions.	required
eta	ndarray	y-coordinates at which to evaluate the basis functions.	required

Returns:

Type	Description
ndarray	np.ndarray: Values of the xy-derivatives of the basis functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_legendre.py

def gradxy(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    This method returns the xy-derivatives of the basis functions at the given (xi, eta) coordinates.

    Args:
        xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
        eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

    Returns:
        np.ndarray: Values of the xy-derivatives of the basis functions.
    """
    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
    grad_test_x = self.dtest_fcn(num_shape_func_in_1d, xi)[0]
    grad_test_y = self.dtest_fcn(num_shape_func_in_1d, eta)[0]
    values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

    for i in range(num_shape_func_in_1d):
        values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
            grad_test_x[i, :] * grad_test_y
        )

    return values

grady(xi, eta)

This method returns the y-derivatives of the basis functions at the given (xi, eta) coordinates.

Parameters:

Name	Type	Description	Default
xi	ndarray	x-coordinates at which to evaluate the basis functions.	required
eta	ndarray	y-coordinates at which to evaluate the basis functions.	required

Returns:

Type	Description
ndarray	np.ndarray: Values of the y-derivatives of the basis functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_legendre.py

def grady(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    This method returns the y-derivatives of the basis functions at the given (xi, eta) coordinates.

    Args:
        xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
        eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

    Returns:
        np.ndarray: Values of the y-derivatives of the basis functions.
    """
    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
    test_x = self.test_fcnx(num_shape_func_in_1d, xi)
    grad_test_y = self.dtest_fcn(num_shape_func_in_1d, eta)[0]
    values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

    for i in range(num_shape_func_in_1d):
        values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
            test_x[i, :] * grad_test_y
        )

    return values

gradyy(xi, eta)

This method returns the yy-derivatives of the basis functions at the given (xi, eta) coordinates.

Parameters:

Name	Type	Description	Default
xi	ndarray	x-coordinates at which to evaluate the basis functions.	required
eta	ndarray	y-coordinates at which to evaluate the basis functions.	required

Returns:

Type	Description
ndarray	np.ndarray: Values of the yy-derivatives of the basis functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_legendre.py

def gradyy(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    This method returns the yy-derivatives of the basis functions at the given (xi, eta) coordinates.

    Args:
        xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
        eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

    Returns:
        np.ndarray: Values of the yy-derivatives of the basis functions.
    """
    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
    test_x = self.test_fcnx(num_shape_func_in_1d, xi)
    grad_grad_y = self.dtest_fcn(num_shape_func_in_1d, eta)[1]
    values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

    for i in range(num_shape_func_in_1d):
        values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
            test_x[i, :] * grad_grad_y
        )

    return values

jacobi_wrapper(n, a, b, x)

Evaluate the Jacobi polynomial of degree n with parameters a and b at the given points x.

Parameters:

Name	Type	Description	Default
n	int	Degree of the Jacobi polynomial.	required
a	int	First parameter of the Jacobi polynomial.	required
b	int	Second parameter of the Jacobi polynomial.	required
x	ndarray	Points at which to evaluate the Jacobi polynomial.	required

Returns:

Type	Description
ndarray	np.ndarray: Values of the Jacobi polynomial at the given points.

Source code in scirex/core/sciml/fe/basis_2d_qn_legendre.py

def jacobi_wrapper(self, n: int, a: int, b: int, x: np.ndarray) -> np.ndarray:
    """
    Evaluate the Jacobi polynomial of degree n with parameters a and b at the given points x.

    Args:
        n (int): Degree of the Jacobi polynomial.
        a (int): First parameter of the Jacobi polynomial.
        b (int): Second parameter of the Jacobi polynomial.
        x (np.ndarray): Points at which to evaluate the Jacobi polynomial.

    Returns:
        np.ndarray: Values of the Jacobi polynomial at the given points.
    """
    x = np.array(x, dtype=np.float64)
    return jacobi(n, a, b)(x)

test_fcnx(n_test, x)

Compute the x-component of the test functions for a given number of test functions and x-coordinates.

Parameters:

Name	Type	Description	Default
n_test	int	Number of test functions.	required
x	ndarray	x-coordinates at which to evaluate the test functions.	required

Returns:

Type	Description
ndarray	np.ndarray: Values of the x-component of the test functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_legendre.py

def test_fcnx(self, n_test: int, x: np.ndarray) -> np.ndarray:
    """
    Compute the x-component of the test functions for a given number of test functions and x-coordinates.

    Args:
        n_test (int): Number of test functions.
        x (np.ndarray): x-coordinates at which to evaluate the test functions.

    Returns:
        np.ndarray: Values of the x-component of the test functions.
    """
    test_total = []
    for n in range(1, n_test + 1):
        test = self.jacobi_wrapper(n + 1, 0, 0, x) - self.jacobi_wrapper(
            n - 1, 0, 0, x
        )
        test_total.append(test)
    return np.asarray(test_total, np.float64)

test_fcny(n_test, y)

Compute the y-component of the test functions for a given number of test functions and y-coordinates.

Parameters:

Name	Type	Description	Default
n_test	int	Number of test functions.	required
y	ndarray	y-coordinates at which to evaluate the test functions.	required

Returns:

Type	Description
ndarray	np.ndarray: Values of the y-component of the test functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_legendre.py

def test_fcny(self, n_test: int, y: np.ndarray) -> np.ndarray:
    """
    Compute the y-component of the test functions for a given number of test functions and y-coordinates.

    Args:
        n_test (int): Number of test functions.
        y (np.ndarray): y-coordinates at which to evaluate the test functions.

    Returns:
        np.ndarray: Values of the y-component of the test functions.
    """
    test_total = []
    for n in range(1, n_test + 1):
        test = self.jacobi_wrapper(n + 1, 0, 0, y) - self.jacobi_wrapper(
            n - 1, 0, 0, y
        )
        test_total.append(test)
    return np.asarray(test_total, np.float64)

value(xi, eta)

This method returns the values of the basis functions at the given (xi, eta) coordinates.

Parameters:

Name	Type	Description	Default
xi	ndarray	x-coordinates at which to evaluate the basis functions.	required
eta	ndarray	y-coordinates at which to evaluate the basis functions.	required

Returns:

Type	Description
ndarray	np.ndarray: Values of the basis functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_legendre.py

def value(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    This method returns the values of the basis functions at the given (xi, eta) coordinates.

    Args:
        xi (np.ndarray): x-coordinates at which to evaluate the basis functions.
        eta (np.ndarray): y-coordinates at which to evaluate the basis functions.

    Returns:
        np.ndarray: Values of the basis functions.
    """
    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))
    test_x = self.test_fcnx(num_shape_func_in_1d, xi)
    test_y = self.test_fcny(num_shape_func_in_1d, eta)
    values = np.zeros((self.num_shape_functions, len(xi)), dtype=np.float64)

    for i in range(num_shape_func_in_1d):
        values[num_shape_func_in_1d * i : num_shape_func_in_1d * (i + 1), :] = (
            test_x[i, :] * test_y
        )

    return values

Basis2DQNLegendreSpecial

Bases: BasisFunction2D

A specialized implementation of two-dimensional basis functions using Legendre polynomials for Q1 elements.

This class provides a complete implementation for computing basis functions and their derivatives in two dimensions. The basis functions are constructed using a special formulation based on differences of consecutive Legendre polynomials.

The class inherits from BasisFunction2D and implements all required methods for computing function values and derivatives up to second order.

Attributes:

Name	Type	Description
num_shape_functions	int	Total number of shape functions in the 2D element. Must be a perfect square as it represents tensor product of 1D functions.

Methods:

Name	Description
test_fcn	Computes test functions using Legendre polynomial differences
test_grad_fcn	Computes first derivatives of test functions
test_grad_grad_fcn	Computes second derivatives of test functions
value	Computes values of all basis functions
gradx	Computes x-derivatives of all basis functions
grady	Computes y-derivatives of all basis functions
gradxx	Computes second x-derivatives of all basis functions
gradyy	Computes second y-derivatives of all basis functions
gradxy	Computes mixed xy-derivatives of all basis functions

Implementation Details

• Basis functions are constructed using differences of consecutive Legendre polynomials
• Test functions are created using Pn+1(x) - Pn-1(x) where Pn is the nth Legendre polynomial
• All computations maintain numerical precision using numpy arrays
• Efficient vectorized operations for multiple point evaluations
• Tensor product construction for 2D basis functions

Example

basis = Basis2DQNLegendreSpecial(num_shape_functions=16)  # Creates 4x4 basis functions
xi = np.linspace(-1, 1, 100)
eta = np.linspace(-1, 1, 100)
values = basis.value(xi, eta)
x_derivatives = basis.gradx(xi, eta)

Source code in scirex/core/sciml/fe/basis_2d_qn_legendre_special.py

class Basis2DQNLegendreSpecial(BasisFunction2D):
    """
    A specialized implementation of two-dimensional basis functions using Legendre polynomials for Q1 elements.

    This class provides a complete implementation for computing basis functions and their derivatives
    in two dimensions. The basis functions are constructed using a special formulation based on
    differences of consecutive Legendre polynomials.

    The class inherits from BasisFunction2D and implements all required methods for computing
    function values and derivatives up to second order.

    Attributes:
        num_shape_functions (int): Total number of shape functions in the 2D element.
            Must be a perfect square as it represents tensor product of 1D functions.

    Methods:
        test_fcn(n_test, x): Computes test functions using Legendre polynomial differences
        test_grad_fcn(n_test, x): Computes first derivatives of test functions
        test_grad_grad_fcn(n_test, x): Computes second derivatives of test functions
        value(xi, eta): Computes values of all basis functions
        gradx(xi, eta): Computes x-derivatives of all basis functions
        grady(xi, eta): Computes y-derivatives of all basis functions
        gradxx(xi, eta): Computes second x-derivatives of all basis functions
        gradyy(xi, eta): Computes second y-derivatives of all basis functions
        gradxy(xi, eta): Computes mixed xy-derivatives of all basis functions

    Implementation Details:
        - Basis functions are constructed using differences of consecutive Legendre polynomials
        - Test functions are created using Pn+1(x) - Pn-1(x) where Pn is the nth Legendre polynomial
        - All computations maintain numerical precision using numpy arrays
        - Efficient vectorized operations for multiple point evaluations
        - Tensor product construction for 2D basis functions

    Example:
        ```python
        basis = Basis2DQNLegendreSpecial(num_shape_functions=16)  # Creates 4x4 basis functions
        xi = np.linspace(-1, 1, 100)
        eta = np.linspace(-1, 1, 100)
        values = basis.value(xi, eta)
        x_derivatives = basis.gradx(xi, eta)
        ```
    """

    def __init__(self, num_shape_functions: int):
        super().__init__(num_shape_functions)

    def test_fcn(self, n_test: int, x: np.ndarray) -> np.ndarray:
        """
        Calculate the test function values for a given number of tests and input values.

        Args:
            n_test (int): The number of test functions to calculate.
            x (np.ndarray): The input values at which to evaluate the test functions.

        Returns:
            np.ndarray: An array containing the results of the test functions at the given input values.
        """
        test_total = []
        for n in range(1, n_test + 1):
            obj1 = legendre(n + 1)
            obj2 = legendre(n - 1)
            test = obj1(x) - obj2(x)
            test_total.append(test)
        return np.asarray(test_total)

    def test_grad_fcn(self, n_test: int, x: np.ndarray) -> np.ndarray:
        """
        Calculate the gradient of the test function at a given point.

        Args:
            n_test (int): The number of test cases to evaluate.
            x (np.ndarray): The input value at which to evaluate the function.

        Returns:
            np.ndarray: An array containing the results of the test cases.
        """
        test_total = []
        for n in range(1, n_test + 1):
            obj1 = legendre(n + 1).deriv()
            obj2 = legendre(n - 1).deriv()
            test = obj1(x) - obj2(x)
            test_total.append(test)
        return np.asarray(test_total)

    def test_grad_grad_fcn(self, n_test: int, x: np.ndarray) -> np.ndarray:
        """
        Calculate the gradient of the second derivative of a function using Legendre polynomials.

        Args:
            n_test (int): The number of test cases to evaluate.
            x (np.ndarray): The input value at which to evaluate the function.

        Returns:
            np.ndarray: An array containing the results of the test cases.
        """
        test_total = []
        for n in range(1, n_test + 1):
            obj1 = legendre(n + 1).deriv(2)
            obj2 = legendre(n - 1).deriv(2)
            test = obj1(x) - obj2(x)

            test_total.append(test)
        return np.asarray(test_total)

    def value(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        This method returns the values of the basis functions at the given (xi, eta) coordinates.

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

        Returns:
            np.ndarray: The values of the basis functions.
        """
        values = np.zeros((self.num_shape_functions, len(xi)))

        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))

        test_function_x = self.test_fcn(num_shape_func_in_1d, xi)
        test_function_y = self.test_fcn(num_shape_func_in_1d, eta)

        # Generate an outer product of the test functions to generate the basis functions
        for i in range(num_shape_func_in_1d):
            values[i * num_shape_func_in_1d : (i + 1) * num_shape_func_in_1d, :] = (
                test_function_x[i, :] * test_function_y
            )

        return values

    def gradx(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        This method returns the x-derivatives of the basis functions at the given (xi, eta) coordinates.

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

        Returns:
            np.ndarray: The x-derivatives of the basis functions.
        """
        values = np.zeros((self.num_shape_functions, len(xi)))

        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))

        test_function_grad_x = self.test_grad_fcn(num_shape_func_in_1d, xi)
        test_function_y = self.test_fcn(num_shape_func_in_1d, eta)

        # Generate an outer product of the test functions to generate the basis functions
        for i in range(num_shape_func_in_1d):
            values[i * num_shape_func_in_1d : (i + 1) * num_shape_func_in_1d, :] = (
                test_function_grad_x[i, :] * test_function_y
            )

        return values

    def grady(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        This method returns the y-derivatives of the basis functions at the given (xi, eta) coordinates.

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

        Returns:
            np.ndarray: The y-derivatives of the basis functions.
        """
        values = np.zeros((self.num_shape_functions, len(xi)))

        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))

        test_function_x = self.test_fcn(num_shape_func_in_1d, xi)
        test_function_grad_y = self.test_grad_fcn(num_shape_func_in_1d, eta)

        # Generate an outer product of the test functions to generate the basis functions
        for i in range(num_shape_func_in_1d):
            values[i * num_shape_func_in_1d : (i + 1) * num_shape_func_in_1d, :] = (
                test_function_x[i, :] * test_function_grad_y
            )

        return values

    def gradxx(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        This method returns the xx-derivatives of the basis functions at the given (xi, eta) coordinates.

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

        Returns:
            np.ndarray: The xx-derivatives of the basis functions.
        """
        values = np.zeros((self.num_shape_functions, len(xi)))

        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))

        test_function_grad_grad_x = self.test_grad_grad_fcn(num_shape_func_in_1d, xi)
        test_function_y = self.test_fcn(num_shape_func_in_1d, eta)

        # Generate an outer product of the test functions to generate the basis functions
        for i in range(num_shape_func_in_1d):
            values[i * num_shape_func_in_1d : (i + 1) * num_shape_func_in_1d, :] = (
                test_function_grad_grad_x[i, :] * test_function_y
            )

        return values

    def gradxy(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        This method returns the xy-derivatives of the basis functions at the given (xi, eta) coordinates.

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

        Returns:
            np.ndarray: The xy-derivatives of the basis functions.
        """
        values = np.zeros((self.num_shape_functions, len(xi)))

        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))

        test_function_grad_x = self.test_grad_fcn(num_shape_func_in_1d, xi)
        test_function_grad_y = self.test_grad_fcn(num_shape_func_in_1d, eta)

        # Generate an outer product of the test functions to generate the basis functions
        for i in range(num_shape_func_in_1d):
            values[i * num_shape_func_in_1d : (i + 1) * num_shape_func_in_1d, :] = (
                test_function_grad_x[i, :] * test_function_grad_y
            )

        return values

    def gradyy(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        This method returns the yy-derivatives of the basis functions at the given (xi, eta) coordinates.

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

        Returns:
            np.ndarray: The yy-derivatives of the basis functions.
        """
        values = np.zeros((self.num_shape_functions, len(xi)))

        num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))

        test_function_x = self.test_fcn(num_shape_func_in_1d, xi)
        test_function_grad_grad_y = self.test_grad_grad_fcn(num_shape_func_in_1d, eta)

        # Generate an outer product of the test functions to generate the basis functions
        for i in range(num_shape_func_in_1d):
            values[i * num_shape_func_in_1d : (i + 1) * num_shape_func_in_1d, :] = (
                test_function_x[i, :] * test_function_grad_grad_y
            )

        return values

gradx(xi, eta)

This method returns the x-derivatives of the basis functions at the given (xi, eta) coordinates.

Parameters:

Name	Type	Description	Default
xi	ndarray	The xi coordinates.	required
eta	ndarray	The eta coordinates.	required

Returns:

Type	Description
ndarray	np.ndarray: The x-derivatives of the basis functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_legendre_special.py

def gradx(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    This method returns the x-derivatives of the basis functions at the given (xi, eta) coordinates.

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

    Returns:
        np.ndarray: The x-derivatives of the basis functions.
    """
    values = np.zeros((self.num_shape_functions, len(xi)))

    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))

    test_function_grad_x = self.test_grad_fcn(num_shape_func_in_1d, xi)
    test_function_y = self.test_fcn(num_shape_func_in_1d, eta)

    # Generate an outer product of the test functions to generate the basis functions
    for i in range(num_shape_func_in_1d):
        values[i * num_shape_func_in_1d : (i + 1) * num_shape_func_in_1d, :] = (
            test_function_grad_x[i, :] * test_function_y
        )

    return values

gradxx(xi, eta)

This method returns the xx-derivatives of the basis functions at the given (xi, eta) 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 xx-derivatives of the basis functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_legendre_special.py

def gradxx(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    This method returns the xx-derivatives of the basis functions at the given (xi, eta) coordinates.

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

    Returns:
        np.ndarray: The xx-derivatives of the basis functions.
    """
    values = np.zeros((self.num_shape_functions, len(xi)))

    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))

    test_function_grad_grad_x = self.test_grad_grad_fcn(num_shape_func_in_1d, xi)
    test_function_y = self.test_fcn(num_shape_func_in_1d, eta)

    # Generate an outer product of the test functions to generate the basis functions
    for i in range(num_shape_func_in_1d):
        values[i * num_shape_func_in_1d : (i + 1) * num_shape_func_in_1d, :] = (
            test_function_grad_grad_x[i, :] * test_function_y
        )

    return values

gradxy(xi, eta)

This method returns the xy-derivatives of the basis functions at the given (xi, eta) 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 xy-derivatives of the basis functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_legendre_special.py

def gradxy(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    This method returns the xy-derivatives of the basis functions at the given (xi, eta) coordinates.

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

    Returns:
        np.ndarray: The xy-derivatives of the basis functions.
    """
    values = np.zeros((self.num_shape_functions, len(xi)))

    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))

    test_function_grad_x = self.test_grad_fcn(num_shape_func_in_1d, xi)
    test_function_grad_y = self.test_grad_fcn(num_shape_func_in_1d, eta)

    # Generate an outer product of the test functions to generate the basis functions
    for i in range(num_shape_func_in_1d):
        values[i * num_shape_func_in_1d : (i + 1) * num_shape_func_in_1d, :] = (
            test_function_grad_x[i, :] * test_function_grad_y
        )

    return values

grady(xi, eta)

This method returns the y-derivatives of the basis functions at the given (xi, eta) coordinates.

Parameters:

Name	Type	Description	Default
xi	ndarray	The xi coordinates.	required
eta	ndarray	The eta coordinates.	required

Returns:

Type	Description
ndarray	np.ndarray: The y-derivatives of the basis functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_legendre_special.py

def grady(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    This method returns the y-derivatives of the basis functions at the given (xi, eta) coordinates.

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

    Returns:
        np.ndarray: The y-derivatives of the basis functions.
    """
    values = np.zeros((self.num_shape_functions, len(xi)))

    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))

    test_function_x = self.test_fcn(num_shape_func_in_1d, xi)
    test_function_grad_y = self.test_grad_fcn(num_shape_func_in_1d, eta)

    # Generate an outer product of the test functions to generate the basis functions
    for i in range(num_shape_func_in_1d):
        values[i * num_shape_func_in_1d : (i + 1) * num_shape_func_in_1d, :] = (
            test_function_x[i, :] * test_function_grad_y
        )

    return values

gradyy(xi, eta)

This method returns the yy-derivatives of the basis functions at the given (xi, eta) coordinates.

Parameters:

Name	Type	Description	Default
xi	ndarray	The xi coordinates.	required
eta	ndarray	The eta coordinates.	required

Returns:

Type	Description
ndarray	np.ndarray: The yy-derivatives of the basis functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_legendre_special.py

def gradyy(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    This method returns the yy-derivatives of the basis functions at the given (xi, eta) coordinates.

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

    Returns:
        np.ndarray: The yy-derivatives of the basis functions.
    """
    values = np.zeros((self.num_shape_functions, len(xi)))

    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))

    test_function_x = self.test_fcn(num_shape_func_in_1d, xi)
    test_function_grad_grad_y = self.test_grad_grad_fcn(num_shape_func_in_1d, eta)

    # Generate an outer product of the test functions to generate the basis functions
    for i in range(num_shape_func_in_1d):
        values[i * num_shape_func_in_1d : (i + 1) * num_shape_func_in_1d, :] = (
            test_function_x[i, :] * test_function_grad_grad_y
        )

    return values

test_fcn(n_test, x)

Calculate the test function values for a given number of tests and input values.

Parameters:

Name	Type	Description	Default
n_test	int	The number of test functions to calculate.	required
x	ndarray	The input values at which to evaluate the test functions.	required

Returns:

Type	Description
ndarray	np.ndarray: An array containing the results of the test functions at the given input values.

Source code in scirex/core/sciml/fe/basis_2d_qn_legendre_special.py

def test_fcn(self, n_test: int, x: np.ndarray) -> np.ndarray:
    """
    Calculate the test function values for a given number of tests and input values.

    Args:
        n_test (int): The number of test functions to calculate.
        x (np.ndarray): The input values at which to evaluate the test functions.

    Returns:
        np.ndarray: An array containing the results of the test functions at the given input values.
    """
    test_total = []
    for n in range(1, n_test + 1):
        obj1 = legendre(n + 1)
        obj2 = legendre(n - 1)
        test = obj1(x) - obj2(x)
        test_total.append(test)
    return np.asarray(test_total)

test_grad_fcn(n_test, x)

Calculate the gradient of the test function at a given point.

Parameters:

Name	Type	Description	Default
n_test	int	The number of test cases to evaluate.	required
x	ndarray	The input value at which to evaluate the function.	required

Returns:

Type	Description
ndarray	np.ndarray: An array containing the results of the test cases.

Source code in scirex/core/sciml/fe/basis_2d_qn_legendre_special.py

def test_grad_fcn(self, n_test: int, x: np.ndarray) -> np.ndarray:
    """
    Calculate the gradient of the test function at a given point.

    Args:
        n_test (int): The number of test cases to evaluate.
        x (np.ndarray): The input value at which to evaluate the function.

    Returns:
        np.ndarray: An array containing the results of the test cases.
    """
    test_total = []
    for n in range(1, n_test + 1):
        obj1 = legendre(n + 1).deriv()
        obj2 = legendre(n - 1).deriv()
        test = obj1(x) - obj2(x)
        test_total.append(test)
    return np.asarray(test_total)

test_grad_grad_fcn(n_test, x)

Calculate the gradient of the second derivative of a function using Legendre polynomials.

Parameters:

Name	Type	Description	Default
n_test	int	The number of test cases to evaluate.	required
x	ndarray	The input value at which to evaluate the function.	required

Returns:

Type	Description
ndarray	np.ndarray: An array containing the results of the test cases.

Source code in scirex/core/sciml/fe/basis_2d_qn_legendre_special.py

def test_grad_grad_fcn(self, n_test: int, x: np.ndarray) -> np.ndarray:
    """
    Calculate the gradient of the second derivative of a function using Legendre polynomials.

    Args:
        n_test (int): The number of test cases to evaluate.
        x (np.ndarray): The input value at which to evaluate the function.

    Returns:
        np.ndarray: An array containing the results of the test cases.
    """
    test_total = []
    for n in range(1, n_test + 1):
        obj1 = legendre(n + 1).deriv(2)
        obj2 = legendre(n - 1).deriv(2)
        test = obj1(x) - obj2(x)

        test_total.append(test)
    return np.asarray(test_total)

value(xi, eta)

This method returns the values of the basis functions at the given (xi, eta) coordinates.

Parameters:

Name	Type	Description	Default
xi	ndarray	The xi coordinates.	required
eta	ndarray	The eta coordinates.	required

Returns:

Type	Description
ndarray	np.ndarray: The values of the basis functions.

Source code in scirex/core/sciml/fe/basis_2d_qn_legendre_special.py

def value(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    This method returns the values of the basis functions at the given (xi, eta) coordinates.

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

    Returns:
        np.ndarray: The values of the basis functions.
    """
    values = np.zeros((self.num_shape_functions, len(xi)))

    num_shape_func_in_1d = int(np.sqrt(self.num_shape_functions))

    test_function_x = self.test_fcn(num_shape_func_in_1d, xi)
    test_function_y = self.test_fcn(num_shape_func_in_1d, eta)

    # Generate an outer product of the test functions to generate the basis functions
    for i in range(num_shape_func_in_1d):
        values[i * num_shape_func_in_1d : (i + 1) * num_shape_func_in_1d, :] = (
            test_function_x[i, :] * test_function_y
        )

    return values

BasisFunction2D

An abstract base class defining the interface for two-dimensional finite element basis functions.

This class serves as a template for implementing various types of 2D basis functions (Legendre, Jacobi, Chebyshev, etc.) used in finite element computations. It defines the required methods for function evaluation and derivatives.

Attributes:

Name	Type	Description
num_shape_functions	int	Number of shape functions in the element. Typically a perfect square for tensor-product bases.

Methods:

Name	Description
value	Evaluates basis functions at given reference coordinates Args: xi (float): First reference coordinate eta (float): Second reference coordinate Returns: float: Values of basis functions at (xi, eta)
gradx	Computes x-derivatives at reference coordinates Args: xi (float): First reference coordinate eta (float): Second reference coordinate Returns: float: Values of x-derivatives at (xi, eta)
grady	Computes y-derivatives at reference coordinates Args: xi (float): First reference coordinate eta (float): Second reference coordinate Returns: float: Values of y-derivatives at (xi, eta)
gradxx	Computes second x-derivatives at reference coordinates Args: xi (float): First reference coordinate eta (float): Second reference coordinate Returns: float: Values of second x-derivatives at (xi, eta)
gradxy	Computes mixed derivatives at reference coordinates Args: xi (float): First reference coordinate eta (float): Second reference coordinate Returns: float: Values of mixed derivatives at (xi, eta)
gradyy	Computes second y-derivatives at reference coordinates Args: xi (float): First reference coordinate eta (float): Second reference coordinate Returns: float: Values of second y-derivatives at (xi, eta)

Notes

• All coordinate inputs (xi, eta) should be in the reference element range
• Subclasses must implement all abstract methods
• Used as base class for specific polynomial implementations:
  • Legendre polynomials (normal and special variants)
  • Jacobi polynomials
  • Chebyshev polynomials

Source code in scirex/core/sciml/fe/basis_function_2d.py

class BasisFunction2D:
    """
    An abstract base class defining the interface for two-dimensional finite element basis functions.

    This class serves as a template for implementing various types of 2D basis functions
    (Legendre, Jacobi, Chebyshev, etc.) used in finite element computations. It defines
    the required methods for function evaluation and derivatives.

    Attributes:
        num_shape_functions (int): Number of shape functions in the element.
            Typically a perfect square for tensor-product bases.

    Methods:
        value(xi, eta): Evaluates basis functions at given reference coordinates
            Args:
                xi (float): First reference coordinate
                eta (float): Second reference coordinate
            Returns:
                float: Values of basis functions at (xi, eta)

        gradx(xi, eta): Computes x-derivatives at reference coordinates
            Args:
                xi (float): First reference coordinate
                eta (float): Second reference coordinate
            Returns:
                float: Values of x-derivatives at (xi, eta)

        grady(xi, eta): Computes y-derivatives at reference coordinates
            Args:
                xi (float): First reference coordinate
                eta (float): Second reference coordinate
            Returns:
                float: Values of y-derivatives at (xi, eta)

        gradxx(xi, eta): Computes second x-derivatives at reference coordinates
            Args:
                xi (float): First reference coordinate
                eta (float): Second reference coordinate
            Returns:
                float: Values of second x-derivatives at (xi, eta)

        gradxy(xi, eta): Computes mixed derivatives at reference coordinates
            Args:
                xi (float): First reference coordinate
                eta (float): Second reference coordinate
            Returns:
                float: Values of mixed derivatives at (xi, eta)

        gradyy(xi, eta): Computes second y-derivatives at reference coordinates
            Args:
                xi (float): First reference coordinate
                eta (float): Second reference coordinate
            Returns:
                float: Values of second y-derivatives at (xi, eta)

    Notes:
        - All coordinate inputs (xi, eta) should be in the reference element range
        - Subclasses must implement all abstract methods
        - Used as base class for specific polynomial implementations:
            - Legendre polynomials (normal and special variants)
            - Jacobi polynomials
            - Chebyshev polynomials
    """

    def __init__(self, num_shape_functions):
        self.num_shape_functions = num_shape_functions

    @abstractmethod
    def value(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        Evaluates the basis function at the given xi and eta coordinates.

        Args:
            xi (float): The xi coordinate.
            eta (float): The eta coordinate.

        Returns:
            float: The value of the basis function at ( xi, eta).
        """
        pass

    @abstractmethod
    def gradx(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        Computes the partial derivative of the basis function with respect to xi.

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

        Returns:
            np.ndarray: The partial derivative of the basis function with respect to xi.
        """
        pass

    @abstractmethod
    def grady(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        Computes the partial derivative of the basis function with respect to eta.

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

        Returns:
            np.ndarray: The partial derivative of the basis function with respect to eta.
        """
        pass

    @abstractmethod
    def gradxx(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        Computes the second partial derivative of the basis function with respect to xi.

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

        Returns:
            np.ndarray: The second partial derivative of the basis function with respect to xi.
        """
        pass

    @abstractmethod
    def gradxy(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        Computes the mixed partial derivative of the basis function with respect to xi and eta.

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

        Returns:
            np.ndarray: The mixed partial derivative of the basis function with respect to xi and eta.
        """
        pass

    @abstractmethod
    def gradyy(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        Computes the second partial derivative of the basis function with respect to eta.

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

        Returns:
            np.ndarray: The second partial derivative of the basis function with respect to eta.
        """
        pass

gradx(xi, eta) abstractmethod

Computes the partial derivative of the basis function with respect to xi.

Parameters:

Name	Type	Description	Default
xi	ndarray	The xi coordinate.	required
eta	ndarray	The eta coordinate.	required

Returns:

Type	Description
ndarray	np.ndarray: The partial derivative of the basis function with respect to xi.

Source code in scirex/core/sciml/fe/basis_function_2d.py

@abstractmethod
def gradx(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    Computes the partial derivative of the basis function with respect to xi.

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

    Returns:
        np.ndarray: The partial derivative of the basis function with respect to xi.
    """
    pass

gradxx(xi, eta) abstractmethod

Computes the second partial derivative of the basis function with respect to xi.

Parameters:

Name	Type	Description	Default
xi	ndarray	The xi coordinate.	required
eta	ndarray	The eta coordinate.	required

Returns:

Type	Description
ndarray	np.ndarray: The second partial derivative of the basis function with respect to xi.

Source code in scirex/core/sciml/fe/basis_function_2d.py

@abstractmethod
def gradxx(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    Computes the second partial derivative of the basis function with respect to xi.

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

    Returns:
        np.ndarray: The second partial derivative of the basis function with respect to xi.
    """
    pass

gradxy(xi, eta) abstractmethod

Computes the mixed partial derivative of the basis function with respect to xi and eta.

Parameters:

Name	Type	Description	Default
xi	ndarray	The xi coordinate.	required
eta	ndarray	The eta coordinate.	required

Returns:

Type	Description
ndarray	np.ndarray: The mixed partial derivative of the basis function with respect to xi and eta.

Source code in scirex/core/sciml/fe/basis_function_2d.py

@abstractmethod
def gradxy(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    Computes the mixed partial derivative of the basis function with respect to xi and eta.

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

    Returns:
        np.ndarray: The mixed partial derivative of the basis function with respect to xi and eta.
    """
    pass

grady(xi, eta) abstractmethod

Computes the partial derivative of the basis function with respect to eta.

Parameters:

Name	Type	Description	Default
xi	ndarray	The xi coordinate.	required
eta	ndarray	The eta coordinate.	required

Returns:

Type	Description
ndarray	np.ndarray: The partial derivative of the basis function with respect to eta.

Source code in scirex/core/sciml/fe/basis_function_2d.py

@abstractmethod
def grady(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    Computes the partial derivative of the basis function with respect to eta.

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

    Returns:
        np.ndarray: The partial derivative of the basis function with respect to eta.
    """
    pass

gradyy(xi, eta) abstractmethod

Computes the second partial derivative of the basis function with respect to eta.

Parameters:

Name	Type	Description	Default
xi	ndarray	The xi coordinate.	required
eta	ndarray	The eta coordinate.	required

Returns:

Type	Description
ndarray	np.ndarray: The second partial derivative of the basis function with respect to eta.

Source code in scirex/core/sciml/fe/basis_function_2d.py

@abstractmethod
def gradyy(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    Computes the second partial derivative of the basis function with respect to eta.

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

    Returns:
        np.ndarray: The second partial derivative of the basis function with respect to eta.
    """
    pass

value(xi, eta) abstractmethod

Evaluates the basis function at the given xi and eta coordinates.

Parameters:

Name	Type	Description	Default
xi	float	The xi coordinate.	required
eta	float	The eta coordinate.	required

Returns:

Name	Type	Description
float	ndarray	The value of the basis function at ( xi, eta).

Source code in scirex/core/sciml/fe/basis_function_2d.py

@abstractmethod
def value(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    Evaluates the basis function at the given xi and eta coordinates.

    Args:
        xi (float): The xi coordinate.
        eta (float): The eta coordinate.

    Returns:
        float: The value of the basis function at ( xi, eta).
    """
    pass

FE2DSetupMain

Main configuration class for 2D finite element analysis setup.

This class handles the configuration and initialization of finite element analysis components, including basis functions, quadrature rules, and geometric transformations.

Attributes:

Name	Type	Description
cell_type	str	Type of finite element ('quadrilateral')
fe_order	int	Order of finite element approximation (1 < order < 1e3)
fe_type	str	Type of basis functions ('legendre', 'legendre_special', 'chebyshev_2', 'jacobi_plain')
quad_order	int	Order of quadrature rule (>= 2)
quad_type	str	Type of quadrature formula
n_nodes	int	Number of nodes in the element

Example

setup = FE2DSetupMain( ... cell_type='quadrilateral', ... fe_order=2, ... fe_type='legendre', ... quad_order=3, ... quad_type='gauss' ... ) basis = setup.assign_basis_function() weights, xi, eta = setup.assign_quadrature_rules()

Notes

• Supports only quadrilateral elements currently
• Validates all input parameters for correctness
• Provides different polynomial basis options
• Handles both affine and bilinear transformations
• Quadrature order must be >= 3 for accuracy

Source code in scirex/core/sciml/fe/fe2d_setup_main.py

class FE2DSetupMain:
    """
    Main configuration class for 2D finite element analysis setup.

    This class handles the configuration and initialization of finite element
    analysis components, including basis functions, quadrature rules, and
    geometric transformations.

    Attributes:
        cell_type (str): Type of finite element ('quadrilateral')
        fe_order (int): Order of finite element approximation (1 < order < 1e3)
        fe_type (str): Type of basis functions
            ('legendre', 'legendre_special', 'chebyshev_2', 'jacobi_plain')
        quad_order (int): Order of quadrature rule (>= 2)
        quad_type (str): Type of quadrature formula
        n_nodes (int): Number of nodes in the element

    Example:
        >>> setup = FE2DSetupMain(
        ...     cell_type='quadrilateral',
        ...     fe_order=2,
        ...     fe_type='legendre',
        ...     quad_order=3,
        ...     quad_type='gauss'
        ... )
        >>> basis = setup.assign_basis_function()
        >>> weights, xi, eta = setup.assign_quadrature_rules()

    Notes:
        - Supports only quadrilateral elements currently
        - Validates all input parameters for correctness
        - Provides different polynomial basis options
        - Handles both affine and bilinear transformations
        - Quadrature order must be >= 3 for accuracy
    """

    def __init__(
        self,
        cell_type: str,
        fe_order: int,
        fe_type: str,
        quad_order: int,
        quad_type: str,
    ):
        """
        Constructor for the FE2DSetupMain class.

        Args:
            cell_type (str): Type of finite element ('quadrilateral')
            fe_order (int): Order of finite element approximation (1 < order < 1e3)
            fe_type (str): Type of basis functions
                ('legendre', 'legendre_special', 'chebyshev_2', 'jacobi_plain')
            quad_order (int): Order of quadrature rule (>= 2)
            quad_type (str): Type of quadrature formula

        Raises:
            None

        Returns:
            None
        """
        self.cell_type = cell_type
        self.fe_order = fe_order
        self.fe_type = fe_type
        self.quad_order = quad_order
        self.quad_type = quad_type

        self.assign_basis_function()

    def assign_basis_function(self) -> BasisFunction2D:
        """
        Assigns the basis function based on the cell type and the fe_order.

        Args:
            None

        Returns:
            BasisFunction2D: The basis function object for the given configuration.

        Raises:
            ValueError: If the fe order is invalid or the cell type is invalid.
        """
        # check for fe order lower bound and higher bound
        if self.fe_order <= 1 or self.fe_order >= 1e3:
            print(
                f"Invalid fe order {self.fe_order} in {self.__class__.__name__} from {__name__}."
            )
            raise ValueError("fe order should be greater than 1 and less than 1e4.")

        if self.cell_type == "quadrilateral":
            self.n_nodes = 4

            # --- LEGENDRE --- #
            if self.fe_type == "legendre" or self.fe_type == "jacobi":
                # jacobi is added for backward compatibility with prev pushes
                # generally, jacobi is referred to as Legendre basis on previous iterations
                return Basis2DQNLegendre(self.fe_order**2)

            elif self.fe_type == "legendre_special":
                return Basis2DQNLegendreSpecial(self.fe_order**2)

            # ----- CHEBYSHEV ---- #
            elif self.fe_type == "chebyshev_2":
                return Basis2DQNChebyshev2(self.fe_order**2)

            # ----- PLain jacobi ---- #
            elif self.fe_type == "jacobi_plain":
                return Basis2DQNJacobi(self.fe_order**2)

            else:
                print(
                    f"Invalid fe order {self.fe_order} in {self.__class__.__name__} from {__name__}."
                )
                raise ValueError(
                    'fe order should be one of the : "legendre" , "jacobi", "legendre_special", "chebyshev_2", "jacobi_plain"'
                )

        print(
            f"Invalid cell type {self.cell_type} in {self.__class__.__name__} from {__name__}."
        )

    def assign_quadrature_rules(self):
        """
        Assigns the quadrature rule based on the quad_order.

        Args:
            None

        Returns:
            tuple: The quadrature weights, xi and eta values in a numpy array format.

        Raises:
            ValueError: If the quad_order is invalid
            ValueError: If the cell type is invalid
            ValueError: If the quad_order is not between 1 and 9999
        """
        if self.cell_type == "quadrilateral":
            if self.quad_order < 3:
                raise ValueError("Quad order should be greater than 2.")
            elif self.quad_order >= 2 and self.quad_order <= 9999:
                weights, xi, eta = Quadratureformulas_Quad2D(
                    self.quad_order, self.quad_type
                ).get_quad_values()
                return weights, xi, eta
            else:
                print(
                    f"Invalid quad order {self.quad_order} in {self.__class__.__name__} from {__name__}."
                )
                raise ValueError("Quad order should be between 1 and 9999.")

        raise ValueError(
            f"Invalid cell type {self.cell_type} in {self.__class__.__name__} from {__name__}."
        )

    def assign_fe_transformation(
        self, fe_transformation_type: str, cell_coordinates: np.ndarray
    ) -> FETransforamtion2D:
        """
        Assigns the fe transformation based on the cell type.

        Args:
            fe_transformation_type (str): Type of fe transformation ('affine', 'bilinear')
            cell_coordinates (np.ndarray): The cell coordinates

        Returns:
            FETransforamtion2D: The fe transformation object for the given configuration.

        Raises:
            ValueError: If the cell type is invalid
            ValueError: If the fe transformation type is invalid
        """
        if self.cell_type == "quadrilateral":
            if fe_transformation_type == "affine":
                return QuadAffin(cell_coordinates)
            elif fe_transformation_type == "bilinear":
                return QuadBilinear(cell_coordinates)
            else:
                raise ValueError(
                    f"Invalid fe transformation type {fe_transformation_type} in {self.__class__.__name__} from {__name__}."
                )

        else:
            raise ValueError(
                f"Invalid cell type {self.cell_type} in {self.__class__.__name__} from {__name__}."
            )

__init__(cell_type, fe_order, fe_type, quad_order, quad_type)

Constructor for the FE2DSetupMain class.

Parameters:

Name	Type	Description	Default
cell_type	str	Type of finite element ('quadrilateral')	required
fe_order	int	Order of finite element approximation (1 < order < 1e3)	required
fe_type	str	Type of basis functions ('legendre', 'legendre_special', 'chebyshev_2', 'jacobi_plain')	required
quad_order	int	Order of quadrature rule (>= 2)	required
quad_type	str	Type of quadrature formula	required

Returns:

Type	Description
	None

Source code in scirex/core/sciml/fe/fe2d_setup_main.py

def __init__(
    self,
    cell_type: str,
    fe_order: int,
    fe_type: str,
    quad_order: int,
    quad_type: str,
):
    """
    Constructor for the FE2DSetupMain class.

    Args:
        cell_type (str): Type of finite element ('quadrilateral')
        fe_order (int): Order of finite element approximation (1 < order < 1e3)
        fe_type (str): Type of basis functions
            ('legendre', 'legendre_special', 'chebyshev_2', 'jacobi_plain')
        quad_order (int): Order of quadrature rule (>= 2)
        quad_type (str): Type of quadrature formula

    Raises:
        None

    Returns:
        None
    """
    self.cell_type = cell_type
    self.fe_order = fe_order
    self.fe_type = fe_type
    self.quad_order = quad_order
    self.quad_type = quad_type

    self.assign_basis_function()

assign_basis_function()

Assigns the basis function based on the cell type and the fe_order.

Returns:

Name	Type	Description
BasisFunction2D	BasisFunction2D	The basis function object for the given configuration.

Raises:

Type	Description
ValueError	If the fe order is invalid or the cell type is invalid.

Source code in scirex/core/sciml/fe/fe2d_setup_main.py

def assign_basis_function(self) -> BasisFunction2D:
    """
    Assigns the basis function based on the cell type and the fe_order.

    Args:
        None

    Returns:
        BasisFunction2D: The basis function object for the given configuration.

    Raises:
        ValueError: If the fe order is invalid or the cell type is invalid.
    """
    # check for fe order lower bound and higher bound
    if self.fe_order <= 1 or self.fe_order >= 1e3:
        print(
            f"Invalid fe order {self.fe_order} in {self.__class__.__name__} from {__name__}."
        )
        raise ValueError("fe order should be greater than 1 and less than 1e4.")

    if self.cell_type == "quadrilateral":
        self.n_nodes = 4

        # --- LEGENDRE --- #
        if self.fe_type == "legendre" or self.fe_type == "jacobi":
            # jacobi is added for backward compatibility with prev pushes
            # generally, jacobi is referred to as Legendre basis on previous iterations
            return Basis2DQNLegendre(self.fe_order**2)

        elif self.fe_type == "legendre_special":
            return Basis2DQNLegendreSpecial(self.fe_order**2)

        # ----- CHEBYSHEV ---- #
        elif self.fe_type == "chebyshev_2":
            return Basis2DQNChebyshev2(self.fe_order**2)

        # ----- PLain jacobi ---- #
        elif self.fe_type == "jacobi_plain":
            return Basis2DQNJacobi(self.fe_order**2)

        else:
            print(
                f"Invalid fe order {self.fe_order} in {self.__class__.__name__} from {__name__}."
            )
            raise ValueError(
                'fe order should be one of the : "legendre" , "jacobi", "legendre_special", "chebyshev_2", "jacobi_plain"'
            )

    print(
        f"Invalid cell type {self.cell_type} in {self.__class__.__name__} from {__name__}."
    )

assign_fe_transformation(fe_transformation_type, cell_coordinates)

Assigns the fe transformation based on the cell type.

Parameters:

Name	Type	Description	Default
fe_transformation_type	str	Type of fe transformation ('affine', 'bilinear')	required
cell_coordinates	ndarray	The cell coordinates	required

Returns:

Name	Type	Description
FETransforamtion2D	FETransforamtion2D	The fe transformation object for the given configuration.

Raises:

Type	Description
ValueError	If the cell type is invalid
ValueError	If the fe transformation type is invalid

Source code in scirex/core/sciml/fe/fe2d_setup_main.py

def assign_fe_transformation(
    self, fe_transformation_type: str, cell_coordinates: np.ndarray
) -> FETransforamtion2D:
    """
    Assigns the fe transformation based on the cell type.

    Args:
        fe_transformation_type (str): Type of fe transformation ('affine', 'bilinear')
        cell_coordinates (np.ndarray): The cell coordinates

    Returns:
        FETransforamtion2D: The fe transformation object for the given configuration.

    Raises:
        ValueError: If the cell type is invalid
        ValueError: If the fe transformation type is invalid
    """
    if self.cell_type == "quadrilateral":
        if fe_transformation_type == "affine":
            return QuadAffin(cell_coordinates)
        elif fe_transformation_type == "bilinear":
            return QuadBilinear(cell_coordinates)
        else:
            raise ValueError(
                f"Invalid fe transformation type {fe_transformation_type} in {self.__class__.__name__} from {__name__}."
            )

    else:
        raise ValueError(
            f"Invalid cell type {self.cell_type} in {self.__class__.__name__} from {__name__}."
        )

assign_quadrature_rules()

Assigns the quadrature rule based on the quad_order.

Returns:

Name	Type	Description
tuple		The quadrature weights, xi and eta values in a numpy array format.

Raises:

Type	Description
ValueError	If the quad_order is invalid
ValueError	If the cell type is invalid
ValueError	If the quad_order is not between 1 and 9999

Source code in scirex/core/sciml/fe/fe2d_setup_main.py

def assign_quadrature_rules(self):
    """
    Assigns the quadrature rule based on the quad_order.

    Args:
        None

    Returns:
        tuple: The quadrature weights, xi and eta values in a numpy array format.

    Raises:
        ValueError: If the quad_order is invalid
        ValueError: If the cell type is invalid
        ValueError: If the quad_order is not between 1 and 9999
    """
    if self.cell_type == "quadrilateral":
        if self.quad_order < 3:
            raise ValueError("Quad order should be greater than 2.")
        elif self.quad_order >= 2 and self.quad_order <= 9999:
            weights, xi, eta = Quadratureformulas_Quad2D(
                self.quad_order, self.quad_type
            ).get_quad_values()
            return weights, xi, eta
        else:
            print(
                f"Invalid quad order {self.quad_order} in {self.__class__.__name__} from {__name__}."
            )
            raise ValueError("Quad order should be between 1 and 9999.")

    raise ValueError(
        f"Invalid cell type {self.cell_type} in {self.__class__.__name__} from {__name__}."
    )

FE2D_Cell

A class for managing finite element computations at the cell level.

This class handles the storage and computation of finite element values, including basis functions, quadrature rules, and transformations for a single cell in a 2D mesh.

Attributes:

Name	Type	Description
cell_coordinates	ndarray	Physical coordinates of the cell vertices
cell_type	str	Type of the cell (e.g., 'quad', 'triangle')
fe_order	int	Order of the finite element approximation
fe_type	str	Type of finite element basis
quad_order	int	Order of quadrature rule
quad_type	str	Type of quadrature formula
fe_transformation	str	Type of geometric transformation
forcing_function	callable	Source term function
basis_function	BasisFunction2D	Basis function implementation
quad_xi	ndarray	Xi coordinates of quadrature points
quad_eta	ndarray	Eta coordinates of quadrature points
quad_weight	ndarray	Quadrature weights
jacobian	ndarray	Transformation Jacobian
basis_at_quad	ndarray	Basis values at quadrature points
basis_gradx_at_quad	ndarray	X-derivatives at quadrature points
basis_grady_at_quad	ndarray	Y-derivatives at quadrature points
quad_actual_coordinates	ndarray	Physical quadrature point coordinates

Example

coords = np.array([[0,0], [1,0], [1,1], [0,1]]) cell = FE2D_Cell( ... cell_coordinates=coords, ... cell_type='quad', ... fe_order=2, ... fe_type='lagrange', ... quad_order=3, ... quad_type='gauss', ... fe_transformation_type='bilinear', ... forcing_function=lambda x, y: x*y ... ) cell.basis_at_quad # Get basis values at quadrature points

Notes

• All gradient and derivative values are stored in the reference domain
• Jacobian and quadrature weights are combined for efficiency
• Forcing function values are typically computed in the fespace class
• Supports multiple types of transformations and element types

Source code in scirex/core/sciml/fe/fe2d_cell.py

class FE2D_Cell:
    """
    A class for managing finite element computations at the cell level.

    This class handles the storage and computation of finite element values,
    including basis functions, quadrature rules, and transformations for a
    single cell in a 2D mesh.

    Attributes:
        cell_coordinates (np.ndarray): Physical coordinates of the cell vertices
        cell_type (str): Type of the cell (e.g., 'quad', 'triangle')
        fe_order (int): Order of the finite element approximation
        fe_type (str): Type of finite element basis
        quad_order (int): Order of quadrature rule
        quad_type (str): Type of quadrature formula
        fe_transformation (str): Type of geometric transformation
        forcing_function (callable): Source term function
        basis_function (BasisFunction2D): Basis function implementation
        quad_xi (np.ndarray): Xi coordinates of quadrature points
        quad_eta (np.ndarray): Eta coordinates of quadrature points
        quad_weight (np.ndarray): Quadrature weights
        jacobian (np.ndarray): Transformation Jacobian
        basis_at_quad (np.ndarray): Basis values at quadrature points
        basis_gradx_at_quad (np.ndarray): X-derivatives at quadrature points
        basis_grady_at_quad (np.ndarray): Y-derivatives at quadrature points
        quad_actual_coordinates (np.ndarray): Physical quadrature point coordinates

    Example:
        >>> coords = np.array([[0,0], [1,0], [1,1], [0,1]])
        >>> cell = FE2D_Cell(
        ...     cell_coordinates=coords,
        ...     cell_type='quad',
        ...     fe_order=2,
        ...     fe_type='lagrange',
        ...     quad_order=3,
        ...     quad_type='gauss',
        ...     fe_transformation_type='bilinear',
        ...     forcing_function=lambda x, y: x*y
        ... )
        >>> cell.basis_at_quad  # Get basis values at quadrature points

    Notes:
        - All gradient and derivative values are stored in the reference domain
        - Jacobian and quadrature weights are combined for efficiency
        - Forcing function values are typically computed in the fespace class
        - Supports multiple types of transformations and element types
    """

    def __init__(
        self,
        cell_coordinates: np.ndarray,
        cell_type: str,
        fe_order: int,
        fe_type: str,
        quad_order: int,
        quad_type: str,
        fe_transformation_type: str,
        forcing_function,
    ):
        """
        Constructor for the FE2D_Cell class.

        Args:
            cell_coordinates (np.ndarray): Physical coordinates of the cell vertices
            cell_type (str): Type of the cell (e.g., 'quad', 'triangle')
            fe_order (int): Order of the finite element approximation
            fe_type (str): Type of finite element basis
            quad_order (int): Order of quadrature rule
            quad_type (str): Type of quadrature formula
            fe_transformation_type (str): Type of geometric transformation
            forcing_function (callable): Source term function

        Returns:
            None
        """
        self.cell_coordinates = cell_coordinates
        self.cell_type = cell_type
        self.fe_order = fe_order
        self.fe_type = fe_type
        self.quad_order = quad_order
        self.quad_type = quad_type
        self.fe_transformation = fe_transformation_type
        self.forcing_function = forcing_function

        # Basis function Class
        self.basis_function = None

        # Quadrature Values
        self.quad_xi = None
        self.quad_eta = None
        self.quad_weight = None
        self.jacobian = None
        self.mult = None

        # FE Values
        self.basis_at_quad = None
        self.basis_gradx_at_quad = None
        self.basis_grady_at_quad = None
        self.basis_gradxy_at_quad = None
        self.basis_gradxx_at_quad = None
        self.basis_gradyy_at_quad = None

        # Quadrature Coordinates
        self.quad_actual_coordinates = None

        # Forcing function values at the quadrature points
        self.forcing_at_quad = None

        # FE Transformation Class
        self.fetransformation = None

        # get instance of the FE_setup class
        self.fe_setup = FE2DSetupMain(
            cell_type=self.cell_type,
            fe_order=self.fe_order,
            fe_type=self.fe_type,
            quad_order=self.quad_order,
            quad_type=self.quad_type,
        )

        # Call the function to assign the basis function
        self.assign_basis_function()

        # Assign the quadrature points and weights
        self.assign_quadrature()

        # Assign the FE Transformation
        self.assign_fe_transformation()

        # calculate mult -> quadrature weights * Jacobian
        self.assign_quad_weights_and_jacobian()

        # Calculate the basis function values at the quadrature points
        self.assign_basis_values_at_quadrature_points()

        # calculate the actual coordinates of the quadrature points
        self.assign_quadrature_coordinates()

        # Calculate the forcing function values at the actual quadrature points
        # NOTE : The function is just for printing the shape of the force matrix, the
        # actual calculation is performed on the fespace class
        self.assign_forcing_term(self.forcing_function)

        # # print the values
        # print("============================================================================")
        # print("Cell Co-ord : ", self.cell_coordinates)
        # print("Basis function values at the quadrature points: \n", self.basis_at_quad / self.mult)
        # print("Basis function gradx at the quadrature points: \n", self.basis_gradx_at_quad)
        # print("Basis function grady at the quadrature points: \n", self.basis_grady_at_quad)
        # print("Forcing function values at the quadrature points: \n", self.forcing_at_quad)

        # grad_x = np.array([5,6,7,8])
        # grad_y = np.array([1,2,3,4])

        # pde = np.matmul(self.basis_gradx_at_quad, grad_x.reshape(-1,1)) + np.matmul(self.basis_grady_at_quad, grad_y.reshape(-1,1))
        # print("PDE values at the quadrature points: \n", pde)

    def assign_basis_function(self) -> BasisFunction2D:
        """
        Assigns the basis function class based on the cell type and the FE order.

        Args:
            None

        Returns:
            BasisFunction2D: The basis function class for the given cell type and FE order.
        """
        self.basis_function = self.fe_setup.assign_basis_function()

    def assign_quadrature(self) -> None:
        """
        Assigns the quadrature points and weights based on the cell type and the quadrature order.

        Args:
            None

        Returns:
            None
        """
        self.quad_weight, self.quad_xi, self.quad_eta = (
            self.fe_setup.assign_quadrature_rules()
        )

    def assign_fe_transformation(self) -> None:
        """
        Assigns the FE Transformation class based on the cell type and the FE order.

        This method assigns the appropriate FE Transformation class based on the cell type and the FE order.
        It sets the cell coordinates for the FE Transformation and obtains the Jacobian of the transformation.

        Args:
            None

        Returns:
            None
        """
        self.fetransformation = self.fe_setup.assign_fe_transformation(
            self.fe_transformation, self.cell_coordinates
        )
        # Sets cell co-ordinates for the FE Transformation
        self.fetransformation.set_cell()

        # obtains the Jacobian of the transformation
        self.jacobian = self.fetransformation.get_jacobian(
            self.quad_xi, self.quad_eta
        ).reshape(-1, 1)

    def assign_basis_values_at_quadrature_points(self) -> None:
        """
        Assigns the basis function values at the quadrature points.

        This method calculates the values of the basis functions and their gradients at the quadrature points.
        The basis function values are stored in `self.basis_at_quad`, while the gradients are stored in
        `self.basis_gradx_at_quad`, `self.basis_grady_at_quad`, `self.basis_gradxy_at_quad`,
        `self.basis_gradxx_at_quad`, and `self.basis_gradyy_at_quad`.

        Args:
            None

        Returns:
            None
        """
        self.basis_at_quad = []
        self.basis_gradx_at_quad = []
        self.basis_grady_at_quad = []
        self.basis_gradxy_at_quad = []
        self.basis_gradxx_at_quad = []
        self.basis_gradyy_at_quad = []

        self.basis_at_quad = self.basis_function.value(self.quad_xi, self.quad_eta)

        # For Gradients we need to perform a transformation to the original cell
        grad_x_ref = self.basis_function.gradx(self.quad_xi, self.quad_eta)
        grad_y_ref = self.basis_function.grady(self.quad_xi, self.quad_eta)

        grad_x_orig, grad_y_orig = self.fetransformation.get_orig_from_ref_derivative(
            grad_x_ref, grad_y_ref, self.quad_xi, self.quad_eta
        )

        self.basis_gradx_at_quad = grad_x_orig
        self.basis_grady_at_quad = grad_y_orig

        self.basis_gradx_at_quad_ref = grad_x_ref
        self.basis_grady_at_quad_ref = grad_y_ref

        # get the double derivatives of the basis functions ( ref co-ordinates )
        grad_xx_ref = self.basis_function.gradxx(self.quad_xi, self.quad_eta)
        grad_xy_ref = self.basis_function.gradxy(self.quad_xi, self.quad_eta)
        grad_yy_ref = self.basis_function.gradyy(self.quad_xi, self.quad_eta)

        # get the double derivatives of the basis functions ( orig co-ordinates )
        grad_xx_orig, grad_xy_orig, grad_yy_orig = (
            self.fetransformation.get_orig_from_ref_second_derivative(
                grad_xx_ref, grad_xy_ref, grad_yy_ref, self.quad_xi, self.quad_eta
            )
        )

        # = the value
        self.basis_gradxy_at_quad = grad_xy_orig
        self.basis_gradxx_at_quad = grad_xx_orig
        self.basis_gradyy_at_quad = grad_yy_orig

        # Multiply each row with the quadrature weights
        # Basis at Quad - n_test * N_quad
        self.basis_at_quad = self.basis_at_quad * self.mult
        self.basis_gradx_at_quad = self.basis_gradx_at_quad * self.mult
        self.basis_grady_at_quad = self.basis_grady_at_quad * self.mult
        self.basis_gradxy_at_quad = self.basis_gradxy_at_quad * self.mult
        self.basis_gradxx_at_quad = self.basis_gradxx_at_quad * self.mult
        self.basis_gradyy_at_quad = self.basis_gradyy_at_quad * self.mult

    def assign_quad_weights_and_jacobian(self) -> None:
        """
        Assigns the quadrature weights and the Jacobian of the transformation.

        This method calculates and assigns the quadrature weights and the Jacobian of the transformation
        for the current cell. The quadrature weights are multiplied by the flattened Jacobian array
        and stored in the `mult` attribute of the class.

        Args:
            None

        Returns:
            None
        """
        self.mult = self.quad_weight * self.jacobian.flatten()

    def assign_quadrature_coordinates(self) -> None:
        """
        Assigns the actual coordinates of the quadrature points.

        This method calculates the actual coordinates of the quadrature points based on the given Xi and Eta values.
        The Xi and Eta values are obtained from the `quad_xi` and `quad_eta` attributes of the class.
        The calculated coordinates are stored in the `quad_actual_coordinates` attribute as a NumPy array.

        Args:
            None

        Returns:
            None
        """
        actual_co_ord = []
        for xi, eta in zip(self.quad_xi, self.quad_eta):
            actual_co_ord.append(self.fetransformation.get_original_from_ref(xi, eta))

        self.quad_actual_coordinates = np.array(actual_co_ord)

    def assign_forcing_term(self, forcing_function) -> None:
        """
        Assigns the forcing function values at the quadrature points.

        This function computes the values of the forcing function at the quadrature points
        and assigns them to the `forcing_at_quad` attribute of the FE2D_Cell object.

        Args:
            forcing_function (callable): The forcing function to be integrated

        Returns:
            None

        Notes:
            - The final shape of `forcing_at_quad` will be N_shape_functions x 1.
            - This function is for backward compatibility with old code and currently assigns
              the values as zeros. The actual calculation is performed in the fespace class.
        """
        # get number of shape functions
        n_shape_functions = self.basis_function.num_shape_functions

        # Loop over all the basis functions and compute the integral
        f_integral = np.zeros((n_shape_functions, 1), dtype=np.float64)

        # The above code is for backward compatibility with old code. this function will just assign the values as zeros
        # the actual calculation is performed in the fespace class

        # for i in range(n_shape_functions):
        #     val = 0
        #     for q in range(self.basis_at_quad.shape[1]):
        #         x = self.quad_actual_coordinates[q, 0]
        #         y = self.quad_actual_coordinates[q, 1]
        #         # print("f_values[q] = ",f_values[q])

        #         # the JAcobian and the quadrature weights are pre multiplied to the basis functions
        #         val +=  ( self.basis_at_quad[i, q] ) * self.forcing_function(x, y)
        #         # print("val = ", val)

        #     f_integral[i] = val

        self.forcing_at_quad = f_integral

__init__(cell_coordinates, cell_type, fe_order, fe_type, quad_order, quad_type, fe_transformation_type, forcing_function)

Constructor for the FE2D_Cell class.

Parameters:

Name	Type	Description	Default
cell_coordinates	ndarray	Physical coordinates of the cell vertices	required
cell_type	str	Type of the cell (e.g., 'quad', 'triangle')	required
fe_order	int	Order of the finite element approximation	required
fe_type	str	Type of finite element basis	required
quad_order	int	Order of quadrature rule	required
quad_type	str	Type of quadrature formula	required
fe_transformation_type	str	Type of geometric transformation	required
forcing_function	callable	Source term function	required

Returns:

Type	Description
	None

Source code in scirex/core/sciml/fe/fe2d_cell.py

def __init__(
    self,
    cell_coordinates: np.ndarray,
    cell_type: str,
    fe_order: int,
    fe_type: str,
    quad_order: int,
    quad_type: str,
    fe_transformation_type: str,
    forcing_function,
):
    """
    Constructor for the FE2D_Cell class.

    Args:
        cell_coordinates (np.ndarray): Physical coordinates of the cell vertices
        cell_type (str): Type of the cell (e.g., 'quad', 'triangle')
        fe_order (int): Order of the finite element approximation
        fe_type (str): Type of finite element basis
        quad_order (int): Order of quadrature rule
        quad_type (str): Type of quadrature formula
        fe_transformation_type (str): Type of geometric transformation
        forcing_function (callable): Source term function

    Returns:
        None
    """
    self.cell_coordinates = cell_coordinates
    self.cell_type = cell_type
    self.fe_order = fe_order
    self.fe_type = fe_type
    self.quad_order = quad_order
    self.quad_type = quad_type
    self.fe_transformation = fe_transformation_type
    self.forcing_function = forcing_function

    # Basis function Class
    self.basis_function = None

    # Quadrature Values
    self.quad_xi = None
    self.quad_eta = None
    self.quad_weight = None
    self.jacobian = None
    self.mult = None

    # FE Values
    self.basis_at_quad = None
    self.basis_gradx_at_quad = None
    self.basis_grady_at_quad = None
    self.basis_gradxy_at_quad = None
    self.basis_gradxx_at_quad = None
    self.basis_gradyy_at_quad = None

    # Quadrature Coordinates
    self.quad_actual_coordinates = None

    # Forcing function values at the quadrature points
    self.forcing_at_quad = None

    # FE Transformation Class
    self.fetransformation = None

    # get instance of the FE_setup class
    self.fe_setup = FE2DSetupMain(
        cell_type=self.cell_type,
        fe_order=self.fe_order,
        fe_type=self.fe_type,
        quad_order=self.quad_order,
        quad_type=self.quad_type,
    )

    # Call the function to assign the basis function
    self.assign_basis_function()

    # Assign the quadrature points and weights
    self.assign_quadrature()

    # Assign the FE Transformation
    self.assign_fe_transformation()

    # calculate mult -> quadrature weights * Jacobian
    self.assign_quad_weights_and_jacobian()

    # Calculate the basis function values at the quadrature points
    self.assign_basis_values_at_quadrature_points()

    # calculate the actual coordinates of the quadrature points
    self.assign_quadrature_coordinates()

    # Calculate the forcing function values at the actual quadrature points
    # NOTE : The function is just for printing the shape of the force matrix, the
    # actual calculation is performed on the fespace class
    self.assign_forcing_term(self.forcing_function)

assign_basis_function()

Assigns the basis function class based on the cell type and the FE order.

Returns:

Name	Type	Description
BasisFunction2D	BasisFunction2D	The basis function class for the given cell type and FE order.

Source code in scirex/core/sciml/fe/fe2d_cell.py

def assign_basis_function(self) -> BasisFunction2D:
    """
    Assigns the basis function class based on the cell type and the FE order.

    Args:
        None

    Returns:
        BasisFunction2D: The basis function class for the given cell type and FE order.
    """
    self.basis_function = self.fe_setup.assign_basis_function()

assign_basis_values_at_quadrature_points()

Assigns the basis function values at the quadrature points.

This method calculates the values of the basis functions and their gradients at the quadrature points. The basis function values are stored in self.basis_at_quad, while the gradients are stored in self.basis_gradx_at_quad, self.basis_grady_at_quad, self.basis_gradxy_at_quad, self.basis_gradxx_at_quad, and self.basis_gradyy_at_quad.

Returns:

Type	Description
None	None

Source code in scirex/core/sciml/fe/fe2d_cell.py

def assign_basis_values_at_quadrature_points(self) -> None:
    """
    Assigns the basis function values at the quadrature points.

    This method calculates the values of the basis functions and their gradients at the quadrature points.
    The basis function values are stored in `self.basis_at_quad`, while the gradients are stored in
    `self.basis_gradx_at_quad`, `self.basis_grady_at_quad`, `self.basis_gradxy_at_quad`,
    `self.basis_gradxx_at_quad`, and `self.basis_gradyy_at_quad`.

    Args:
        None

    Returns:
        None
    """
    self.basis_at_quad = []
    self.basis_gradx_at_quad = []
    self.basis_grady_at_quad = []
    self.basis_gradxy_at_quad = []
    self.basis_gradxx_at_quad = []
    self.basis_gradyy_at_quad = []

    self.basis_at_quad = self.basis_function.value(self.quad_xi, self.quad_eta)

    # For Gradients we need to perform a transformation to the original cell
    grad_x_ref = self.basis_function.gradx(self.quad_xi, self.quad_eta)
    grad_y_ref = self.basis_function.grady(self.quad_xi, self.quad_eta)

    grad_x_orig, grad_y_orig = self.fetransformation.get_orig_from_ref_derivative(
        grad_x_ref, grad_y_ref, self.quad_xi, self.quad_eta
    )

    self.basis_gradx_at_quad = grad_x_orig
    self.basis_grady_at_quad = grad_y_orig

    self.basis_gradx_at_quad_ref = grad_x_ref
    self.basis_grady_at_quad_ref = grad_y_ref

    # get the double derivatives of the basis functions ( ref co-ordinates )
    grad_xx_ref = self.basis_function.gradxx(self.quad_xi, self.quad_eta)
    grad_xy_ref = self.basis_function.gradxy(self.quad_xi, self.quad_eta)
    grad_yy_ref = self.basis_function.gradyy(self.quad_xi, self.quad_eta)

    # get the double derivatives of the basis functions ( orig co-ordinates )
    grad_xx_orig, grad_xy_orig, grad_yy_orig = (
        self.fetransformation.get_orig_from_ref_second_derivative(
            grad_xx_ref, grad_xy_ref, grad_yy_ref, self.quad_xi, self.quad_eta
        )
    )

    # = the value
    self.basis_gradxy_at_quad = grad_xy_orig
    self.basis_gradxx_at_quad = grad_xx_orig
    self.basis_gradyy_at_quad = grad_yy_orig

    # Multiply each row with the quadrature weights
    # Basis at Quad - n_test * N_quad
    self.basis_at_quad = self.basis_at_quad * self.mult
    self.basis_gradx_at_quad = self.basis_gradx_at_quad * self.mult
    self.basis_grady_at_quad = self.basis_grady_at_quad * self.mult
    self.basis_gradxy_at_quad = self.basis_gradxy_at_quad * self.mult
    self.basis_gradxx_at_quad = self.basis_gradxx_at_quad * self.mult
    self.basis_gradyy_at_quad = self.basis_gradyy_at_quad * self.mult

assign_fe_transformation()

Assigns the FE Transformation class based on the cell type and the FE order.

This method assigns the appropriate FE Transformation class based on the cell type and the FE order. It sets the cell coordinates for the FE Transformation and obtains the Jacobian of the transformation.

Returns:

Type	Description
None	None

Source code in scirex/core/sciml/fe/fe2d_cell.py

def assign_fe_transformation(self) -> None:
    """
    Assigns the FE Transformation class based on the cell type and the FE order.

    This method assigns the appropriate FE Transformation class based on the cell type and the FE order.
    It sets the cell coordinates for the FE Transformation and obtains the Jacobian of the transformation.

    Args:
        None

    Returns:
        None
    """
    self.fetransformation = self.fe_setup.assign_fe_transformation(
        self.fe_transformation, self.cell_coordinates
    )
    # Sets cell co-ordinates for the FE Transformation
    self.fetransformation.set_cell()

    # obtains the Jacobian of the transformation
    self.jacobian = self.fetransformation.get_jacobian(
        self.quad_xi, self.quad_eta
    ).reshape(-1, 1)

assign_forcing_term(forcing_function)

Assigns the forcing function values at the quadrature points.

This function computes the values of the forcing function at the quadrature points and assigns them to the forcing_at_quad attribute of the FE2D_Cell object.

Parameters:

Name	Type	Description	Default
forcing_function	callable	The forcing function to be integrated	required

Returns:

Type	Description
None	None

Notes

• The final shape of forcing_at_quad will be N_shape_functions x 1.
• This function is for backward compatibility with old code and currently assigns the values as zeros. The actual calculation is performed in the fespace class.

Source code in scirex/core/sciml/fe/fe2d_cell.py

def assign_forcing_term(self, forcing_function) -> None:
    """
    Assigns the forcing function values at the quadrature points.

    This function computes the values of the forcing function at the quadrature points
    and assigns them to the `forcing_at_quad` attribute of the FE2D_Cell object.

    Args:
        forcing_function (callable): The forcing function to be integrated

    Returns:
        None

    Notes:
        - The final shape of `forcing_at_quad` will be N_shape_functions x 1.
        - This function is for backward compatibility with old code and currently assigns
          the values as zeros. The actual calculation is performed in the fespace class.
    """
    # get number of shape functions
    n_shape_functions = self.basis_function.num_shape_functions

    # Loop over all the basis functions and compute the integral
    f_integral = np.zeros((n_shape_functions, 1), dtype=np.float64)

    # The above code is for backward compatibility with old code. this function will just assign the values as zeros
    # the actual calculation is performed in the fespace class

    # for i in range(n_shape_functions):
    #     val = 0
    #     for q in range(self.basis_at_quad.shape[1]):
    #         x = self.quad_actual_coordinates[q, 0]
    #         y = self.quad_actual_coordinates[q, 1]
    #         # print("f_values[q] = ",f_values[q])

    #         # the JAcobian and the quadrature weights are pre multiplied to the basis functions
    #         val +=  ( self.basis_at_quad[i, q] ) * self.forcing_function(x, y)
    #         # print("val = ", val)

    #     f_integral[i] = val

    self.forcing_at_quad = f_integral

assign_quad_weights_and_jacobian()

Assigns the quadrature weights and the Jacobian of the transformation.

This method calculates and assigns the quadrature weights and the Jacobian of the transformation for the current cell. The quadrature weights are multiplied by the flattened Jacobian array and stored in the mult attribute of the class.

Returns:

Type	Description
None	None

Source code in scirex/core/sciml/fe/fe2d_cell.py

def assign_quad_weights_and_jacobian(self) -> None:
    """
    Assigns the quadrature weights and the Jacobian of the transformation.

    This method calculates and assigns the quadrature weights and the Jacobian of the transformation
    for the current cell. The quadrature weights are multiplied by the flattened Jacobian array
    and stored in the `mult` attribute of the class.

    Args:
        None

    Returns:
        None
    """
    self.mult = self.quad_weight * self.jacobian.flatten()

assign_quadrature()

Assigns the quadrature points and weights based on the cell type and the quadrature order.

Returns:

Type	Description
None	None

Source code in scirex/core/sciml/fe/fe2d_cell.py

def assign_quadrature(self) -> None:
    """
    Assigns the quadrature points and weights based on the cell type and the quadrature order.

    Args:
        None

    Returns:
        None
    """
    self.quad_weight, self.quad_xi, self.quad_eta = (
        self.fe_setup.assign_quadrature_rules()
    )

assign_quadrature_coordinates()

Assigns the actual coordinates of the quadrature points.

This method calculates the actual coordinates of the quadrature points based on the given Xi and Eta values. The Xi and Eta values are obtained from the quad_xi and quad_eta attributes of the class. The calculated coordinates are stored in the quad_actual_coordinates attribute as a NumPy array.

Returns:

Type	Description
None	None

Source code in scirex/core/sciml/fe/fe2d_cell.py

def assign_quadrature_coordinates(self) -> None:
    """
    Assigns the actual coordinates of the quadrature points.

    This method calculates the actual coordinates of the quadrature points based on the given Xi and Eta values.
    The Xi and Eta values are obtained from the `quad_xi` and `quad_eta` attributes of the class.
    The calculated coordinates are stored in the `quad_actual_coordinates` attribute as a NumPy array.

    Args:
        None

    Returns:
        None
    """
    actual_co_ord = []
    for xi, eta in zip(self.quad_xi, self.quad_eta):
        actual_co_ord.append(self.fetransformation.get_original_from_ref(xi, eta))

    self.quad_actual_coordinates = np.array(actual_co_ord)

FETransforamtion2D

A base class for 2D finite element transformations.

This abstract class defines the interface for mapping between reference and physical coordinates in 2D finite element analysis. Implementations must provide specific transformation rules for different element types.

Methods:

Name	Description
set_cell	Sets the physical coordinates of the element vertices. Must be implemented by derived classes.
get_original_from_ref	Maps coordinates from reference to physical domain. Must be implemented by derived classes.
get_jacobian	Computes the Jacobian matrix of the transformation. Must be implemented by derived classes.

Example

class QuadTransform(FETransformation2D): ... def set_cell(self, vertices): ... self.vertices = vertices ... def get_original_from_ref(self, xi:np.ndarray, eta:np.ndarray) -> np.ndarray: ... # Implementation for quad element ... pass ... def get_jacobian(self, xi: np.ndarray, eta:np.ndarray) -> np.ndarray: ... # Implementation for quad element ... pass

Notes

• Reference domain is typically [-1,1] × [-1,1]
• Transformations must be invertible
• Implementations should handle element distortion
• Jacobian is used for both mapping and integration

Source code in scirex/core/sciml/fe/fe_transformation_2d.py

class FETransforamtion2D:
    """
    A base class for 2D finite element transformations.

    This abstract class defines the interface for mapping between reference and physical
    coordinates in 2D finite element analysis. Implementations must provide specific
    transformation rules for different element types.

    Attributes:
        None

    Methods:
        set_cell():
            Sets the physical coordinates of the element vertices.
            Must be implemented by derived classes.

        get_original_from_ref(xi, eta):
            Maps coordinates from reference to physical domain.
            Must be implemented by derived classes.

        get_jacobian(xi, eta):
            Computes the Jacobian matrix of the transformation.
            Must be implemented by derived classes.

    Example:
        >>> class QuadTransform(FETransformation2D):
        ...     def set_cell(self, vertices):
        ...         self.vertices = vertices
        ...     def get_original_from_ref(self, xi:np.ndarray, eta:np.ndarray) -> np.ndarray:
        ...         # Implementation for quad element
        ...         pass
        ...     def get_jacobian(self, xi: np.ndarray, eta:np.ndarray) -> np.ndarray:
        ...         # Implementation for quad element
        ...         pass

    Notes:
        - Reference domain is typically [-1,1] × [-1,1]
        - Transformations must be invertible
        - Implementations should handle element distortion
        - Jacobian is used for both mapping and integration
    """

    def __init__(self):
        """
        Constructor for the FETransforamtion2D class.
        """

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

        :return: None
        """

    @abstractmethod
    def get_original_from_ref(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        This method returns the original coordinates from the reference coordinates.

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

        Returns:
            np.ndarray: Returns the transformed original coordinates from the reference coordinates.
        """

    @abstractmethod
    def get_jacobian(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
        """
        This method returns the Jacobian of the transformation.

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

        Returns:
            np.ndarray: Returns the Jacobian of the transformation.
        """

__init__()

Constructor for the FETransforamtion2D class.

Source code in scirex/core/sciml/fe/fe_transformation_2d.py

def __init__(self):
    """
    Constructor for the FETransforamtion2D class.
    """

get_jacobian(xi, eta) abstractmethod

This method 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: Returns the Jacobian of the transformation.

Source code in scirex/core/sciml/fe/fe_transformation_2d.py

@abstractmethod
def get_jacobian(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    This method returns the Jacobian of the transformation.

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

    Returns:
        np.ndarray: Returns the Jacobian of the transformation.
    """

get_original_from_ref(xi, eta) abstractmethod

This method 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: Returns the transformed original coordinates from the reference coordinates.

Source code in scirex/core/sciml/fe/fe_transformation_2d.py

@abstractmethod
def get_original_from_ref(self, xi: np.ndarray, eta: np.ndarray) -> np.ndarray:
    """
    This method returns the original coordinates from the reference coordinates.

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

    Returns:
        np.ndarray: Returns the transformed original coordinates from the reference coordinates.
    """

set_cell() abstractmethod

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

:return: None

Source code in scirex/core/sciml/fe/fe_transformation_2d.py

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

    :return: None
    """

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

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

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

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

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

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

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

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

QuadBilinear

Bases: FETransforamtion2D

Implements bilinear transformation for quadrilateral elements.

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

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, xc3)	x-coordinate transformation coefficients
yc0,	(yc1, yc2, yc3)	y-coordinate transformation coefficients
detjk		Determinant of the Jacobian matrix

Example

coords = np.array([[0,0], [1,0], [1.2,1], [0.2,1.1]]) quad = QuadBilinear(coords) ref_point = np.array([0.5, 0.5]) physical_point = quad.get_original_from_ref(*ref_point)

Note

• Implementation assumes counterclockwise vertex ordering
• Second derivatives computation is not fully implemented
• Jacobian is computed point-wise due to non-constant nature of bilinear transformation

References

[1] ParMooN Project: QuadBilineare.C implementation

Source code in scirex/core/sciml/fe/quad_bilinear.py

class QuadBilinear(FETransforamtion2D):
    """
    Implements bilinear transformation for quadrilateral elements.

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

    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, xc3: x-coordinate transformation coefficients
        yc0, yc1, yc2, yc3: y-coordinate transformation coefficients
        detjk: Determinant of the Jacobian matrix

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

    Note:
        - Implementation assumes counterclockwise vertex ordering
        - Second derivatives computation is not fully implemented
        - Jacobian is computed point-wise due to non-constant nature
        of bilinear transformation

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

    def __init__(self, co_ordinates: np.ndarray) -> None:
        """
        Constructor for the QuadBilinear 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.detjk = None  # Jacobian of the transformation

    def set_cell(self):
        """
        Set the cell coordinates, which will be used as intermediate values 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.x0 + self.x1 + self.x2 + self.x3) * 0.25
        self.xc1 = (-self.x0 + self.x1 + self.x2 - self.x3) * 0.25
        self.xc2 = (-self.x0 - self.x1 + self.x2 + self.x3) * 0.25
        self.xc3 = (self.x0 - self.x1 + self.x2 - self.x3) * 0.25

        self.yc0 = (self.y0 + self.y1 + self.y2 + self.y3) * 0.25
        self.yc1 = (-self.y0 + self.y1 + self.y2 - self.y3) * 0.25
        self.yc2 = (-self.y0 - self.y1 + self.y2 + self.y3) * 0.25
        self.yc3 = (self.y0 - self.y1 + self.y2 - self.y3) * 0.25

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

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

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

        return np.array([x, y], dtype=np.float64)

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

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

        Returns:
            np.ndarray: Returns the Jacobian of the transformation.
        """
        self.detjk = abs(
            (self.xc1 + self.xc3 * eta) * (self.yc2 + self.yc3 * xi)
            - (self.xc2 + self.xc3 * xi) * (self.yc1 + self.yc3 * eta)
        )
        return self.detjk

    def get_orig_from_ref_derivative(
        self,
        ref_gradx: np.ndarray,
        ref_grady: np.ndarray,
        xi: np.ndarray,
        eta: np.ndarray,
    ) -> np.ndarray:
        """
        This method returns the derivatives of the original coordinates with respect to the reference coordinates.

        Args:
            ref_gradx (np.ndarray): The derivative of the xi coordinate in the reference element.
            ref_grady (np.ndarray): The derivative of the eta coordinate in the reference element.
            xi (np.ndarray): The xi coordinate.
            eta (np.ndarray): The eta coordinate.

        Returns:
            np.ndarray: The derivatives of the original coordinates [x, y] with respect to the reference coordinates.

        """
        n_test = ref_gradx.shape[0]
        gradx_orig = np.zeros(ref_gradx.shape, dtype=np.float64)
        grady_orig = np.zeros(ref_grady.shape, dtype=np.float64)

        for j in range(n_test):
            Xi = xi
            Eta = eta
            rec_detjk = 1 / (
                (self.xc1 + self.xc3 * Eta) * (self.yc2 + self.yc3 * Xi)
                - (self.xc2 + self.xc3 * Xi) * (self.yc1 + self.yc3 * Eta)
            )
            gradx_orig[j] = (
                (self.yc2 + self.yc3 * Xi) * ref_gradx[j]
                - (self.yc1 + self.yc3 * Eta) * ref_grady[j]
            ) * rec_detjk
            grady_orig[j] = (
                -(self.xc2 + self.xc3 * Xi) * ref_gradx[j]
                + (self.xc1 + self.xc3 * Eta) * ref_grady[j]
            ) * rec_detjk

        return gradx_orig, grady_orig

    def get_orig_from_ref_second_derivative(
        self,
        grad_xx_ref: np.ndarray,
        grad_xy_ref: np.ndarray,
        grad_yy_ref: np.ndarray,
        xi: np.ndarray,
        eta: np.ndarray,
    ):
        """
        This method returns the second derivatives of the original coordinates with respect to the reference coordinates.

        Args:
            grad_xx_ref (np.ndarray): The second derivative of the xi coordinate in the reference element.
            grad_xy_ref (np.ndarray): The second derivative of the xi and eta coordinates in the reference element.
            grad_yy_ref (np.ndarray): The second derivative of the eta coordinate in the reference element.
            xi (np.ndarray): The xi coordinate.
            eta (np.ndarray): The eta coordinate.

        Note:
            Second derivative calculations are not fully implemented in this method. Needs further development.
        """
        # print(" Error : Second Derivative not implemented -- Ignore this error, if second derivative is not required ")
        return grad_xx_ref, grad_xy_ref, grad_yy_ref

__init__(co_ordinates)

Constructor for the QuadBilinear 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_bilinear.py

def __init__(self, co_ordinates: np.ndarray) -> None:
    """
    Constructor for the QuadBilinear 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.detjk = None  # Jacobian of the transformation

get_jacobian(xi, eta)

This method 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: Returns the Jacobian of the transformation.

Source code in scirex/core/sciml/fe/quad_bilinear.py

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

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

    Returns:
        np.ndarray: Returns the Jacobian of the transformation.
    """
    self.detjk = abs(
        (self.xc1 + self.xc3 * eta) * (self.yc2 + self.yc3 * xi)
        - (self.xc2 + self.xc3 * xi) * (self.yc1 + self.yc3 * eta)
    )
    return self.detjk

get_orig_from_ref_derivative(ref_gradx, ref_grady, xi, eta)

This method returns the derivatives of the original coordinates with respect to the reference coordinates.

Parameters:

Name	Type	Description	Default
ref_gradx	ndarray	The derivative of the xi coordinate in the reference element.	required
ref_grady	ndarray	The derivative of the eta coordinate in the reference element.	required
xi	ndarray	The xi coordinate.	required
eta	ndarray	The eta coordinate.	required

Returns:

Type	Description
ndarray	np.ndarray: The derivatives of the original coordinates [x, y] with respect to the reference coordinates.

Source code in scirex/core/sciml/fe/quad_bilinear.py

def get_orig_from_ref_derivative(
    self,
    ref_gradx: np.ndarray,
    ref_grady: np.ndarray,
    xi: np.ndarray,
    eta: np.ndarray,
) -> np.ndarray:
    """
    This method returns the derivatives of the original coordinates with respect to the reference coordinates.

    Args:
        ref_gradx (np.ndarray): The derivative of the xi coordinate in the reference element.
        ref_grady (np.ndarray): The derivative of the eta coordinate in the reference element.
        xi (np.ndarray): The xi coordinate.
        eta (np.ndarray): The eta coordinate.

    Returns:
        np.ndarray: The derivatives of the original coordinates [x, y] with respect to the reference coordinates.

    """
    n_test = ref_gradx.shape[0]
    gradx_orig = np.zeros(ref_gradx.shape, dtype=np.float64)
    grady_orig = np.zeros(ref_grady.shape, dtype=np.float64)

    for j in range(n_test):
        Xi = xi
        Eta = eta
        rec_detjk = 1 / (
            (self.xc1 + self.xc3 * Eta) * (self.yc2 + self.yc3 * Xi)
            - (self.xc2 + self.xc3 * Xi) * (self.yc1 + self.yc3 * Eta)
        )
        gradx_orig[j] = (
            (self.yc2 + self.yc3 * Xi) * ref_gradx[j]
            - (self.yc1 + self.yc3 * Eta) * ref_grady[j]
        ) * rec_detjk
        grady_orig[j] = (
            -(self.xc2 + self.xc3 * Xi) * ref_gradx[j]
            + (self.xc1 + self.xc3 * Eta) * ref_grady[j]
        ) * rec_detjk

    return gradx_orig, grady_orig

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

This method returns the second derivatives of the original coordinates with respect to the reference coordinates.

Parameters:

Name	Type	Description	Default
grad_xx_ref	ndarray	The second derivative of the xi coordinate in the reference element.	required
grad_xy_ref	ndarray	The second derivative of the xi and eta coordinates in the reference element.	required
grad_yy_ref	ndarray	The second derivative of the eta coordinate in the reference element.	required
xi	ndarray	The xi coordinate.	required
eta	ndarray	The eta coordinate.	required

Note

Second derivative calculations are not fully implemented in this method. Needs further development.

Source code in scirex/core/sciml/fe/quad_bilinear.py

def get_orig_from_ref_second_derivative(
    self,
    grad_xx_ref: np.ndarray,
    grad_xy_ref: np.ndarray,
    grad_yy_ref: np.ndarray,
    xi: np.ndarray,
    eta: np.ndarray,
):
    """
    This method returns the second derivatives of the original coordinates with respect to the reference coordinates.

    Args:
        grad_xx_ref (np.ndarray): The second derivative of the xi coordinate in the reference element.
        grad_xy_ref (np.ndarray): The second derivative of the xi and eta coordinates in the reference element.
        grad_yy_ref (np.ndarray): The second derivative of the eta coordinate in the reference element.
        xi (np.ndarray): The xi coordinate.
        eta (np.ndarray): The eta coordinate.

    Note:
        Second derivative calculations are not fully implemented in this method. Needs further development.
    """
    # print(" Error : Second Derivative not implemented -- Ignore this error, if second derivative is not required ")
    return grad_xx_ref, grad_xy_ref, grad_yy_ref

get_original_from_ref(xi, eta)

This method 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: Returns the transformed original coordinates from the reference coordinates.

Source code in scirex/core/sciml/fe/quad_bilinear.py

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

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

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

    return np.array([x, y], dtype=np.float64)

set_cell()

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

Returns:

Type	Description
	None

Source code in scirex/core/sciml/fe/quad_bilinear.py

def set_cell(self):
    """
    Set the cell coordinates, which will be used as intermediate values 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.x0 + self.x1 + self.x2 + self.x3) * 0.25
    self.xc1 = (-self.x0 + self.x1 + self.x2 - self.x3) * 0.25
    self.xc2 = (-self.x0 - self.x1 + self.x2 + self.x3) * 0.25
    self.xc3 = (self.x0 - self.x1 + self.x2 - self.x3) * 0.25

    self.yc0 = (self.y0 + self.y1 + self.y2 + self.y3) * 0.25
    self.yc1 = (-self.y0 + self.y1 + self.y2 - self.y3) * 0.25
    self.yc2 = (-self.y0 - self.y1 + self.y2 + self.y3) * 0.25
    self.yc3 = (self.y0 - self.y1 + self.y2 - self.y3) * 0.25

Quadratureformulas

Abstract base class for numerical quadrature formulas.

This class defines the interface that all quadrature implementations must follow. It provides the basic structure for implementing various quadrature rules while ensuring consistent access to quadrature data.

Attributes:

Name	Type	Description
quad_order		Order of the quadrature rule
quad_type		Type of quadrature (e.g., 'gauss-legendre', 'gauss-jacobi')
num_quad_points		Total number of quadrature points

Example

class MyQuadrature(Quadratureformulas): ... def init(self): ... super().init(quad_order=3, ... quad_type='custom', ... num_quad_points=9) ... def get_quad_values(self): ... # Implementation ... pass ... def get_num_quad_points(self): ... return self.num_quad_points

Note

This is an abstract base class. Concrete implementations must override: - get_quad_values() - get_num_quad_points()

The implementation should ensure proper initialization of: - Quadrature points - Quadrature weights - Number of quadrature points

Source code in scirex/core/sciml/fe/quadratureformulas.py

class Quadratureformulas:
    """Abstract base class for numerical quadrature formulas.

    This class defines the interface that all quadrature implementations must
    follow. It provides the basic structure for implementing various quadrature
    rules while ensuring consistent access to quadrature data.

    Attributes:
        quad_order: Order of the quadrature rule
        quad_type: Type of quadrature (e.g., 'gauss-legendre', 'gauss-jacobi')
        num_quad_points: Total number of quadrature points

    Example:
        >>> class MyQuadrature(Quadratureformulas):
        ...     def __init__(self):
        ...         super().__init__(quad_order=3,
        ...                         quad_type='custom',
        ...                         num_quad_points=9)
        ...     def get_quad_values(self):
        ...         # Implementation
        ...         pass
        ...     def get_num_quad_points(self):
        ...         return self.num_quad_points

    Note:
        This is an abstract base class. Concrete implementations must override:
        - get_quad_values()
        - get_num_quad_points()

        The implementation should ensure proper initialization of:
        - Quadrature points
        - Quadrature weights
        - Number of quadrature points
    """

    def __init__(self, quad_order: int, quad_type: str, num_quad_points: int):
        """
        Constructor for the Quadratureformulas_Quad2D class.

        Args:
            quad_order: Order of quadrature rule
            quad_type: Type of quadrature ('gauss-legendre' or 'gauss-jacobi')
            num_quad_points: Total number of quadrature points

        Returns:
            None
        """
        self.quad_order = quad_order
        self.quad_type = quad_type
        self.num_quad_points = num_quad_points

    @abstractmethod
    def get_quad_values(self):
        """
        Returns the quadrature weights, xi and eta values.

        Args:
            None

        Returns:
            weights: Weights for each quadrature point
            xi: x-coordinates of quadrature points in reference element
            eta: y-coordinates of quadrature points in reference element
        """

    @abstractmethod
    def get_num_quad_points(self):
        """
        Returns the number of quadrature points.

        Args:
            None

        Returns:
            num_quad_points: Total number of quadrature points
        """

__init__(quad_order, quad_type, num_quad_points)

Constructor for the Quadratureformulas_Quad2D class.

Parameters:

Name	Type	Description	Default
quad_order	int	Order of quadrature rule	required
quad_type	str	Type of quadrature ('gauss-legendre' or 'gauss-jacobi')	required
num_quad_points	int	Total number of quadrature points	required

Returns:

Type	Description
	None

Source code in scirex/core/sciml/fe/quadratureformulas.py

def __init__(self, quad_order: int, quad_type: str, num_quad_points: int):
    """
    Constructor for the Quadratureformulas_Quad2D class.

    Args:
        quad_order: Order of quadrature rule
        quad_type: Type of quadrature ('gauss-legendre' or 'gauss-jacobi')
        num_quad_points: Total number of quadrature points

    Returns:
        None
    """
    self.quad_order = quad_order
    self.quad_type = quad_type
    self.num_quad_points = num_quad_points

get_num_quad_points() abstractmethod

Returns the number of quadrature points.

Returns:

Name	Type	Description
num_quad_points		Total number of quadrature points

Source code in scirex/core/sciml/fe/quadratureformulas.py

@abstractmethod
def get_num_quad_points(self):
    """
    Returns the number of quadrature points.

    Args:
        None

    Returns:
        num_quad_points: Total number of quadrature points
    """

get_quad_values() abstractmethod

Returns the quadrature weights, xi and eta values.

Returns:

Name	Type	Description
weights		Weights for each quadrature point
xi		x-coordinates of quadrature points in reference element
eta		y-coordinates of quadrature points in reference element

Source code in scirex/core/sciml/fe/quadratureformulas.py

@abstractmethod
def get_quad_values(self):
    """
    Returns the quadrature weights, xi and eta values.

    Args:
        None

    Returns:
        weights: Weights for each quadrature point
        xi: x-coordinates of quadrature points in reference element
        eta: y-coordinates of quadrature points in reference element
    """

Quadratureformulas_Quad2D

Bases: Quadratureformulas

Implements quadrature formulas for 2D quadrilateral elements.

This class provides methods to compute quadrature points and weights for 2D quadrilateral elements using either Gauss-Legendre or Gauss-Jacobi quadrature schemes. The implementation uses tensor products of 1D rules.

Attributes:

Name	Type	Description
quad_order		Order of quadrature rule
quad_type		Type of quadrature ('gauss-legendre' or 'gauss-jacobi')
num_quad_points		Total number of quadrature points (quad_order^2)
xi_quad		x-coordinates of quadrature points in reference element
eta_quad		y-coordinates of quadrature points in reference element
quad_weights		Weights for each quadrature point

Example

quad = Quadratureformulas_Quad2D(quad_order=3, quad_type='gauss-legendre') weights, xi, eta = quad.get_quad_values() n_points = quad.get_num_quad_points()

Note

• Gauss-Legendre points are optimal for polynomial integrands
• Gauss-Jacobi points include element vertices (useful for certain FEM applications)
• All computations are performed in the reference element [-1,1]×[-1,1]

Source code in scirex/core/sciml/fe/quadratureformulas_quad2d.py

class Quadratureformulas_Quad2D(Quadratureformulas):
    """Implements quadrature formulas for 2D quadrilateral elements.

    This class provides methods to compute quadrature points and weights for
    2D quadrilateral elements using either Gauss-Legendre or Gauss-Jacobi
    quadrature schemes. The implementation uses tensor products of 1D rules.

    Attributes:
        quad_order: Order of quadrature rule
        quad_type: Type of quadrature ('gauss-legendre' or 'gauss-jacobi')
        num_quad_points: Total number of quadrature points (quad_order^2)
        xi_quad: x-coordinates of quadrature points in reference element
        eta_quad: y-coordinates of quadrature points in reference element
        quad_weights: Weights for each quadrature point

    Example:
        >>> quad = Quadratureformulas_Quad2D(quad_order=3, quad_type='gauss-legendre')
        >>> weights, xi, eta = quad.get_quad_values()
        >>> n_points = quad.get_num_quad_points()

    Note:
        - Gauss-Legendre points are optimal for polynomial integrands
        - Gauss-Jacobi points include element vertices (useful for certain FEM applications)
        - All computations are performed in the reference element [-1,1]×[-1,1]

    """

    def __init__(self, quad_order: int, quad_type: str):
        """
        Constructor for the Quadratureformulas_Quad2D class.

        Args:
            quad_order: Order of quadrature rule
            quad_type: Type of quadrature ('gauss-legendre' or 'gauss-jacobi')

        Returns:
            None

        Raises:
            ValueError: If the quadrature type is not supported.
        """
        # initialize the super class
        super().__init__(
            quad_order=quad_order,
            quad_type=quad_type,
            num_quad_points=quad_order * quad_order,
        )

        # Calculate the Gauss-Legendre quadrature points and weights for 1D
        # nodes_1d, weights_1d = roots_jacobi(self.quad_order, 1, 1)

        quad_type = self.quad_type

        if quad_type == "gauss-legendre":
            # Commented out by THIVIN -  to Just use legendre quadrature points as it is
            # if quad_order == 2:
            #     nodes_1d = np.array([-1, 1])
            #     weights_1d = np.array([1, 1])
            # else:
            nodes_1d, weights_1d = np.polynomial.legendre.leggauss(
                quad_order
            )  # Interior points
            # nodes_1d = np.concatenate(([-1, 1], nodes_1d))
            # weights_1d = np.concatenate(([1, 1], weights_1d))

            # Generate the tensor outer product of the nodes
            xi_quad, eta_quad = np.meshgrid(nodes_1d, nodes_1d)
            xi_quad = xi_quad.flatten()
            eta_quad = eta_quad.flatten()

            # Multiply the weights accordingly for 2D
            quad_weights = (weights_1d[:, np.newaxis] * weights_1d).flatten()

            # Assign the values
            self.xi_quad = xi_quad
            self.eta_quad = eta_quad
            self.quad_weights = quad_weights

        elif quad_type == "gauss-jacobi":

            def GaussJacobiWeights(Q: int, a, b):
                [X, W] = roots_jacobi(Q, a, b)
                return [X, W]

            def jacobi_wrapper(n, a, b, x):

                x = np.array(x, dtype=np.float64)

                return jacobi(n, a, b)(x)

            # Weight coefficients
            def GaussLobattoJacobiWeights(Q: int, a, b):
                W = []
                X = roots_jacobi(Q - 2, a + 1, b + 1)[0]
                if a == 0 and b == 0:
                    W = 2 / ((Q - 1) * (Q) * (jacobi_wrapper(Q - 1, 0, 0, X) ** 2))
                    Wl = 2 / ((Q - 1) * (Q) * (jacobi_wrapper(Q - 1, 0, 0, -1) ** 2))
                    Wr = 2 / ((Q - 1) * (Q) * (jacobi_wrapper(Q - 1, 0, 0, 1) ** 2))
                else:
                    W = (
                        2 ** (a + b + 1)
                        * gamma(a + Q)
                        * gamma(b + Q)
                        / (
                            (Q - 1)
                            * gamma(Q)
                            * gamma(a + b + Q + 1)
                            * (jacobi_wrapper(Q - 1, a, b, X) ** 2)
                        )
                    )
                    Wl = (
                        (b + 1)
                        * 2 ** (a + b + 1)
                        * gamma(a + Q)
                        * gamma(b + Q)
                        / (
                            (Q - 1)
                            * gamma(Q)
                            * gamma(a + b + Q + 1)
                            * (jacobi_wrapper(Q - 1, a, b, -1) ** 2)
                        )
                    )
                    Wr = (
                        (a + 1)
                        * 2 ** (a + b + 1)
                        * gamma(a + Q)
                        * gamma(b + Q)
                        / (
                            (Q - 1)
                            * gamma(Q)
                            * gamma(a + b + Q + 1)
                            * (jacobi_wrapper(Q - 1, a, b, 1) ** 2)
                        )
                    )
                W = np.append(W, Wr)
                W = np.append(Wl, W)
                X = np.append(X, 1)
                X = np.append(-1, X)
                return [X, W]

            # get quadrature points and weights in 1D
            x, w = GaussLobattoJacobiWeights(self.quad_order, 0, 0)

            # Generate the tensor outer product of the nodes
            xi_quad, eta_quad = np.meshgrid(x, x)
            xi_quad = xi_quad.flatten()
            eta_quad = eta_quad.flatten()

            # Multiply the weights accordingly for 2D
            quad_weights = (w[:, np.newaxis] * w).flatten()

            # Assign the values
            self.xi_quad = xi_quad
            self.eta_quad = eta_quad
            self.quad_weights = quad_weights

        else:
            print("Supported quadrature types are: gauss-legendre, gauss-jacobi")
            print(
                f"Invalid quadrature type {quad_type} in {self.__class__.__name__} from {__name__}."
            )
            raise ValueError("Quadrature type not supported.")

    def get_quad_values(self):
        """
        Returns the quadrature weights, xi and eta values.

        Args:
            None

        Returns:
            tuple: The quadrature weights, xi and eta values in a numpy array format
        """
        return self.quad_weights, self.xi_quad, self.eta_quad

    def get_num_quad_points(self):
        """
        Returns the number of quadrature points.

        Args:
            None

        Returns:
            int: The number of quadrature points
        """
        return self.num_quad_points

__init__(quad_order, quad_type)

Constructor for the Quadratureformulas_Quad2D class.

Parameters:

Name	Type	Description	Default
quad_order	int	Order of quadrature rule	required
quad_type	str	Type of quadrature ('gauss-legendre' or 'gauss-jacobi')	required

Returns:

Type	Description
	None

Raises:

Type	Description
ValueError	If the quadrature type is not supported.

Source code in scirex/core/sciml/fe/quadratureformulas_quad2d.py

def __init__(self, quad_order: int, quad_type: str):
    """
    Constructor for the Quadratureformulas_Quad2D class.

    Args:
        quad_order: Order of quadrature rule
        quad_type: Type of quadrature ('gauss-legendre' or 'gauss-jacobi')

    Returns:
        None

    Raises:
        ValueError: If the quadrature type is not supported.
    """
    # initialize the super class
    super().__init__(
        quad_order=quad_order,
        quad_type=quad_type,
        num_quad_points=quad_order * quad_order,
    )

    # Calculate the Gauss-Legendre quadrature points and weights for 1D
    # nodes_1d, weights_1d = roots_jacobi(self.quad_order, 1, 1)

    quad_type = self.quad_type

    if quad_type == "gauss-legendre":
        # Commented out by THIVIN -  to Just use legendre quadrature points as it is
        # if quad_order == 2:
        #     nodes_1d = np.array([-1, 1])
        #     weights_1d = np.array([1, 1])
        # else:
        nodes_1d, weights_1d = np.polynomial.legendre.leggauss(
            quad_order
        )  # Interior points
        # nodes_1d = np.concatenate(([-1, 1], nodes_1d))
        # weights_1d = np.concatenate(([1, 1], weights_1d))

        # Generate the tensor outer product of the nodes
        xi_quad, eta_quad = np.meshgrid(nodes_1d, nodes_1d)
        xi_quad = xi_quad.flatten()
        eta_quad = eta_quad.flatten()

        # Multiply the weights accordingly for 2D
        quad_weights = (weights_1d[:, np.newaxis] * weights_1d).flatten()

        # Assign the values
        self.xi_quad = xi_quad
        self.eta_quad = eta_quad
        self.quad_weights = quad_weights

    elif quad_type == "gauss-jacobi":

        def GaussJacobiWeights(Q: int, a, b):
            [X, W] = roots_jacobi(Q, a, b)
            return [X, W]

        def jacobi_wrapper(n, a, b, x):

            x = np.array(x, dtype=np.float64)

            return jacobi(n, a, b)(x)

        # Weight coefficients
        def GaussLobattoJacobiWeights(Q: int, a, b):
            W = []
            X = roots_jacobi(Q - 2, a + 1, b + 1)[0]
            if a == 0 and b == 0:
                W = 2 / ((Q - 1) * (Q) * (jacobi_wrapper(Q - 1, 0, 0, X) ** 2))
                Wl = 2 / ((Q - 1) * (Q) * (jacobi_wrapper(Q - 1, 0, 0, -1) ** 2))
                Wr = 2 / ((Q - 1) * (Q) * (jacobi_wrapper(Q - 1, 0, 0, 1) ** 2))
            else:
                W = (
                    2 ** (a + b + 1)
                    * gamma(a + Q)
                    * gamma(b + Q)
                    / (
                        (Q - 1)
                        * gamma(Q)
                        * gamma(a + b + Q + 1)
                        * (jacobi_wrapper(Q - 1, a, b, X) ** 2)
                    )
                )
                Wl = (
                    (b + 1)
                    * 2 ** (a + b + 1)
                    * gamma(a + Q)
                    * gamma(b + Q)
                    / (
                        (Q - 1)
                        * gamma(Q)
                        * gamma(a + b + Q + 1)
                        * (jacobi_wrapper(Q - 1, a, b, -1) ** 2)
                    )
                )
                Wr = (
                    (a + 1)
                    * 2 ** (a + b + 1)
                    * gamma(a + Q)
                    * gamma(b + Q)
                    / (
                        (Q - 1)
                        * gamma(Q)
                        * gamma(a + b + Q + 1)
                        * (jacobi_wrapper(Q - 1, a, b, 1) ** 2)
                    )
                )
            W = np.append(W, Wr)
            W = np.append(Wl, W)
            X = np.append(X, 1)
            X = np.append(-1, X)
            return [X, W]

        # get quadrature points and weights in 1D
        x, w = GaussLobattoJacobiWeights(self.quad_order, 0, 0)

        # Generate the tensor outer product of the nodes
        xi_quad, eta_quad = np.meshgrid(x, x)
        xi_quad = xi_quad.flatten()
        eta_quad = eta_quad.flatten()

        # Multiply the weights accordingly for 2D
        quad_weights = (w[:, np.newaxis] * w).flatten()

        # Assign the values
        self.xi_quad = xi_quad
        self.eta_quad = eta_quad
        self.quad_weights = quad_weights

    else:
        print("Supported quadrature types are: gauss-legendre, gauss-jacobi")
        print(
            f"Invalid quadrature type {quad_type} in {self.__class__.__name__} from {__name__}."
        )
        raise ValueError("Quadrature type not supported.")

get_num_quad_points()

Returns the number of quadrature points.

Returns:

Name	Type	Description
int		The number of quadrature points

Source code in scirex/core/sciml/fe/quadratureformulas_quad2d.py

def get_num_quad_points(self):
    """
    Returns the number of quadrature points.

    Args:
        None

    Returns:
        int: The number of quadrature points
    """
    return self.num_quad_points

get_quad_values()

Returns the quadrature weights, xi and eta values.

Returns:

Name	Type	Description
tuple		The quadrature weights, xi and eta values in a numpy array format

Source code in scirex/core/sciml/fe/quadratureformulas_quad2d.py

def get_quad_values(self):
    """
    Returns the quadrature weights, xi and eta values.

    Args:
        None

    Returns:
        tuple: The quadrature weights, xi and eta values in a numpy array format
    """
    return self.quad_weights, self.xi_quad, self.eta_quad