Source code for Muscat.FE.Spaces.IPSpaces

# -*- coding: utf-8 -*-
#
# This file is subject to the terms and conditions defined in
# file 'LICENSE.txt', which is part of this source code package.
#
from weakref import WeakKeyDictionary
import numpy as np

import Muscat.Containers.ElementsDescription as ED
from Muscat.FE.IntegrationRules import MeshQuadrature, ElementaryQuadrature
from Muscat.FE.Spaces.FESpaces import FESpace

_cache_:WeakKeyDictionary = WeakKeyDictionary()

def GenerateSpaceForIntegrationPointInterpolation(integrationRule:MeshQuadrature):
    if integrationRule in _cache_:
        return _cache_[integrationRule]

    res = FESpace("FESpaceFor"+integrationRule.name)
    res.hasSharedShapeFunctions = False

    from Muscat.FE.Spaces.SymSpace import SymSpaceBase
    class IntegrationPointSpace(SymSpaceBase):
        def __init__(self,elementIR:ElementaryQuadrature) -> None:
            super().__init__()

            self.geoSupport = ED.geoSupport[elementIR.elementType]
            from sympy.matrices import Matrix
            from sympy import DiracDelta
            integrationPoints  = elementIR.points
            coord  = [self.xi, self.eta, self.phi]

            self.symN = Matrix([ np.prod([DiracDelta(c-x) for x,c in zip(y,coord)]) for y in integrationPoints])
            self.posN = np.array(integrationPoints)
            self.dofAttachments = [("IP",i,None) for i in range(len(integrationPoints)) ]

    for eIR in integrationRule.values():
            res[eIR.elementType] = IntegrationPointSpace(eIR)

    _cache_[integrationRule] = res

    return res


def CheckIntegrity(GUI:bool=False):
    from Muscat.FE.IntegrationRules import LagrangeP1Quadrature

    a = GenerateSpaceForIntegrationPointInterpolation(LagrangeP1Quadrature)
    #print(a)
    a[ED.Triangle_3].GetSpaceEvaluatedAt(LagrangeP1Quadrature[ED.Triangle_3].points,LagrangeP1Quadrature[ED.Triangle_3].weights)
    return "ok"

if __name__ == '__main__':
    print(CheckIntegrity(True))# pragma: no cover