# Source code for mass_matrices_trace

"""

⭕ To access the source code, click on the [source] button at the right
side or click on
:download:[mass_matrices_trace.py]</contents/LIBRARY/ptc/mathischeap_ptc/mass_matrices_trace.py>.
Dependence may exist. In case of error, check import and install required

"""

import numpy as np
from mimetic_basis_polynomials_2d import grid2d, MimeticBasisPolynomials2D
from scipy import sparse as spspa
from coordinate_transformation_surface import extract_surface_coordinate_transformations_of

[docs]class MassMatricesTrace(object):
""" The mass matrices of trace spaces
:math:\\text{TN}_{N}(\\partial\\Omega),
:math:\\text{TE}_{N-1}(\\partial\\Omega) and
:math:\\text{TF}_{N-1}(\\partial\\Omega).

:param bf: A :class:MimeticBasisPolynomials
(not :class:MimeticBasisPolynomials2D) instance.
:type bf: MimeticBasisPolynomials
:param ct: A CoordinateTransformation instance that represents the
mapping :math:\\Phi:\\Omega_{\mathrm{ref}}\\to\\Omega.
:type ct: CoordinateTransformation
:param quad_degree: (default: None) The degree used for the
numerical integral. It should be a list or tuple of three
positive integers. If it is None, a suitable degree will be
obtained from bf.

:example:

>>> from coordinate_transformation import CoordinateTransformation
>>> from coordinate_transformation import Phi, d_Phi
>>> from mimetic_basis_polynomials import MimeticBasisPolynomials
>>> ct = CoordinateTransformation(Phi, d_Phi)
>>> bf = MimeticBasisPolynomials('Lobatto-4', 'Lobatto-3', 'Lobatto-2')
>>> MMT = MassMatricesTrace(bf, ct)
>>> MMT.TF # doctest: +ELLIPSIS
<52x52 sparse matrix of type '<class 'numpy.float64'>'...

"""
assert bf.__class__.__name__ == 'MimeticBasisPolynomials'
assert ct.__class__.__name__ == 'CoordinateTransformation'
self.bf = bf
self.ct = ct
bf_N = bf.degree
quad_degree = [bf_N[i]+2 for i in range(3)]
else:
pass
'WE': np.kron(qw_2, qw_0),
'BF': np.kron(qw_1, qw_0)}
'WE': [qn_0, qn_2],
'BF': [qn_0, qn_1],}
self._rho_tau_ = {'NS':grid2d(qn_1, qn_2),
'WE':grid2d(qn_0, qn_2),
'BF':grid2d(qn_0, qn_1)}
nodes = self.bf.nodes
self.bft = {'NS':MimeticBasisPolynomials2D(nodes[1], nodes[2]),
'WE':MimeticBasisPolynomials2D(nodes[0], nodes[2]),
'BF':MimeticBasisPolynomials2D(nodes[0], nodes[1])}

@property
def TN(self):
"""The mass matrix :math:\\mathbb{M}_{\\text{N}}.

:return: Return a csc_matrix representing the mass matrix.
"""
raise NotImplementedError("Could you code it?")

@property
def TE(self):
"""The mass matrix :math:\\mathbb{M}_{\\text{E}}.

:return: Return a csc_matrix representing the mass matrix.
"""
raise NotImplementedError("Could you code it?")

@property
def TF(self):
"""The mass matrix :math:\\mathbb{M}_{\\text{F}}.

:return: Return a csc_matrix representing the mass matrix.
"""
N, S, W, E, B, F = extract_surface_coordinate_transformations_of(self.ct)

bf = self.bft['NS'].face_polynomials(rho, tau)
# North side
g = N.metric(*self._rho_tau_['NS'])
M_N = np.einsum('im, jm, m -> ij',
optimize='greedy')
# South side
g = S.metric(*self._rho_tau_['NS'])
M_S = np.einsum('im, jm, m -> ij',
optimize='greedy')

bf = self.bft['WE'].face_polynomials(rho, tau)
# West side
g = W.metric(*self._rho_tau_['WE'])
M_W = np.einsum('im, jm, m -> ij',
optimize='greedy')
# East side
g = E.metric(*self._rho_tau_['WE'])
M_E = np.einsum('im, jm, m -> ij',
optimize='greedy')

bf = self.bft['BF'].face_polynomials(rho, tau)
# Back side
g = B.metric(*self._rho_tau_['BF'])
M_B = np.einsum('im, jm, m -> ij',
optimize='greedy')
# Front side
g = F.metric(*self._rho_tau_['BF'])
M_F = np.einsum('im, jm, m -> ij',
optimize='greedy')

M_N = spspa.csc_matrix(M_N)
M_S = spspa.csc_matrix(M_S)
M_W = spspa.csc_matrix(M_W)
M_E = spspa.csc_matrix(M_E)
M_B = spspa.csc_matrix(M_B)
M_F = spspa.csc_matrix(M_F)

M = spspa.bmat([(M_N , None, None, None, None, None),
(None, M_S , None, None, None, None),
(None, None, M_W , None, None, None),
(None, None, None, M_E , None, None),
(None, None, None, None, M_B , None),
(None, None, None, None, None, M_F )], format='csc')

return M

if __name__ == '__main__':
import doctest
doctest.testmod()

from coordinate_transformation import Phi, d_Phi
from coordinate_transformation import CoordinateTransformation
from mimetic_basis_polynomials import MimeticBasisPolynomials
ct = CoordinateTransformation(Phi, d_Phi)
bf = MimeticBasisPolynomials('Lobatto-4', 'Lobatto-3', 'Lobatto-2')
MMT = MassMatricesTrace(bf, ct)
M2 = MMT.TF

print(M2)