Source code for sfepy.discrete.quadratures

"""
`quadrature_tables` are organized as follows::

    quadrature_tables = {
        '<geometry1>' : {
            order1 : QuadraturePoints(args1),
            order2 : QuadraturePoints(args2),
            ...
        },
        '<geometry2>' : {
            order1 : QuadraturePoints(args1),
            order2 : QuadraturePoints(args2),
            ...
        },
        ...
    }

**Note** The order for quadratures on tensor product domains (`'2_4'`,
`'3_8'` geometries) in case of composite Gauss quadratures (products of
1D quadratures) holds for each component separately, so the actual
polynomial order may be much higher (up to `order * dimension`).

Naming conventions in problem description files::

    `<family>_<order>_<dimension>`

Integral 'family' is just an arbitrary name given by user.

Low order quadrature coordinates and weights copied from The Finite Element
Method Displayed by Gouri Dhatt and Gilbert Touzat, Wiley-Interscience
Production, 1984.

The line integral (geometry '1_2') coordinates and weights are from Abramowitz,
M. and Stegun, I.A., Handbook of Mathematical Functions, Dover Publications,
New York, 1972. The triangle (geometry '2_3') coordinates and weights are from
Dunavant, D.A., High Degree Efficient Symmetrical Gaussian Quadrature Rules for
the Triangle, Int. J. Num. Meth. Eng., 21 (1985) pp 1129-1148 - only rules with
points inside the reference triangle are used. The actual values were copied
from PHAML (http://math.nist.gov/phaml/), see also Mitchell, W.F., PHAML User's
Guide, NISTIR 7374, 2006.

Quadrature rules for the quadrilateral (geometry '2_4') and hexahedron
(geometry '3_8') of order higher than 5 are computed as the tensor product of
the line (geometry '1_2') rules.

Quadrature rules for the triangle (geometry '2_3') and tetrahedron (geometry
'3_4') of order higher than 19 and 6, respectively follow A. Grundmann and
H.M. Moeller, Invariant integration formulas for the n-simplex by combinatorial
methods, SIAM J. Numer. Anal.  15 (1978), 282--290. The generating function was
adapted from pytools/hegde codes (http://mathema.tician.de/software/hedge) by
Andreas Kloeckner.
"""
from __future__ import absolute_import
import numpy as nm

from sfepy.base.base import output, assert_, Struct
from sfepy.discrete.simplex_cubature import get_simplex_cubature

simplex_geometries = ['1_2', '2_3', '3_4']
tp_geometries = ['2_4', '3_8']

_msg1 = 'WARNING: quadrature order %s is not available for geometry %s!'
_msg2 = 'WARNING: using %d instead!'

[docs] def get_actual_order(geometry, order): """ Return the actual integration order for given geometry. Parameters ---------- geometry : str The geometry key describing the integration domain, see the keys of `quadrature_tables`. Returns ------- order : int If `order` is in quadrature tables it is this value. Otherwise it is the closest higher order. If no higher order is available, a warning is printed and the highest available order is used. """ table = quadrature_tables[geometry] if order not in table: orders = list(table.keys()) ii = nm.searchsorted(orders, order) if ii >= len(orders): omax = max(orders) output(_msg1 % (order, geometry)) output(_msg2 % omax) order = omax else: order = orders[ii] return order
[docs] class QuadraturePoints(Struct): """ Representation of a set of quadrature points. Parameters ---------- data : array_like The array of shape `(n_point, dim + 1)` of quadrature point coordinates (first `dim` columns) and weights (the last column). coors : array_like, optional Optionally, instead of using `data`, the coordinates and weights can be provided separately - `data` are then ignored. weights : array_like, optional Optionally, instead of using `data`, the coordinates and weights can be provided separately - `data` are then ignored. bounds : (float, float), optional The coordinates and weights should correspond to a reference element in `[0, 1]` x `dim`. Provide the correct bounds if this is not the case. tp_fix : float, optional The value that is used to multiply the tensor product element volume (= 1.0) to get the correct volume. weight_fix : float, optional The value that is used to multiply the weights to get the correct values. symmetric : bool If True, the integral is 1D and the given coordinates and weights are symmetric w.r.t. the centre of bounds; only the non-negative coordinates are given. """
[docs] @staticmethod def from_table(geometry, order): """ Create a new :class:`QuadraturePoints` instance, given reference element geometry name and polynomial order. For tensor product geometries, the polynomial order is the 1D (line) order. """ table = quadrature_tables[geometry] if geometry in simplex_geometries: if order > max_orders[geometry]: oo = order // 2 dim = int(geometry[0]) tp_fix = 0.5 if dim == 2 else 1.0 / 6.0 coors, weights, exact = get_simplex_cubature(oo, dim) qp = QuadraturePoints(None, coors=coors, weights=weights, bounds=(-1.0, 1.0), tp_fix=tp_fix) assert_(exact >= order) else: order = get_actual_order(geometry, order) qp = table[order] qp.order = order else: order1d = order dim = int(geometry[0]) order = dim * order1d if order <= max_orders[geometry]: order = get_actual_order(geometry, order) qp = table[order] qp.order = order else: oo = get_actual_order('1_2', order1d) qp1d = quadrature_tables['1_2'][oo] weights = nm.outer(qp1d.weights, qp1d.weights) nc = qp1d.coors.shape[0] if dim == 3: weights = nm.outer(qp1d.weights, weights) iz, iy, ix = nm.mgrid[0:nc, 0:nc, 0:nc] coors = nm.c_[qp1d.coors[ix.ravel()], qp1d.coors[iy.ravel()], qp1d.coors[iz.ravel()]].copy() else: iy, ix = nm.mgrid[0:nc, 0:nc] coors = nm.c_[qp1d.coors[ix.ravel()], qp1d.coors[iy.ravel()]].copy() weights = weights.ravel() qp = QuadraturePoints(None, coors=coors, weights=weights) qp.order = dim * oo return qp
def __init__(self, data, coors=None, weights=None, bounds=None, tp_fix=1.0, weight_fix=1.0, symmetric=False): if coors is None: data = nm.array(data, dtype=nm.float64, ndmin=2) self.coors = data[:,:-1].copy() self.weights = data[:,-1].copy() elif weights is not None: self.coors = nm.array(coors, dtype=nm.float64, ndmin=2) self.weights = nm.array(weights, dtype=nm.float64) else: raise ValueError('both "coors" and "weights" have to be provided!') self.weights *= weight_fix self.n_point, self.dim = self.coors.shape self.bounds = (0, 1) bbox = nm.array([self.bounds] * self.dim, dtype=nm.float64) self.volume = nm.prod(bbox.sum(axis=1)) * tp_fix if symmetric: isym = 0 if data[0, 0] == 0 else None if bounds is not None: # Transform from given bounds to self.bounds. bbox = nm.array([bounds] * self.dim, dtype=nm.float64) volume = nm.prod(nm.diff(bbox, axis=1)) * tp_fix a, b = bounds c, d = self.bounds c1 = (d - c) / (b - a) c2 = ((b * c) - (a * d)) / (b - a) self.coors = c1 * self.coors + c2 self.weights *= self.volume / volume if symmetric: if self.coors.shape[1] != 1: msg = 'symmetric mode is allowed for 1D integrals only!' raise ValueError(msg) origin = 0.5 * (self.bounds[0] + self.bounds[1]) self.coors = nm.r_[2 * origin - self.coors[:isym:-1], self.coors] self.weights = nm.r_[self.weights[:isym:-1], self.weights]
_QP = QuadraturePoints quadrature_tables = { '0_1' : { 1 : _QP([[0.0, 1.0]]) }, '1_2' : { 1 : _QP([[0.000000000000000e+00, 2.0]], bounds=(-1.0, 1.0), symmetric=True), 3 : _QP([[0.577350269189626e+00, 1.0]], bounds=(-1.0, 1.0), symmetric=True), 5 : _QP([[0.000000000000000e+00, 0.888888888888889e+00], [0.774596669241483e+00, 0.555555555555556e+00]], bounds=(-1.0, 1.0), symmetric=True), 7 : _QP([[0.339981043584856e+00, 0.652145154862546e+00], [0.861136311594053e+00, 0.347854845137454e+00]], bounds=(-1.0, 1.0), symmetric=True), 9 : _QP([[0.000000000000000e+00, 0.568888888888889e+00], [0.538469310105683e+00, 0.478628670499366e+00], [0.906179845938664e+00, 0.236926885056189e+00]], bounds=(-1.0, 1.0), symmetric=True), 11 : _QP([[0.238619186083197e+00, 0.467913934572691e+00], [0.661209386466265e+00, 0.360761573048139e+00], [0.932469514203152e+00, 0.171324492379170e+00]], bounds=(-1.0, 1.0), symmetric=True), 13 : _QP([[0.000000000000000e+00, 0.417959183673469e+00], [0.405845151377397e+00, 0.381830050505119e+00], [0.741531185599394e+00, 0.279705391489277e+00], [0.949107912342759e+00, 0.129484966168870e+00]], bounds=(-1.0, 1.0), symmetric=True), 15 : _QP([[0.183434642495650e+00, 0.362683783378362e+00], [0.525532409916329e+00, 0.313706645877887e+00], [0.796666477413627e+00, 0.222381034453374e+00], [0.960289856497536e+00, 0.101228536290376e+00]], bounds=(-1.0, 1.0), symmetric=True), 17 : _QP([[0.000000000000000e+00, 0.330239355001260e+00], [0.324253423403809e+00, 0.312347077040003e+00], [0.613371432700590e+00, 0.260610696402935e+00], [0.836031107326636e+00, 0.180648160694857e+00], [0.968160239507626e+00, 0.081274388361574e+00]], bounds=(-1.0, 1.0), symmetric=True), 19 : _QP([[0.148874338981631e+00, 0.295524224714753e+00], [0.433395394129247e+00, 0.269266719309996e+00], [0.679409568299024e+00, 0.219086362515982e+00], [0.865063366688985e+00, 0.149451349150581e+00], [0.973906528517172e+00, 0.066671344308688e+00]], bounds=(-1.0, 1.0), symmetric=True), 23 : _QP([[0.125233408511469e+00, 0.249147045813403e+00], [0.367831498998180e+00, 0.233492536538355e+00], [0.587317954286617e+00, 0.203167426723066e+00], [0.769902674194305e+00, 0.160078328543346e+00], [0.904117256370475e+00, 0.106939325995318e+00], [0.981560634246719e+00, 0.047175336386512e+00]], bounds=(-1.0, 1.0), symmetric=True), 31 : _QP([[0.095012509837637440185e+00, 0.189450610455068496285e+00], [0.281603550779258913230e+00, 0.182603415044923588867e+00], [0.458016777657227386342e+00, 0.169156519395002538189e+00], [0.617876244402643748447e+00, 0.149595988816576732081e+00], [0.755404408355003033895e+00, 0.124628971255533872052e+00], [0.865631202387831743880e+00, 0.095158511682492784810e+00], [0.944575023073232576078e+00, 0.062253523938647892863e+00], [0.989400934991649932596e+00, 0.027152459411754094852e+00]], bounds=(-1.0, 1.0), symmetric=True), 39 : _QP([[0.076526521133497333755e+00, 0.152753387130725850698e+00], [0.227785851141645078080e+00, 0.149172986472603746788e+00], [0.373706088715419560673e+00, 0.142096109318382051329e+00], [0.510867001950827098004e+00, 0.131688638449176626898e+00], [0.636053680726515025453e+00, 0.118194531961518417312e+00], [0.746331906460150792614e+00, 0.101930119817240435037e+00], [0.839116971822218823395e+00, 0.083276741576704748725e+00], [0.912234428251325905868e+00, 0.062672048334109063570e+00], [0.963971927277913791268e+00, 0.040601429800386941331e+00], [0.993128599185094924786e+00, 0.017614007139152118312e+00]], bounds=(-1.0, 1.0), symmetric=True), 47 : _QP([[0.064056892862605626085e+00, 0.127938195346752156974e+00], [0.191118867473616309159e+00, 0.125837456346828296121e+00], [0.315042679696163374387e+00, 0.121670472927803391204e+00], [0.433793507626045138487e+00, 0.115505668053725601353e+00], [0.545421471388839535658e+00, 0.107444270115965634783e+00], [0.648093651936975569252e+00, 0.097618652104113888270e+00], [0.740124191578554364244e+00, 0.086190161531953275917e+00], [0.820001985973902921954e+00, 0.073346481411080305734e+00], [0.886415527004401034213e+00, 0.059298584915436780746e+00], [0.938274552002732758524e+00, 0.044277438817419806169e+00], [0.974728555971309498198e+00, 0.028531388628933663181e+00], [0.995187219997021360180e+00, 0.012341229799987199547e+00]], bounds=(-1.0, 1.0), symmetric=True), }, '2_3' : { 1 : _QP([[1.0/3.0, 1.0/3.0, 0.5]], tp_fix=0.5), 2 : _QP([[1.0/6.0, 1.0/6.0, 1.0/6.0], [2.0/3.0, 1.0/6.0, 1.0/6.0], [1.0/6.0, 2.0/3.0, 1.0/6.0]], tp_fix=0.5), 3 : _QP([[1.0/3.0, 1.0/3.0,-27.0/96.0], [1.0/5.0, 1.0/5.0, 25.0/96.0], [3.0/5.0, 1.0/5.0, 25.0/96.0], [1.0/5.0, 3.0/5.0, 25.0/96.0]], tp_fix=0.5), 4 : _QP([[0.445948490915965e+00, 0.445948490915965e+00, 0.223381589678011e+00], [0.108103018168070e+00, 0.445948490915965e+00, 0.223381589678011e+00], [0.445948490915965e+00, 0.108103018168070e+00, 0.223381589678011e+00], [0.091576213509771e+00, 0.091576213509771e+00, 0.109951743655322e+00], [0.816847572980459e+00, 0.091576213509771e+00, 0.109951743655322e+00], [0.091576213509771e+00, 0.816847572980459e+00, 0.109951743655322e+00]], tp_fix=0.5, weight_fix=0.5), 5 : _QP([[0.333333333333333e+00, 0.333333333333333e+00, 0.225000000000000e+00], [0.470142064105115e+00, 0.470142064105115e+00, 0.132394152788506e+00], [0.059715871789770e+00, 0.470142064105115e+00, 0.132394152788506e+00], [0.470142064105115e+00, 0.059715871789770e+00, 0.132394152788506e+00], [0.101286507323456e+00, 0.101286507323456e+00, 0.125939180544827e+00], [0.797426985353087e+00, 0.101286507323456e+00, 0.125939180544827e+00], [0.101286507323456e+00, 0.797426985353087e+00, 0.125939180544827e+00]], tp_fix=0.5, weight_fix=0.5), 6 : _QP([[0.249286745170910e+00, 0.249286745170910e+00, 0.116786275726379e+00], [0.501426509658179e+00, 0.249286745170910e+00, 0.116786275726379e+00], [0.249286745170910e+00, 0.501426509658179e+00, 0.116786275726379e+00], [0.063089014491502e+00, 0.063089014491502e+00, 0.050844906370207e+00], [0.873821971016996e+00, 0.063089014491502e+00, 0.050844906370207e+00], [0.063089014491502e+00, 0.873821971016996e+00, 0.050844906370207e+00], [0.310352451033784e+00, 0.636502499121399e+00, 0.082851075618374e+00], [0.636502499121399e+00, 0.310352451033784e+00, 0.082851075618374e+00], [0.053145049844817e+00, 0.636502499121399e+00, 0.082851075618374e+00], [0.636502499121399e+00, 0.053145049844817e+00, 0.082851075618374e+00], [0.310352451033784e+00, 0.053145049844817e+00, 0.082851075618374e+00], [0.053145049844817e+00, 0.310352451033784e+00, 0.082851075618374e+00]], tp_fix=0.5, weight_fix=0.5), 7 : _QP([[0.333333333333333e+00, 0.333333333333333e+00,-0.149570044467682e+00], [0.260345966079040e+00, 0.260345966079040e+00, 0.175615257433208e+00], [0.479308067841920e+00, 0.260345966079040e+00, 0.175615257433208e+00], [0.260345966079040e+00, 0.479308067841920e+00, 0.175615257433208e+00], [0.065130102902216e+00, 0.065130102902216e+00, 0.053347235608838e+00], [0.869739794195568e+00, 0.065130102902216e+00, 0.053347235608838e+00], [0.065130102902216e+00, 0.869739794195568e+00, 0.053347235608838e+00], [0.312865496004874e+00, 0.638444188569810e+00, 0.077113760890257e+00], [0.638444188569810e+00, 0.312865496004874e+00, 0.077113760890257e+00], [0.048690315425316e+00, 0.638444188569810e+00, 0.077113760890257e+00], [0.638444188569810e+00, 0.048690315425316e+00, 0.077113760890257e+00], [0.312865496004874e+00, 0.048690315425316e+00, 0.077113760890257e+00], [0.048690315425316e+00, 0.312865496004874e+00, 0.077113760890257e+00]], tp_fix=0.5, weight_fix=0.5), 8 : _QP([[0.333333333333333e+00, 0.333333333333333e+00, 0.144315607677787e+00], [0.459292588292723e+00, 0.459292588292723e+00, 0.095091634267285e+00], [0.081414823414554e+00, 0.459292588292723e+00, 0.095091634267285e+00], [0.459292588292723e+00, 0.081414823414554e+00, 0.095091634267285e+00], [0.170569307751760e+00, 0.170569307751760e+00, 0.103217370534718e+00], [0.658861384496480e+00, 0.170569307751760e+00, 0.103217370534718e+00], [0.170569307751760e+00, 0.658861384496480e+00, 0.103217370534718e+00], [0.050547228317031e+00, 0.050547228317031e+00, 0.032458497623198e+00], [0.898905543365938e+00, 0.050547228317031e+00, 0.032458497623198e+00], [0.050547228317031e+00, 0.898905543365938e+00, 0.032458497623198e+00], [0.263112829634638e+00, 0.728492392955404e+00, 0.027230314174435e+00], [0.728492392955404e+00, 0.263112829634638e+00, 0.027230314174435e+00], [0.008394777409958e+00, 0.728492392955404e+00, 0.027230314174435e+00], [0.728492392955404e+00, 0.008394777409958e+00, 0.027230314174435e+00], [0.263112829634638e+00, 0.008394777409958e+00, 0.027230314174435e+00], [0.008394777409958e+00, 0.263112829634638e+00, 0.027230314174435e+00]], tp_fix=0.5, weight_fix=0.5), 9 : _QP([[0.333333333333333e+00, 0.333333333333333e+00, 0.097135796282799e+00], [0.489682519198738e+00, 0.489682519198738e+00, 0.031334700227139e+00], [0.020634961602525e+00, 0.489682519198738e+00, 0.031334700227139e+00], [0.489682519198738e+00, 0.020634961602525e+00, 0.031334700227139e+00], [0.437089591492937e+00, 0.437089591492937e+00, 0.077827541004774e+00], [0.125820817014127e+00, 0.437089591492937e+00, 0.077827541004774e+00], [0.437089591492937e+00, 0.125820817014127e+00, 0.077827541004774e+00], [0.188203535619033e+00, 0.188203535619033e+00, 0.079647738927210e+00], [0.623592928761935e+00, 0.188203535619033e+00, 0.079647738927210e+00], [0.188203535619033e+00, 0.623592928761935e+00, 0.079647738927210e+00], [0.044729513394453e+00, 0.044729513394453e+00, 0.025577675658698e+00], [0.910540973211095e+00, 0.044729513394453e+00, 0.025577675658698e+00], [0.044729513394453e+00, 0.910540973211095e+00, 0.025577675658698e+00], [0.221962989160766e+00, 0.741198598784498e+00, 0.043283539377289e+00], [0.741198598784498e+00, 0.221962989160766e+00, 0.043283539377289e+00], [0.036838412054736e+00, 0.741198598784498e+00, 0.043283539377289e+00], [0.741198598784498e+00, 0.036838412054736e+00, 0.043283539377289e+00], [0.221962989160766e+00, 0.036838412054736e+00, 0.043283539377289e+00], [0.036838412054736e+00, 0.221962989160766e+00, 0.043283539377289e+00]], tp_fix=0.5, weight_fix=0.5), 10 : _QP([[0.333333333333333e+00, 0.333333333333333e+00, 0.908179903827540e-01], [0.485577633383657e+00, 0.485577633383657e+00, 0.367259577564670e-01], [0.288447332326850e-01, 0.485577633383657e+00, 0.367259577564670e-01], [0.485577633383657e+00, 0.288447332326850e-01, 0.367259577564670e-01], [0.109481575485037e+00, 0.109481575485037e+00, 0.453210594355280e-01], [0.781036849029926e+00, 0.109481575485037e+00, 0.453210594355280e-01], [0.109481575485037e+00, 0.781036849029926e+00, 0.453210594355280e-01], [0.307939838764121e+00, 0.550352941820999e+00, 0.727579168454200e-01], [0.550352941820999e+00, 0.307939838764121e+00, 0.727579168454200e-01], [0.141707219414880e+00, 0.550352941820999e+00, 0.727579168454200e-01], [0.550352941820999e+00, 0.141707219414880e+00, 0.727579168454200e-01], [0.307939838764121e+00, 0.141707219414880e+00, 0.727579168454200e-01], [0.141707219414880e+00, 0.307939838764121e+00, 0.727579168454200e-01], [0.246672560639903e+00, 0.728323904597411e+00, 0.283272425310570e-01], [0.728323904597411e+00, 0.246672560639903e+00, 0.283272425310570e-01], [0.250035347626860e-01, 0.728323904597411e+00, 0.283272425310570e-01], [0.728323904597411e+00, 0.250035347626860e-01, 0.283272425310570e-01], [0.246672560639903e+00, 0.250035347626860e-01, 0.283272425310570e-01], [0.250035347626860e-01, 0.246672560639903e+00, 0.283272425310570e-01], [0.668032510122000e-01, 0.923655933587500e+00, 0.942166696373300e-02], [0.923655933587500e+00, 0.668032510122000e-01, 0.942166696373300e-02], [0.954081540029900e-02, 0.923655933587500e+00, 0.942166696373300e-02], [0.923655933587500e+00, 0.954081540029900e-02, 0.942166696373300e-02], [0.668032510122000e-01, 0.954081540029900e-02, 0.942166696373300e-02], [0.954081540029900e-02, 0.668032510122000e-01, 0.942166696373300e-02]], tp_fix=0.5, weight_fix=0.5), 12 : _QP([[0.488217389773805e+00, 0.488217389773805e+00, 0.257310664404550e-01], [0.235652204523900e-01, 0.488217389773805e+00, 0.257310664404550e-01], [0.488217389773805e+00, 0.235652204523900e-01, 0.257310664404550e-01], [0.439724392294460e+00, 0.439724392294460e+00, 0.436925445380380e-01], [0.120551215411079e+00, 0.439724392294460e+00, 0.436925445380380e-01], [0.439724392294460e+00, 0.120551215411079e+00, 0.436925445380380e-01], [0.271210385012116e+00, 0.271210385012116e+00, 0.628582242178850e-01], [0.457579229975768e+00, 0.271210385012116e+00, 0.628582242178850e-01], [0.271210385012116e+00, 0.457579229975768e+00, 0.628582242178850e-01], [0.127576145541586e+00, 0.127576145541586e+00, 0.347961129307090e-01], [0.744847708916828e+00, 0.127576145541586e+00, 0.347961129307090e-01], [0.127576145541586e+00, 0.744847708916828e+00, 0.347961129307090e-01], [0.213173504532100e-01, 0.213173504532100e-01, 0.616626105155900e-02], [0.957365299093579e+00, 0.213173504532100e-01, 0.616626105155900e-02], [0.213173504532100e-01, 0.957365299093579e+00, 0.616626105155900e-02], [0.275713269685514e+00, 0.608943235779788e+00, 0.403715577663810e-01], [0.608943235779788e+00, 0.275713269685514e+00, 0.403715577663810e-01], [0.115343494534698e+00, 0.608943235779788e+00, 0.403715577663810e-01], [0.608943235779788e+00, 0.115343494534698e+00, 0.403715577663810e-01], [0.275713269685514e+00, 0.115343494534698e+00, 0.403715577663810e-01], [0.115343494534698e+00, 0.275713269685514e+00, 0.403715577663810e-01], [0.281325580989940e+00, 0.695836086787803e+00, 0.223567732023030e-01], [0.695836086787803e+00, 0.281325580989940e+00, 0.223567732023030e-01], [0.228383322222570e-01, 0.695836086787803e+00, 0.223567732023030e-01], [0.695836086787803e+00, 0.228383322222570e-01, 0.223567732023030e-01], [0.281325580989940e+00, 0.228383322222570e-01, 0.223567732023030e-01], [0.228383322222570e-01, 0.281325580989940e+00, 0.223567732023030e-01], [0.116251915907597e+00, 0.858014033544073e+00, 0.173162311086590e-01], [0.858014033544073e+00, 0.116251915907597e+00, 0.173162311086590e-01], [0.257340505483300e-01, 0.858014033544073e+00, 0.173162311086590e-01], [0.858014033544073e+00, 0.257340505483300e-01, 0.173162311086590e-01], [0.116251915907597e+00, 0.257340505483300e-01, 0.173162311086590e-01], [0.257340505483300e-01, 0.116251915907597e+00, 0.173162311086590e-01]], tp_fix=0.5, weight_fix=0.5), 13 : _QP([[0.333333333333333e+00, 0.333333333333333e+00, 0.525209234008020e-01], [0.495048184939705e+00, 0.495048184939705e+00, 0.112801452093300e-01], [0.990363012059100e-02, 0.495048184939705e+00, 0.112801452093300e-01], [0.495048184939705e+00, 0.990363012059100e-02, 0.112801452093300e-01], [0.468716635109574e+00, 0.468716635109574e+00, 0.314235183624540e-01], [0.625667297808520e-01, 0.468716635109574e+00, 0.314235183624540e-01], [0.468716635109574e+00, 0.625667297808520e-01, 0.314235183624540e-01], [0.414521336801277e+00, 0.414521336801277e+00, 0.470725025041940e-01], [0.170957326397447e+00, 0.414521336801277e+00, 0.470725025041940e-01], [0.414521336801277e+00, 0.170957326397447e+00, 0.470725025041940e-01], [0.229399572042831e+00, 0.229399572042831e+00, 0.473635865363550e-01], [0.541200855914337e+00, 0.229399572042831e+00, 0.473635865363550e-01], [0.229399572042831e+00, 0.541200855914337e+00, 0.473635865363550e-01], [0.114424495196330e+00, 0.114424495196330e+00, 0.311675290457940e-01], [0.771151009607340e+00, 0.114424495196330e+00, 0.311675290457940e-01], [0.114424495196330e+00, 0.771151009607340e+00, 0.311675290457940e-01], [0.248113913634590e-01, 0.248113913634590e-01, 0.797577146507400e-02], [0.950377217273082e+00, 0.248113913634590e-01, 0.797577146507400e-02], [0.248113913634590e-01, 0.950377217273082e+00, 0.797577146507400e-02], [0.268794997058761e+00, 0.636351174561660e+00, 0.368484027287320e-01], [0.636351174561660e+00, 0.268794997058761e+00, 0.368484027287320e-01], [0.948538283795790e-01, 0.636351174561660e+00, 0.368484027287320e-01], [0.636351174561660e+00, 0.948538283795790e-01, 0.368484027287320e-01], [0.268794997058761e+00, 0.948538283795790e-01, 0.368484027287320e-01], [0.948538283795790e-01, 0.268794997058761e+00, 0.368484027287320e-01], [0.291730066734288e+00, 0.690169159986905e+00, 0.174014633038220e-01], [0.690169159986905e+00, 0.291730066734288e+00, 0.174014633038220e-01], [0.181007732788070e-01, 0.690169159986905e+00, 0.174014633038220e-01], [0.690169159986905e+00, 0.181007732788070e-01, 0.174014633038220e-01], [0.291730066734288e+00, 0.181007732788070e-01, 0.174014633038220e-01], [0.181007732788070e-01, 0.291730066734288e+00, 0.174014633038220e-01], [0.126357385491669e+00, 0.851409537834241e+00, 0.155217868390450e-01], [0.851409537834241e+00, 0.126357385491669e+00, 0.155217868390450e-01], [0.222330766740900e-01, 0.851409537834241e+00, 0.155217868390450e-01], [0.851409537834241e+00, 0.222330766740900e-01, 0.155217868390450e-01], [0.126357385491669e+00, 0.222330766740900e-01, 0.155217868390450e-01], [0.222330766740900e-01, 0.126357385491669e+00, 0.155217868390450e-01]], tp_fix=0.5, weight_fix=0.5), 14 : _QP([[0.488963910362179e+00, 0.488963910362179e+00, 0.218835813694290e-01], [0.220721792756430e-01, 0.488963910362179e+00, 0.218835813694290e-01], [0.488963910362179e+00, 0.220721792756430e-01, 0.218835813694290e-01], [0.417644719340454e+00, 0.417644719340454e+00, 0.327883535441250e-01], [0.164710561319092e+00, 0.417644719340454e+00, 0.327883535441250e-01], [0.417644719340454e+00, 0.164710561319092e+00, 0.327883535441250e-01], [0.273477528308839e+00, 0.273477528308839e+00, 0.517741045072920e-01], [0.453044943382323e+00, 0.273477528308839e+00, 0.517741045072920e-01], [0.273477528308839e+00, 0.453044943382323e+00, 0.517741045072920e-01], [0.177205532412543e+00, 0.177205532412543e+00, 0.421625887369930e-01], [0.645588935174913e+00, 0.177205532412543e+00, 0.421625887369930e-01], [0.177205532412543e+00, 0.645588935174913e+00, 0.421625887369930e-01], [0.617998830908730e-01, 0.617998830908730e-01, 0.144336996697770e-01], [0.876400233818255e+00, 0.617998830908730e-01, 0.144336996697770e-01], [0.617998830908730e-01, 0.876400233818255e+00, 0.144336996697770e-01], [0.193909612487010e-01, 0.193909612487010e-01, 0.492340360240000e-02], [0.961218077502598e+00, 0.193909612487010e-01, 0.492340360240000e-02], [0.193909612487010e-01, 0.961218077502598e+00, 0.492340360240000e-02], [0.172266687821356e+00, 0.770608554774996e+00, 0.246657532125640e-01], [0.770608554774996e+00, 0.172266687821356e+00, 0.246657532125640e-01], [0.571247574036480e-01, 0.770608554774996e+00, 0.246657532125640e-01], [0.770608554774996e+00, 0.571247574036480e-01, 0.246657532125640e-01], [0.172266687821356e+00, 0.571247574036480e-01, 0.246657532125640e-01], [0.571247574036480e-01, 0.172266687821356e+00, 0.246657532125640e-01], [0.336861459796345e+00, 0.570222290846683e+00, 0.385715107870610e-01], [0.570222290846683e+00, 0.336861459796345e+00, 0.385715107870610e-01], [0.929162493569720e-01, 0.570222290846683e+00, 0.385715107870610e-01], [0.570222290846683e+00, 0.929162493569720e-01, 0.385715107870610e-01], [0.336861459796345e+00, 0.929162493569720e-01, 0.385715107870610e-01], [0.929162493569720e-01, 0.336861459796345e+00, 0.385715107870610e-01], [0.298372882136258e+00, 0.686980167808088e+00, 0.144363081135340e-01], [0.686980167808088e+00, 0.298372882136258e+00, 0.144363081135340e-01], [0.146469500556540e-01, 0.686980167808088e+00, 0.144363081135340e-01], [0.686980167808088e+00, 0.146469500556540e-01, 0.144363081135340e-01], [0.298372882136258e+00, 0.146469500556540e-01, 0.144363081135340e-01], [0.146469500556540e-01, 0.298372882136258e+00, 0.144363081135340e-01], [0.118974497696957e+00, 0.879757171370171e+00, 0.501022883850100e-02], [0.879757171370171e+00, 0.118974497696957e+00, 0.501022883850100e-02], [0.126833093287200e-02, 0.879757171370171e+00, 0.501022883850100e-02], [0.879757171370171e+00, 0.126833093287200e-02, 0.501022883850100e-02], [0.118974497696957e+00, 0.126833093287200e-02, 0.501022883850100e-02], [0.126833093287200e-02, 0.118974497696957e+00, 0.501022883850100e-02]], tp_fix=0.5, weight_fix=0.5), 17 : _QP([[0.333333333333333e+00, 0.333333333333333e+00, 0.334371992908030e-01], [0.497170540556774e+00, 0.497170540556774e+00, 0.509341544050700e-02], [0.565891888645200e-02, 0.497170540556774e+00, 0.509341544050700e-02], [0.497170540556774e+00, 0.565891888645200e-02, 0.509341544050700e-02], [0.482176322624625e+00, 0.482176322624625e+00, 0.146708645276380e-01], [0.356473547507510e-01, 0.482176322624625e+00, 0.146708645276380e-01], [0.482176322624625e+00, 0.356473547507510e-01, 0.146708645276380e-01], [0.450239969020782e+00, 0.450239969020782e+00, 0.243508783536720e-01], [0.995200619584370e-01, 0.450239969020782e+00, 0.243508783536720e-01], [0.450239969020782e+00, 0.995200619584370e-01, 0.243508783536720e-01], [0.400266239377397e+00, 0.400266239377397e+00, 0.311075508689690e-01], [0.199467521245206e+00, 0.400266239377397e+00, 0.311075508689690e-01], [0.400266239377397e+00, 0.199467521245206e+00, 0.311075508689690e-01], [0.252141267970953e+00, 0.252141267970953e+00, 0.312571112186200e-01], [0.495717464058095e+00, 0.252141267970953e+00, 0.312571112186200e-01], [0.252141267970953e+00, 0.495717464058095e+00, 0.312571112186200e-01], [0.162047004658461e+00, 0.162047004658461e+00, 0.248156543396650e-01], [0.675905990683077e+00, 0.162047004658461e+00, 0.248156543396650e-01], [0.162047004658461e+00, 0.675905990683077e+00, 0.248156543396650e-01], [0.758758822607460e-01, 0.758758822607460e-01, 0.140560730705570e-01], [0.848248235478508e+00, 0.758758822607460e-01, 0.140560730705570e-01], [0.758758822607460e-01, 0.848248235478508e+00, 0.140560730705570e-01], [0.156547269678220e-01, 0.156547269678220e-01, 0.319467617377900e-02], [0.968690546064356e+00, 0.156547269678220e-01, 0.319467617377900e-02], [0.156547269678220e-01, 0.968690546064356e+00, 0.319467617377900e-02], [0.334319867363658e+00, 0.655493203809423e+00, 0.811965531899300e-02], [0.655493203809423e+00, 0.334319867363658e+00, 0.811965531899300e-02], [0.101869288269190e-01, 0.655493203809423e+00, 0.811965531899300e-02], [0.655493203809423e+00, 0.101869288269190e-01, 0.811965531899300e-02], [0.334319867363658e+00, 0.101869288269190e-01, 0.811965531899300e-02], [0.101869288269190e-01, 0.334319867363658e+00, 0.811965531899300e-02], [0.292221537796944e+00, 0.572337590532020e+00, 0.268057422831630e-01], [0.572337590532020e+00, 0.292221537796944e+00, 0.268057422831630e-01], [0.135440871671036e+00, 0.572337590532020e+00, 0.268057422831630e-01], [0.572337590532020e+00, 0.135440871671036e+00, 0.268057422831630e-01], [0.292221537796944e+00, 0.135440871671036e+00, 0.268057422831630e-01], [0.135440871671036e+00, 0.292221537796944e+00, 0.268057422831630e-01], [0.319574885423190e+00, 0.626001190286228e+00, 0.184599932108220e-01], [0.626001190286228e+00, 0.319574885423190e+00, 0.184599932108220e-01], [0.544239242905830e-01, 0.626001190286228e+00, 0.184599932108220e-01], [0.626001190286228e+00, 0.544239242905830e-01, 0.184599932108220e-01], [0.319574885423190e+00, 0.544239242905830e-01, 0.184599932108220e-01], [0.544239242905830e-01, 0.319574885423190e+00, 0.184599932108220e-01], [0.190704224192292e+00, 0.796427214974071e+00, 0.847686853432800e-02], [0.796427214974071e+00, 0.190704224192292e+00, 0.847686853432800e-02], [0.128685608336370e-01, 0.796427214974071e+00, 0.847686853432800e-02], [0.796427214974071e+00, 0.128685608336370e-01, 0.847686853432800e-02], [0.190704224192292e+00, 0.128685608336370e-01, 0.847686853432800e-02], [0.128685608336370e-01, 0.190704224192292e+00, 0.847686853432800e-02], [0.180483211648746e+00, 0.752351005937729e+00, 0.182927967700250e-01], [0.752351005937729e+00, 0.180483211648746e+00, 0.182927967700250e-01], [0.671657824135240e-01, 0.752351005937729e+00, 0.182927967700250e-01], [0.752351005937729e+00, 0.671657824135240e-01, 0.182927967700250e-01], [0.180483211648746e+00, 0.671657824135240e-01, 0.182927967700250e-01], [0.671657824135240e-01, 0.180483211648746e+00, 0.182927967700250e-01], [0.807113136795640e-01, 0.904625504095608e+00, 0.666563200416500e-02], [0.904625504095608e+00, 0.807113136795640e-01, 0.666563200416500e-02], [0.146631822248280e-01, 0.904625504095608e+00, 0.666563200416500e-02], [0.904625504095608e+00, 0.146631822248280e-01, 0.666563200416500e-02], [0.807113136795640e-01, 0.146631822248280e-01, 0.666563200416500e-02], [0.146631822248280e-01, 0.807113136795640e-01, 0.666563200416500e-02]], tp_fix=0.5, weight_fix=0.5), 19 : _QP([[0.333333333333333e+00, 0.333333333333333e+00, 0.329063313889190e-01], [0.489609987073006e+00, 0.489609987073006e+00, 0.103307318912720e-01], [0.207800258539870e-01, 0.489609987073006e+00, 0.103307318912720e-01], [0.489609987073006e+00, 0.207800258539870e-01, 0.103307318912720e-01], [0.454536892697893e+00, 0.454536892697893e+00, 0.223872472630160e-01], [0.909262146042150e-01, 0.454536892697893e+00, 0.223872472630160e-01], [0.454536892697893e+00, 0.909262146042150e-01, 0.223872472630160e-01], [0.401416680649431e+00, 0.401416680649431e+00, 0.302661258694680e-01], [0.197166638701138e+00, 0.401416680649431e+00, 0.302661258694680e-01], [0.401416680649431e+00, 0.197166638701138e+00, 0.302661258694680e-01], [0.255551654403098e+00, 0.255551654403098e+00, 0.304909678021980e-01], [0.488896691193805e+00, 0.255551654403098e+00, 0.304909678021980e-01], [0.255551654403098e+00, 0.488896691193805e+00, 0.304909678021980e-01], [0.177077942152130e+00, 0.177077942152130e+00, 0.241592127416410e-01], [0.645844115695741e+00, 0.177077942152130e+00, 0.241592127416410e-01], [0.177077942152130e+00, 0.645844115695741e+00, 0.241592127416410e-01], [0.110061053227952e+00, 0.110061053227952e+00, 0.160508035868010e-01], [0.779877893544096e+00, 0.110061053227952e+00, 0.160508035868010e-01], [0.110061053227952e+00, 0.779877893544096e+00, 0.160508035868010e-01], [0.555286242518400e-01, 0.555286242518400e-01, 0.808458026178400e-02], [0.888942751496321e+00, 0.555286242518400e-01, 0.808458026178400e-02], [0.555286242518400e-01, 0.888942751496321e+00, 0.808458026178400e-02], [0.126218637772290e-01, 0.126218637772290e-01, 0.207936202748500e-02], [0.974756272445543e+00, 0.126218637772290e-01, 0.207936202748500e-02], [0.126218637772290e-01, 0.974756272445543e+00, 0.207936202748500e-02], [0.395754787356943e+00, 0.600633794794645e+00, 0.388487690498100e-02], [0.600633794794645e+00, 0.395754787356943e+00, 0.388487690498100e-02], [0.361141784841200e-02, 0.600633794794645e+00, 0.388487690498100e-02], [0.600633794794645e+00, 0.361141784841200e-02, 0.388487690498100e-02], [0.395754787356943e+00, 0.361141784841200e-02, 0.388487690498100e-02], [0.361141784841200e-02, 0.395754787356943e+00, 0.388487690498100e-02], [0.307929983880436e+00, 0.557603261588784e+00, 0.255741606120220e-01], [0.557603261588784e+00, 0.307929983880436e+00, 0.255741606120220e-01], [0.134466754530780e+00, 0.557603261588784e+00, 0.255741606120220e-01], [0.557603261588784e+00, 0.134466754530780e+00, 0.255741606120220e-01], [0.307929983880436e+00, 0.134466754530780e+00, 0.255741606120220e-01], [0.134466754530780e+00, 0.307929983880436e+00, 0.255741606120220e-01], [0.264566948406520e+00, 0.720987025817365e+00, 0.888090357333800e-02], [0.720987025817365e+00, 0.264566948406520e+00, 0.888090357333800e-02], [0.144460257761150e-01, 0.720987025817365e+00, 0.888090357333800e-02], [0.720987025817365e+00, 0.144460257761150e-01, 0.888090357333800e-02], [0.264566948406520e+00, 0.144460257761150e-01, 0.888090357333800e-02], [0.144460257761150e-01, 0.264566948406520e+00, 0.888090357333800e-02], [0.358539352205951e+00, 0.594527068955871e+00, 0.161245467617310e-01], [0.594527068955871e+00, 0.358539352205951e+00, 0.161245467617310e-01], [0.469335788381780e-01, 0.594527068955871e+00, 0.161245467617310e-01], [0.594527068955871e+00, 0.469335788381780e-01, 0.161245467617310e-01], [0.358539352205951e+00, 0.469335788381780e-01, 0.161245467617310e-01], [0.469335788381780e-01, 0.358539352205951e+00, 0.161245467617310e-01], [0.157807405968595e+00, 0.839331473680839e+00, 0.249194181749100e-02], [0.839331473680839e+00, 0.157807405968595e+00, 0.249194181749100e-02], [0.286112035056700e-02, 0.839331473680839e+00, 0.249194181749100e-02], [0.839331473680839e+00, 0.286112035056700e-02, 0.249194181749100e-02], [0.157807405968595e+00, 0.286112035056700e-02, 0.249194181749100e-02], [0.286112035056700e-02, 0.157807405968595e+00, 0.249194181749100e-02], [0.750505969759110e-01, 0.701087978926173e+00, 0.182428401189510e-01], [0.701087978926173e+00, 0.750505969759110e-01, 0.182428401189510e-01], [0.223861424097916e+00, 0.701087978926173e+00, 0.182428401189510e-01], [0.701087978926173e+00, 0.223861424097916e+00, 0.182428401189510e-01], [0.750505969759110e-01, 0.223861424097916e+00, 0.182428401189510e-01], [0.223861424097916e+00, 0.750505969759110e-01, 0.182428401189510e-01], [0.142421601113383e+00, 0.822931324069857e+00, 0.102585637361990e-01], [0.822931324069857e+00, 0.142421601113383e+00, 0.102585637361990e-01], [0.346470748167600e-01, 0.822931324069857e+00, 0.102585637361990e-01], [0.822931324069857e+00, 0.346470748167600e-01, 0.102585637361990e-01], [0.142421601113383e+00, 0.346470748167600e-01, 0.102585637361990e-01], [0.346470748167600e-01, 0.142421601113383e+00, 0.102585637361990e-01], [0.654946280829380e-01, 0.924344252620784e+00, 0.379992885530200e-02], [0.924344252620784e+00, 0.654946280829380e-01, 0.379992885530200e-02], [0.101611192962780e-01, 0.924344252620784e+00, 0.379992885530200e-02], [0.924344252620784e+00, 0.101611192962780e-01, 0.379992885530200e-02], [0.654946280829380e-01, 0.101611192962780e-01, 0.379992885530200e-02], [0.101611192962780e-01, 0.654946280829380e-01, 0.379992885530200e-02]], tp_fix=0.5, weight_fix=0.5), }, '2_4' : { 2 : _QP([[ nm.sqrt(2.0/3.0), 0.0 , 4.0/3.0], [-1/nm.sqrt(6) , 1/nm.sqrt(2), 4.0/3.0], [-1/nm.sqrt(6) ,-1/nm.sqrt(2), 4.0/3.0]], bounds=(-1.0, 1.0)), 3 : _QP([[-1/nm.sqrt(3),-1/nm.sqrt(3), 1.0], [ 1/nm.sqrt(3),-1/nm.sqrt(3), 1.0], [ 1/nm.sqrt(3), 1/nm.sqrt(3), 1.0], [-1/nm.sqrt(3), 1/nm.sqrt(3), 1.0]], bounds=(-1.0, 1.0)), 5 : _QP([[ nm.sqrt(7.0/15.0), 0.0 , 0.816326530612245], [-nm.sqrt(7.0/15.0), 0.0 , 0.816326530612245], [ 0.0 , nm.sqrt(7.0/15.0), 0.816326530612245], [ 0.0 ,-nm.sqrt(7.0/15.0), 0.816326530612245], [ 0.881917103688197, 0.881917103688197, 0.183673469387755], [ 0.881917103688197,-0.881917103688197, 0.183673469387755], [-0.881917103688197, 0.881917103688197, 0.183673469387755], [-0.881917103688197,-0.881917103688197, 0.183673469387755]], bounds=(-1.0, 1.0)), }, '3_4' : { 1 : _QP([[ 1.0/4.0, 1.0/4.0, 1.0/4.0, 1.0/6.0]], tp_fix=1.0/6.0), 2 : _QP([[ (5-nm.sqrt(5))/20 , (5-nm.sqrt(5))/20 , (5-nm.sqrt(5))/20 , 1.0/24.0], [ (5-nm.sqrt(5))/20 , (5-nm.sqrt(5))/20 , (5+3*nm.sqrt(5))/20, 1.0/24.0], [ (5-nm.sqrt(5))/20 , (5+3*nm.sqrt(5))/20, (5-nm.sqrt(5))/20 , 1.0/24.0], [ (5+3*nm.sqrt(5))/20, (5-nm.sqrt(5))/20 , (5-nm.sqrt(5))/20 , 1.0/24.0]], tp_fix=1.0/6.0), 3 : _QP([[ 1.0/4.0, 1.0/4.0, 1.0/4.0,-2.0/15.0], [ 1.0/6.0, 1.0/6.0, 1.0/6.0, 3.0/40.0], [ 1.0/6.0, 1.0/6.0, 1.0/2.0, 3.0/40.0], [ 1.0/6.0, 1.0/2.0, 1.0/6.0, 3.0/40.0], [ 1.0/2.0, 1.0/6.0, 1.0/6.0, 3.0/40.0]], tp_fix=1.0/6.0), 4 : _QP([[-0.5000000000000000, -0.5000000000000000, -0.5000000000000000, -0.1052444444444440], [-0.8571428571428570, -0.8571428571428570, -0.8571428571428570, 0.0609777777777780], [-0.8571428571428570, -0.8571428571428570, 0.5714285714285710, 0.0609777777777780], [-0.8571428571428570, 0.5714285714285710, -0.8571428571428570, 0.0609777777777780], [ 0.5714285714285710, -0.8571428571428570, -0.8571428571428570, 0.0609777777777780], [-0.2011928476664020, -0.2011928476664020, -0.7988071523335980, 0.1991111111111110], [-0.2011928476664020, -0.7988071523335980, -0.2011928476664020, 0.1991111111111110], [-0.7988071523335980, -0.2011928476664020, -0.2011928476664020, 0.1991111111111110], [-0.2011928476664020, -0.7988071523335980, -0.7988071523335980, 0.1991111111111110], [-0.7988071523335980, -0.2011928476664020, -0.7988071523335980, 0.1991111111111110], [-0.7988071523335980, -0.7988071523335980, -0.2011928476664020, 0.1991111111111110]], bounds=(-1.0, 1.0), tp_fix=1.0/6.0), 6 : _QP([[-0.5707942574816960, -0.5707942574816960, -0.5707942574816960, 0.0532303336775570], [-0.2876172275549120, -0.5707942574816960, -0.5707942574816960, 0.0532303336775570], [-0.5707942574816960, -0.2876172275549120, -0.5707942574816960, 0.0532303336775570], [-0.5707942574816960, -0.5707942574816960, -0.2876172275549120, 0.0532303336775570], [-0.9186520829307770, -0.9186520829307770, -0.9186520829307770, 0.0134362814070940], [0.7559562487923320, -0.9186520829307770, -0.9186520829307770, 0.0134362814070940], [-0.9186520829307770, 0.7559562487923320, -0.9186520829307770, 0.0134362814070940], [-0.9186520829307770, -0.9186520829307770, 0.7559562487923320, 0.0134362814070940], [-0.3553242197154490, -0.3553242197154490, -0.3553242197154490, 0.0738095753915400], [-0.9340273408536530, -0.3553242197154490, -0.3553242197154490, 0.0738095753915400], [-0.3553242197154490, -0.9340273408536530, -0.3553242197154490, 0.0738095753915400], [-0.3553242197154490, -0.3553242197154490, -0.9340273408536530, 0.0738095753915400], [-0.8726779962499650, -0.8726779962499650, -0.4606553370833680, 0.0642857142857140], [-0.8726779962499650, -0.4606553370833680, -0.8726779962499650, 0.0642857142857140], [-0.8726779962499650, -0.8726779962499650, 0.2060113295832980, 0.0642857142857140], [-0.8726779962499650, 0.2060113295832980, -0.8726779962499650, 0.0642857142857140], [-0.8726779962499650, -0.4606553370833680, 0.2060113295832980, 0.0642857142857140], [-0.8726779962499650, 0.2060113295832980, -0.4606553370833680, 0.0642857142857140], [-0.4606553370833680, -0.8726779962499650, -0.8726779962499650, 0.0642857142857140], [-0.4606553370833680, -0.8726779962499650, 0.2060113295832980, 0.0642857142857140], [-0.4606553370833680, 0.2060113295832980, -0.8726779962499650, 0.0642857142857140], [0.2060113295832980, -0.8726779962499650, -0.4606553370833680, 0.0642857142857140], [0.2060113295832980, -0.8726779962499650, -0.8726779962499650, 0.0642857142857140], [0.2060113295832980, -0.4606553370833680, -0.8726779962499650, 0.0642857142857140]], bounds=(-1.0, 1.0), tp_fix=1.0/6.0), }, '3_6': { 3: _QP([[1.0/3.0, 1.0/3.0, 0.5, 0.5]], tp_fix=1/2), 6: _QP([[1.0/6.0, 1.0/6.0, 0.2113248654051871, 1.0/12.0], [2.0/3.0, 1.0/6.0, 0.2113248654051871, 1.0/12.0], [1.0/6.0, 2.0/3.0, 0.2113248654051871, 1.0/12.0], [1.0/6.0, 1.0/6.0, 0.7886751345948129, 1.0/12.0], [2.0/3.0, 1.0/6.0, 0.7886751345948129, 1.0/12.0], [1.0/6.0, 2.0/3.0, 0.7886751345948129, 1.0/12.0]], tp_fix=1/2), 9: _QP([[0.1081030181680702, 0.4459484909159649, 0.1127016653792583, 0.0310252207886127], [0.4459484909159649, 0.1081030181680702, 0.1127016653792583, 0.0310252207886127], [0.4459484909159649, 0.4459484909159649, 0.1127016653792583, 0.0310252207886127], [0.8168475729804585, 0.0915762135097707, 0.1127016653792583, 0.0152710755076836], [0.0915762135097707, 0.8168475729804585, 0.1127016653792583, 0.0152710755076836], [0.0915762135097707, 0.0915762135097707, 0.1127016653792583, 0.0152710755076836], [0.1081030181680702, 0.4459484909159649, 0.5000000000000000, 0.0496403532617803], [0.4459484909159649, 0.1081030181680702, 0.5000000000000000, 0.0496403532617803], [0.4459484909159649, 0.4459484909159649, 0.5000000000000000, 0.0496403532617803], [0.8168475729804585, 0.0915762135097707, 0.5000000000000000, 0.0244337208122937], [0.0915762135097707, 0.8168475729804585, 0.5000000000000000, 0.0244337208122937], [0.0915762135097707, 0.0915762135097707, 0.5000000000000000, 0.0244337208122937], [0.1081030181680702, 0.4459484909159649, 0.8872983346207417, 0.0310252207886127], [0.4459484909159649, 0.1081030181680702, 0.8872983346207417, 0.0310252207886127], [0.4459484909159649, 0.4459484909159649, 0.8872983346207417, 0.0310252207886127], [0.8168475729804585, 0.0915762135097707, 0.8872983346207417, 0.0152710755076836], [0.0915762135097707, 0.8168475729804585, 0.8872983346207417, 0.0152710755076836], [0.0915762135097707, 0.0915762135097707, 0.8872983346207417, 0.0152710755076836]], tp_fix=1/2), }, '3_8' : { 2 : _QP([[ 0.0 , nm.sqrt(2.0/3.0),-1/nm.sqrt(3), 2.0], [ 0.0 ,-nm.sqrt(2.0/3.0),-1/nm.sqrt(3), 2.0], [ nm.sqrt(2.0/3.0), 0.0 , 1/nm.sqrt(3), 2.0], [-nm.sqrt(2.0/3.0), 0.0 , 1/nm.sqrt(3), 2.0]], bounds=(-1.0, 1.0)), 3 : _QP([[-1.0, 0.0, 0.0, 4.0/3.0], [ 1.0, 0.0, 0.0, 4.0/3.0], [ 0.0,-1.0, 0.0, 4.0/3.0], [ 0.0, 1.0, 0.0, 4.0/3.0], [ 0.0, 0.0,-1.0, 4.0/3.0], [ 0.0, 0.0, 1.0, 4.0/3.0]], bounds=(-1.0, 1.0)), 5 : _QP([[-nm.sqrt(19.0/30.0), 0.0 , 0.0 , 320.0/361.0], [ nm.sqrt(19.0/30.0), 0.0 , 0.0 , 320.0/361.0], [ 0.0 ,-nm.sqrt(19.0/30.0), 0.0 , 320.0/361.0], [ 0.0 , nm.sqrt(19.0/30.0), 0.0 , 320.0/361.0], [ 0.0 , 0.0 ,-nm.sqrt(19.0/30.0), 320.0/361.0], [ 0.0 , 0.0 , nm.sqrt(19.0/30.0), 320.0/361.0], [ nm.sqrt(19.0/33.0), nm.sqrt(19.0/33.0), nm.sqrt(19.0/33.0), 121.0/361.0], [ nm.sqrt(19.0/33.0), nm.sqrt(19.0/33.0),-nm.sqrt(19.0/33.0), 121.0/361.0], [ nm.sqrt(19.0/33.0),-nm.sqrt(19.0/33.0), nm.sqrt(19.0/33.0), 121.0/361.0], [ nm.sqrt(19.0/33.0),-nm.sqrt(19.0/33.0),-nm.sqrt(19.0/33.0), 121.0/361.0], [-nm.sqrt(19.0/33.0), nm.sqrt(19.0/33.0), nm.sqrt(19.0/33.0), 121.0/361.0], [-nm.sqrt(19.0/33.0), nm.sqrt(19.0/33.0),-nm.sqrt(19.0/33.0), 121.0/361.0], [-nm.sqrt(19.0/33.0),-nm.sqrt(19.0/33.0), nm.sqrt(19.0/33.0), 121.0/361.0], [-nm.sqrt(19.0/33.0),-nm.sqrt(19.0/33.0),-nm.sqrt(19.0/33.0), 121.0/361.0]], bounds=(-1.0, 1.0)), }, } del _QP def _get_max_orders(): max_orders = {} for key, table in quadrature_tables.items(): orders = list(table.keys()) max_orders[key] = max(orders) return max_orders max_orders = _get_max_orders()