Source code for FIAT.bell

# Copyright (C) 2018 Robert C. Kirby
#
# This file is part of FIAT (https://www.fenicsproject.org)
#
# SPDX-License-Identifier:    LGPL-3.0-or-later

# This is not quite Bell, but is 21-dofs and includes 3 extra constraint
# functionals.  The first 18 basis functions are the reference element
# bfs, but the extra three are used in the transformation theory.

from FIAT import finite_element, polynomial_set, dual_set, functional
from FIAT.reference_element import TRIANGLE, ufc_simplex


[docs] class BellDualSet(dual_set.DualSet): def __init__(self, ref_el): entity_ids = {} nodes = [] cur = 0 # make nodes by getting points # need to do this dimension-by-dimension, facet-by-facet top = ref_el.get_topology() verts = ref_el.get_vertices() sd = ref_el.get_spatial_dimension() if ref_el.get_shape() != TRIANGLE: raise ValueError("Bell only defined on triangles") pd = functional.PointDerivative # get jet at each vertex entity_ids[0] = {} for v in sorted(top[0]): nodes.append(functional.PointEvaluation(ref_el, verts[v])) # first derivatives for i in range(sd): alpha = [0] * sd alpha[i] = 1 nodes.append(pd(ref_el, verts[v], alpha)) # second derivatives alphas = [[2, 0], [1, 1], [0, 2]] for alpha in alphas: nodes.append(pd(ref_el, verts[v], alpha)) entity_ids[0][v] = list(range(cur, cur + 6)) cur += 6 # we need an edge quadrature rule for the moment from FIAT.quadrature_schemes import create_quadrature from FIAT.jacobi import eval_jacobi rline = ufc_simplex(1) q1d = create_quadrature(rline, 8) q1dpts = q1d.get_points() leg4_at_qpts = eval_jacobi(0, 0, 4, 2.0*q1dpts - 1) imond = functional.IntegralMomentOfNormalDerivative entity_ids[1] = {} for e in sorted(top[1]): entity_ids[1][e] = [18+e] nodes.append(imond(ref_el, e, q1d, leg4_at_qpts)) entity_ids[2] = {0: []} super(BellDualSet, self).__init__(nodes, ref_el, entity_ids)
[docs] class Bell(finite_element.CiarletElement): """The Bell finite element.""" def __init__(self, ref_el): poly_set = polynomial_set.ONPolynomialSet(ref_el, 5) dual = BellDualSet(ref_el) super(Bell, self).__init__(poly_set, dual, 5)