Skip to content

Commit

Permalink
Got the permutations rihgt finally
Browse files Browse the repository at this point in the history
Signed-off-by: Umberto Zerbinati <[email protected]>
  • Loading branch information
Umberto Zerbinati committed Dec 14, 2024
1 parent f7dafda commit 4bed5ff
Showing 1 changed file with 20 additions and 14 deletions.
34 changes: 20 additions & 14 deletions FIAT/crouzeix_raviart.py
Original file line number Diff line number Diff line change
Expand Up @@ -10,6 +10,7 @@
# Last changed: 2010-01-28

from FIAT import finite_element, polynomial_set, dual_set, functional
from FIAT.check_format_variant import parse_lagrange_variant


def _initialize_entity_ids(topology):
Expand Down Expand Up @@ -62,27 +63,26 @@ def __init__(self, cell, degree):
entity_ids = _initialize_entity_ids(topology)
nodes = [None for i in range(10)]
# Construct nodes and entity_ids
entity_permutations = {0: {0: {0: []}, 1: {0: []}, 2: {0: []}},
1: {0: {0: [0, 1, 2], 1: [2, 1, 0]}, 1: {0: [0, 1, 2], 1: [2, 1, 0]}, 2: {0: [0, 1, 2], 1: [2, 1, 0]}},
2: {0: {0: [0], 1: [0], 2: [0], 3: [0], 4: [0], 5: [0]}}
}
dof_counter = 0
for i in topology[d - 1]:

# Construct midpoint
x = cell.make_points(d-1,i,4)
print(x)
x = cell.make_points(d-1,i,4, variant="equi")
# Degree of freedom number i is evaluation at midpoint
nodes[dof_counter] = functional.PointEvaluation(cell, x[0])
nodes[dof_counter+1] = functional.PointEvaluation(cell, x[1])
nodes[dof_counter] = functional.PointEvaluation(cell, x[1])
nodes[dof_counter+1] = functional.PointEvaluation(cell, x[0])
nodes[dof_counter+2] = functional.PointEvaluation(cell, x[2])
entity_ids[d - 1][i] += [dof_counter]
entity_ids[d - 1][i] += [dof_counter+1]
entity_ids[d - 1][i] += [dof_counter+2]
dof_counter += 3
x = cell.make_points(d,0,3)
x = cell.make_points(d,0,3, variant="equi")
nodes[-1] = functional.PointEvaluation(cell,x[0])
entity_ids[2][0] += [9]
print("N",nodes)
print("E",entity_ids)
# Initialize super-class
super(CrouzeixRaviartThreeDualSet, self).__init__(nodes, cell, entity_ids)
super(CrouzeixRaviartThreeDualSet, self).__init__(nodes, cell, entity_ids, entity_permutations)



Expand All @@ -94,7 +94,7 @@ class CrouzeixRaviart(finite_element.CiarletElement):
Dual basis: Evaluation at facet midpoints
"""

def __init__(self, cell, degree):
def __init__(self, cell, degree,variant=None):

# Crouzeix Raviart is only defined for polynomial degree == 1
if not (degree in [1,3]):
Expand All @@ -106,7 +106,13 @@ def __init__(self, cell, degree):
dual = CrouzeixRaviartDualSet(cell, 1)
super(CrouzeixRaviart, self).__init__(space, dual, 1)
elif degree == 3:
space = polynomial_set.ONPolynomialSet(cell, 3)
dual = CrouzeixRaviartThreeDualSet(cell, 3)
super(CrouzeixRaviart, self).__init__(space, dual, 3)
ref_el = cell
splitting, point_variant = parse_lagrange_variant(variant)
if splitting is not None:
ref_el = splitting(ref_el)
poly_variant = "bubble" if ref_el.is_macrocell() else None
space = polynomial_set.ONPolynomialSet(ref_el, degree, variant=poly_variant)
formdegree = 0 # 0-form
dual = CrouzeixRaviartThreeDualSet(cell, degree)
super(CrouzeixRaviart, self).__init__(space, dual, 3, formdegree)

0 comments on commit 4bed5ff

Please sign in to comment.