Source code for pysofe.quadrature.gaussian

"""
Provides the data structures for Gauss-Legendre quadrature on simplicial domains.
"""

# IMPORTS
import numpy as np

from .base import QuadRule


[docs]class GaussQuadSimp(QuadRule): """ Provides Gauss-Legendre quadrature points and weights for standard simplicial domains in several spatial dimensions. Parameters ---------- order : int The polynomial order up to which the numerical integration should be exact dimension : int The spatial dimension up to which to provide quadrature data """ def __init__(self, order, dimension): QuadRule.__init__(self, order, dimension) def _set_data(self): """ Generates the quadrature points and weights. """ if self.dimension >= 0: # 0-dimensional data (just for conveinience) self._points[0] = np.empty((0,0)) self._weights[0] = np.asarray([1.]) if self.dimension >= 1: # 1-dimensional data are the Gauss-Legendre points and weights points, weights = np.polynomial.legendre.leggauss(self.order+1) # the points and weights are defined for the interval [-1,1] # so we need to transform them to our reference interval [0,1] self._points[1] = np.atleast_2d(0.5 * (points + 1.)) self._weights[1] = 0.5 * weights if self.dimension >= 2: # 2-dimensional points and weights are symmetrical gaussian quadrature # rules of D.A. Dunavant # (tabulated for orders 1-10) if self.order > 10: raise ValueError("Quadrature order not supported, yet!") points_weights = np.asarray(dunavant_points_weights[self.order-1]) points = points_weights.take((0,1), axis=1).T weights = points_weights.take(2, axis=1) # weights are normalized to the triangle area weights *= 0.5 self._points[2] = points self._weights[2] = weights # # DEPRECATED # #------------ # points, weights = np.polynomial.legendre.leggauss(self.order+1) # points = 0.5 * (points + 1.) # weights = 0.5 * weights # # Duffy Transform # W = np.outer(weights, weights).ravel(order='F') # P = np.vstack([np.tile(points, (np.size(weights, 0),)), # np.repeat(points, repeats=np.size(weights, 0), axis=0)]) # self._points[2] = np.vstack([P[0] * (1 - P[1]), P[1]]) # self._weights[2] = (W * (1 - P[1])) if self.dimension >= 3: # 3-dimensional points and weights are not yet supported # for orders higher than 2 if self.order == 1: self._points[3] = np.array([[0.25], [0.25], [0.25]]) self._weights[3] = np.array([1/6.]) elif self.order == 2: a1 = 0.138196601150000 a2 = 0.585410196600000 self._points[3] = np.array([[a1, a2, a1, a1], [a1, a1, a2, a1], [a1, a1, a1, a2]]) self._weights[3] = 1/24. * np.array([1., 1., 1., 1.]) else: msg = "Quadrature points and weights on tetrahedron not supported for order {} (>2)!" raise NotImplementedError(msg.format(self.order))
#-------------------------------------- # TABLE OF DUNAVANT'S QUADRATURE RULES #-------------------------------------- dunavant_points_weights = np.array([ [[0.33333333333333, 0.33333333333333, 1.00000000000000]], [[0.16666666666667, 0.16666666666667, 0.33333333333333], [0.16666666666667, 0.66666666666667, 0.33333333333333], [0.66666666666667, 0.16666666666667, 0.33333333333333]], [[0.33333333333333, 0.33333333333333, -0.56250000000000], [0.20000000000000, 0.20000000000000, 0.52083333333333], [0.20000000000000, 0.60000000000000, 0.52083333333333], [0.60000000000000, 0.20000000000000, 0.52083333333333]], [[0.44594849091597, 0.44594849091597, 0.22338158967801], [0.44594849091597, 0.10810301816807, 0.22338158967801], [0.10810301816807, 0.44594849091597, 0.22338158967801], [0.09157621350977, 0.09157621350977, 0.10995174365532], [0.09157621350977, 0.81684757298046, 0.10995174365532], [0.81684757298046, 0.09157621350977, 0.10995174365532]], [[0.33333333333333, 0.33333333333333, 0.22500000000000], [0.47014206410511, 0.47014206410511, 0.13239415278851], [0.47014206410511, 0.05971587178977, 0.13239415278851], [0.05971587178977, 0.47014206410511, 0.13239415278851], [0.10128650732346, 0.10128650732346, 0.12593918054483], [0.10128650732346, 0.79742698535309, 0.12593918054483], [0.79742698535309, 0.10128650732346, 0.12593918054483]], [[0.24928674517091, 0.24928674517091, 0.11678627572638], [0.24928674517091, 0.50142650965818, 0.11678627572638], [0.50142650965818, 0.24928674517091, 0.11678627572638], [0.06308901449150, 0.06308901449150, 0.05084490637021], [0.06308901449150, 0.87382197101700, 0.05084490637021], [0.87382197101700, 0.06308901449150, 0.05084490637021], [0.31035245103378, 0.63650249912140, 0.08285107561837], [0.63650249912140, 0.05314504984482, 0.08285107561837], [0.05314504984482, 0.31035245103378, 0.08285107561837], [0.63650249912140, 0.31035245103378, 0.08285107561837], [0.31035245103378, 0.05314504984482, 0.08285107561837], [0.05314504984482, 0.63650249912140, 0.08285107561837]], [[0.33333333333333, 0.33333333333333, -0.14957004446768], [0.26034596607904, 0.26034596607904, 0.17561525743321], [0.26034596607904, 0.47930806784192, 0.17561525743321], [0.47930806784192, 0.26034596607904, 0.17561525743321], [0.06513010290222, 0.06513010290222, 0.05334723560884], [0.06513010290222, 0.86973979419557, 0.05334723560884], [0.86973979419557, 0.06513010290222, 0.05334723560884], [0.31286549600487, 0.63844418856981, 0.07711376089026], [0.63844418856981, 0.04869031542532, 0.07711376089026], [0.04869031542532, 0.31286549600487, 0.07711376089026], [0.63844418856981, 0.31286549600487, 0.07711376089026], [0.31286549600487, 0.04869031542532, 0.07711376089026], [0.04869031542532, 0.63844418856981, 0.07711376089026]], [[0.33333333333333, 0.33333333333333, 0.14431560767779], [0.45929258829272, 0.45929258829272, 0.09509163426728], [0.45929258829272, 0.08141482341455, 0.09509163426728], [0.08141482341455, 0.45929258829272, 0.09509163426728], [0.17056930775176, 0.17056930775176, 0.10321737053472], [0.17056930775176, 0.65886138449648, 0.10321737053472], [0.65886138449648, 0.17056930775176, 0.10321737053472], [0.05054722831703, 0.05054722831703, 0.03245849762320], [0.05054722831703, 0.89890554336594, 0.03245849762320], [0.89890554336594, 0.05054722831703, 0.03245849762320], [0.26311282963464, 0.72849239295540, 0.02723031417443], [0.72849239295540, 0.00839477740996, 0.02723031417443], [0.00839477740996, 0.26311282963464, 0.02723031417443], [0.72849239295540, 0.26311282963464, 0.02723031417443], [0.26311282963464, 0.00839477740996, 0.02723031417443], [0.00839477740996, 0.72849239295540, 0.02723031417443]], [[0.33333333333333, 0.33333333333333, 0.09713579628280], [0.48968251919874, 0.48968251919874, 0.03133470022714], [0.48968251919874, 0.02063496160252, 0.03133470022714], [0.02063496160252, 0.48968251919874, 0.03133470022714], [0.43708959149294, 0.43708959149294, 0.07782754100477], [0.43708959149294, 0.12582081701413, 0.07782754100477], [0.12582081701413, 0.43708959149294, 0.07782754100477], [0.18820353561903, 0.18820353561903, 0.07964773892721], [0.18820353561903, 0.62359292876193, 0.07964773892721], [0.62359292876193, 0.18820353561903, 0.07964773892721], [0.04472951339445, 0.04472951339445, 0.02557767565870], [0.04472951339445, 0.91054097321109, 0.02557767565870], [0.91054097321109, 0.04472951339445, 0.02557767565870], [0.22196298916077, 0.74119859878450, 0.04328353937729], [0.74119859878450, 0.03683841205474, 0.04328353937729], [0.03683841205474, 0.22196298916077, 0.04328353937729], [0.74119859878450, 0.22196298916077, 0.04328353937729], [0.22196298916077, 0.03683841205474, 0.04328353937729], [0.03683841205474, 0.74119859878450, 0.04328353937729]], [[0.33333333333333, 0.33333333333333, 0.09081799038275], [0.48557763338366, 0.48557763338366, 0.03672595775647], [0.48557763338366, 0.02884473323269, 0.03672595775647], [0.02884473323269, 0.48557763338366, 0.03672595775647], [0.10948157548504, 0.10948157548504, 0.04532105943553], [0.10948157548504, 0.78103684902993, 0.04532105943553], [0.78103684902993, 0.10948157548504, 0.04532105943553], [0.30793983876412, 0.55035294182100, 0.07275791684542], [0.55035294182100, 0.14170721941488, 0.07275791684542], [0.14170721941488, 0.30793983876412, 0.07275791684542], [0.55035294182100, 0.30793983876412, 0.07275791684542], [0.30793983876412, 0.14170721941488, 0.07275791684542], [0.14170721941488, 0.55035294182100, 0.07275791684542], [0.24667256063990, 0.72832390459741, 0.02832724253106], [0.72832390459741, 0.02500353476269, 0.02832724253106], [0.02500353476269, 0.24667256063990, 0.02832724253106], [0.72832390459741, 0.24667256063990, 0.02832724253106], [0.24667256063990, 0.02500353476269, 0.02832724253106], [0.02500353476269, 0.72832390459741, 0.02832724253106], [0.06680325101220, 0.92365593358750, 0.00942166696373], [0.92365593358750, 0.00954081540030, 0.00942166696373], [0.00954081540030, 0.06680325101220, 0.00942166696373], [0.92365593358750, 0.06680325101220, 0.00942166696373], [0.06680325101220, 0.00954081540030, 0.00942166696373], [0.00954081540030, 0.92365593358750, 0.00942166696373]] ])