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 edited 9 May 2005
9
import shapes, functional, polynomial, points, functionalset
11
def PhiK( shape , n , U ):
12
"""The set PhiK defined by Brezzi & Fortin of vector-valued
13
polynomials of degree n with vanishing divergence and
14
vanishing normal components around the boundary."""
15
U = polynomial.OrthogonalPolynomialArraySet( shape , n )
16
# need Phi and U to share a common base, so I have to
17
# take degree n and then throw away some to get degree n-1
18
Phi = polynomial.OrthogonalPolynomialSet( shape , n )
19
Phi_lower = Phi[:shapes.polynomial_dimension( shape , n-1 )]
20
DCPE = functional.DirectionalComponentPointEvaluation
21
d = shapes.dimension( shape )
22
pts_per_edge = [ [ x \
23
for x in points.make_points( shape , \
27
for i in shapes.entity_range( shape , d-1 ) ]
28
pts_flat = reduce( lambda a,b:a+b , pts_per_edge )
30
ls = [ functional.IntegralMomentOfDivergence( U , phi ) \
31
for phi in Phi_lower ]
34
nrml = shapes.normals[ shape ][ i ]
35
ls_cur = [ DCPE(U,nrml,pt) for pt in pts_per_edge[i] ]
37
fset = functionalset.FunctionalSet(U,ls)
39
return polynomial.ConstrainedPolynomialSet( fset )