fiat team mailing list archive

[Andy Terrel <aterrel@xxxxxxxxxxxx>]

Hello,

Here is a patch for Crouzeix-Raviart.  It fixes a bug with the 
shapes.dim versus shapes.dimension and adds a VectorCrouzeixRaviart 
element pretty much following VectorLagrange class.  I needed this for 
some stokes stuff with CrouzeixRaviart-P0 elements.


Thanks,
Andy

--

====================
Andy Terrel
Computer Science Dept
University of Chicago
aterrel@xxxxxxxxxxxx
---------------------

In mathematics you don't understand things. 
You just get used to them.
                      - Johann Von Neumann

diff -u --new-file --recursive FIAT-0.3.0/FIAT/CrouzeixRaviart.py FIAT-0.3.0-mod/FIAT/CrouzeixRaviart.py
--- FIAT-0.3.0/FIAT/CrouzeixRaviart.py	2006-11-30 16:52:36.000000000 -0600
+++ FIAT-0.3.0-mod/FIAT/CrouzeixRaviart.py	2007-01-13 22:24:47.000000000 -0600
@@ -15,9 +15,10 @@
     def __init__( self , shape , U ):
         # in d dimensions, evaluate at midpoints of (d-1)-dimensional
         # entities
-        d = shapes.dims[ shape ]
+        d = shapes.dimension(shape)
         pts = [ pt for i in shapes.entity_range(shape,d-1) \
                 for pt in shapes.make_points( shape , d-1 , i , d ) ]
+	self.pts = pts
         ls = functional.make_point_evaluations( U , pts )
         entity_ids = {}
         for i in range(d-1):
@@ -40,6 +41,42 @@
         polynomial.FiniteElement.__init__( self , Udual , U )
         return
 
+class VectorCrouzeixRaviartDual( dualbasis.DualBasis ):
+    """Dual basis for Crouzeix-Raviart element (linears continuous at
+    boundary midpoints)"""
+    def __init__( self , shape , U ):
+        # in d dimensions, evaluate at midpoints of (d-1)-dimensional
+        # entities
+	nc = U.tensor_shape()[0]
+        d = shapes.dimension(shape)
+        pts = [ pt for i in shapes.entity_range(shape,d-1) \
+                for pt in shapes.make_points( shape , d-1 , i , d ) ]
+	self.pts = pts
+        ls = [ functional.ComponentPointEvaluation( U , c, pt ) \
+	       for c in range(nc) \
+	       for pt in pts ]
+
+        entity_ids = {}
+        for i in range(d-1):
+            entity_ids[i] = {}
+            for j in shapes.entity_range(shape,i):
+                entity_ids[i][j] = []
+        entity_ids[d-1] = {}
+        for i in shapes.entity_range(shape,d-1):
+            entity_ids[d-1][i] = [i]
+        entity_ids[d] = { 0: [] }
+
+        fset = functionalset.FunctionalSet( U , ls )
+        
+        dualbasis.DualBasis.__init__( self , fset , entity_ids,nc )
+
+class VectorCrouzeixRaviart( polynomial.FiniteElement ):
+    def __init__( self , shape , nc = None):
+        U = polynomial.OrthogonalPolynomialArraySet( shape , 1, nc )
+        Udual = VectorCrouzeixRaviartDual( shape , U )
+        polynomial.FiniteElement.__init__( self , Udual , U )
+        return
+