"""
`Node`, `edge` and `face` polynomials in the 2d reference domain
:math:`\Pi_{\mathrm{ref}}=[-1,1]^2`.
.. _grid2d:
======
grid2d
======
Let ``A=(a1, a2)``, ``B=(b1,b2,b3)``. Then
:ref:`grid2d` (``A``, ``B``) refers to a sequence of coordinates,
``(a1, b1)``, ``(a2, b1)``, ``(a1, b2)``, ``(a2, b2)``,
``(a1, b3)``, ``(a2, b3)``. Namely, we
first do a ``meshgrid`` (``A``, ``B``), then put the coordinates
into a sequence one by one picking along ``A`` direction firstly, ``B``
direction secondly and finally ``C`` direction. Also see :func:`grid2d`.
⭕ To access the source code, click on the [source] button at the right
side or click on
:download:`[mimetic_basis_polynomials_2d.py]</contents/LIBRARY/ptc/mathischeap_ptc/mimetic_basis_polynomials_2d.py>`.
Dependence may exist. In case of error, check import and install required
packages or download required scripts. © mathischeap.com
"""
import numpy as np
from quadrature import Lobatto
from Lagrange_and_edge_polynomials import Lagrange_polynomials as Lp1
from Lagrange_and_edge_polynomials import edge_polynomials as ep1
[docs]
class MimeticBasisPolynomials2D(object):
"""A wrapper of basis node, edge, face polynomials in the
2d reference domain :math:`\Pi_{\mathrm{ref}}=[-1,1]^2`.
:param nodes_rho: The nodes on which the 1D mimetic polynomials are
built along the first axis (:math:`\\varrho`).
:param nodes_tau: The nodes on which the 1D mimetic polynomials are
built along the second axis (:math:`\\tau`).
:type nodes_rho: 1d np.array
:type nodes_tau: 1d np.array
:example:
>>> bf = MimeticBasisPolynomials2D('Lobatto-3', 'Lobatto-3')
>>> bf.degree # N = nodes_rho = nodes_tau = 3
[3, 3]
>>> rho = np.linspace(-1, 1, 5)
>>> tau = np.linspace(-1, 1, 6)
>>> NP = bf.node_polynomials(rho, tau)
>>> NP.shape # 4^2=16 node polynomials evaluated at 5*6=30 points
(16, 30)
>>> EP_rho, EP_tau = bf.edge_polynomials(rho, tau)
>>> EP_rho.shape # 3*4=12 edge polynomials
(12, 30)
>>> FP = bf.face_polynomials(rho, tau)
>>> FP.shape # 3*3=9 face polynomials
(9, 30)
"""
def __init__(self, nodes_rho, nodes_tau):
self.nodes = [self._parse_nodes_(nodes_rho),
self._parse_nodes_(nodes_tau)]
self.degree = [len(self.nodes[i])-1 for i in range(2)]
@staticmethod
def _parse_nodes_(nodes):
"""Parse the nodes. You can customize your own input shortcut
through this method.
"""
if isinstance(nodes, str):
if nodes[:7] == 'Lobatto': # for example, nodes = "Lobatto-5"
_, nodes_degree = nodes.split('-')
nodes_degree = int(nodes_degree)
nodes, _ = Lobatto(nodes_degree)
else:
pass
return nodes
[docs]
def node_polynomials(self, rho, tau):
"""Evaluate the node polynomials at ``grid2d(rho, tau)``.
:param rho: :ref:`grid2d` (``rho``, ``tau``) is the grid
to evaluate the polynomials.
:param tau: :ref:`grid2d` (``rho``, ``tau``)is the grid
to evaluate the polynomials.
:type rho: 1d np.array
:type tau: 1d np.array
"""
nodes_rho, nodes_tau = self.nodes
np_rho = Lp1(nodes_rho, rho)
np_tau = Lp1(nodes_tau, tau)
return np.kron(np_tau, np_rho)
[docs]
def edge_polynomials(self, rho, tau):
"""Evaluate the edge polynomials at ``grid2d(rho, tau)``.
:param rho: :ref:`grid2d` (``rho``, ``tau``) is the grid
to evaluate the polynomials.
:param tau: :ref:`grid2d` (``rho``, ``tau``)is the grid
to evaluate the polynomials.
:type rho: 1d np.array
:type tau: 1d np.array
"""
nodes_rho, nodes_tau = self.nodes
np_rho = Lp1(nodes_rho, rho)
np_tau = Lp1(nodes_tau, tau)
ep_rho = ep1(nodes_rho, rho)
ep_tau = ep1(nodes_tau, tau)
return np.kron(np_tau, ep_rho), np.kron(ep_tau, np_rho)
[docs]
def face_polynomials(self, rho, tau):
"""Evaluate the face polynomials at ``grid2d(rho, tau)``.
:param rho: :ref:`grid2d` (``rho``, ``tau``) is the grid
to evaluate the polynomials.
:param tau: :ref:`grid2d` (``rho``, ``tau``)is the grid
to evaluate the polynomials.
:type rho: 1d np.array
:type tau: 1d np.array
"""
nodes_rho, nodes_tau = self.nodes
ep_rho = ep1(nodes_rho, rho)
ep_tau = ep1(nodes_tau, tau)
return np.kron(ep_tau, ep_rho)
[docs]
def grid2d(A, B):
"""The function to compute the grid of two sets of nodes.
:param A: The first (along :math:`\\varrho`) set of nodes.
:param B: The second (along :math:`\\tau`) set of nodes.
:type A: 1d data object
:type B: 1d data object
:return: A tuple of three outputs:
#. (1d np.array) The grided :math:`\\varrho` coordinates.
#. (1d np.array) The grided :math:`\\tau` coordinates.
:example:
>>> A = np.array([1, 2])
>>> B = np.array([3, 4, 5])
>>> x, y = grid2d(A, B)
>>> D = np.vstack((x,y)).T
>>> D
array([[1, 3],
[2, 3],
[1, 4],
[2, 4],
[1, 5],
[2, 5]])
"""
x, y = np.meshgrid(A, B, indexing='ij')
x = x.ravel('F')
y = y.ravel('F')
return x, y
if __name__ == '__main__':
import doctest
doctest.testmod()
