1
# Written by Robert C. Kirby
2
# Copyright 2005 by The University of Chicago
3
# Distributed under the LGPL license
4
# This work is partially supported by the US Department of Energy
5
# under award number DE-FG02-04ER25650
7
# Last modified 16 May 2005
10
"""This module defines the standard Lagrange finite element over
11
any shape for which an appropriate expansion basis is defined. The
12
module simply provides the class Lagrange, which takes a shape and a
16
import dualbasis, functional, points, shapes, polynomial, functionalset
18
class LagrangeDual( dualbasis.DualBasis ):
19
"""Dual basis for the plain vanilla Lagrange finite element."""
20
def __init__( self , shape , n , U ):
21
# Nodes are point evaluation on the lattice, which here is
22
# ordered by topological dimension and entity within that
23
# dimension. For triangles we have vertex points plus
24
# n - 1 points per edge. We then have max(0,(n-2)(n-1)/2)
25
# points on the interior of the triangle, or that many points
26
# per face of tetrahedra plus max(0,(n-2)(n-1)n/6 points
28
pts = [ [ points.make_points(shape,d,e,n) \
29
for e in shapes.entity_range(shape,d) ] \
30
for d in shapes.dimension_range(shape) ]
32
pts_flat = reduce( lambda a,b:a+b , \
33
reduce( lambda a,b:a+b , pts ) )
35
# get node location for each node i by self.pts[i]
38
ls = functional.make_point_evaluations( U , pts_flat )
40
# entity_ids is a dictionary whose keys are the topological
41
# dimensions (0,1,2,3) and values are dictionaries mapping
42
# the ids of entities of that dimension to the list of
43
# ids of nodes associated with that entity.
44
# see FIAT.base.dualbasis for description of entity_ids
47
for d in shapes.dimension_range( shape ):
49
for v in shapes.entity_range( shape, d ):
50
num_nods = len( pts[d][v] )
51
entity_ids[d][v] = range(id_cur,id_cur + num_nods)
54
dualbasis.DualBasis.__init__( self , \
55
functionalset.FunctionalSet( U , ls ) , \
60
class Lagrange( polynomial.FiniteElement ):
61
def __init__( self , shape , n ):
64
"Lagrange elements are only defined for n >= 1"
65
U = polynomial.OrthogonalPolynomialSet( shape , n )
66
Udual = LagrangeDual( shape , n , U )
67
polynomial.FiniteElement.__init__( self , \
70
class VectorLagrangeDual( dualbasis.DualBasis ):
71
def __init__( self , shape , n , U ):
72
space_dim = shapes.dimension( shape )
73
nc = U.tensor_shape()[0]
74
pts = [ [ points.make_points(shape,d,e,n) \
75
for e in shapes.entity_range(shape,d) ] \
76
for d in shapes.dimension_range(shape) ]
78
pts_flat = reduce( lambda a,b:a+b , \
79
reduce( lambda a,b:a+b , pts ) )
83
ls = [ functional.ComponentPointEvaluation( U , c , pt ) \
88
# see FIAT.base.dualbasis for description of entity_ids
91
for d in shapes.dimension_range( shape ):
93
for v in shapes.entity_range( shape, d ):
94
num_nods = len( pts[d][v] )
95
entity_ids[d][v] = range(id_cur,id_cur + num_nods)
98
fset = functionalset.FunctionalSet( U , ls )
100
dualbasis.DualBasis.__init__( self , fset , entity_ids , nc )
103
class VectorLagrange( polynomial.FiniteElement ):
104
def __init__( self , shape , n , nc = None ):
106
raise RuntimeError, \
107
"Lagrange elements are only defined for n >= 1"
108
U = polynomial.OrthogonalPolynomialArraySet( shape , n , nc )
109
Udual = VectorLagrangeDual( shape , n , U )
110
polynomial.FiniteElement.__init__( self , Udual , U )