crouzeixraviart.py
Go to the documentation of this file.00001
00002
00003 from swiginac import *
00004 from SyFi import *
00005
00006 initSyFi(3)
00007
00008 x = cvar.x; y = cvar.y; z = cvar.z
00009
00010
00011 class CrouzeixRaviart:
00012 """
00013 Python implementation of the Crouzeix-Raviart element.
00014 The corresponding C++ implementation is in the
00015 file CrouzeixRaviart.cpp.
00016 """
00017
00018 def __init__(self, polygon):
00019 """ Constructor """
00020 self.Ns = []
00021 self.dofs = []
00022 self.polygon = polygon
00023 self.compute_basis_functions()
00024
00025 def compute_basis_functions(self):
00026 """
00027 Compute the basis functions and degrees of freedom
00028 and put them in Ns and dofs, respectively.
00029 """
00030 polspace = bernstein(1,triangle,"a")
00031 N = polspace[0]
00032 variables = polspace[1]
00033
00034 for i in range(0,3):
00035 line = triangle.line(i)
00036 dofi = line.integrate(N)
00037 self.dofs.append(dofi)
00038
00039 for i in range(0,3):
00040 equations = []
00041 for j in range(0,3):
00042 equations.append(relational(self.dofs[j], dirac(i,j)))
00043 sub = lsolve(equations, variables)
00044 Ni = N.subs(sub)
00045 self.Ns.append(Ni);
00046
00047 def N(self,i): return self.Ns[i]
00048 def dof(self,i): return self.dofs[i]
00049 def nbf(self): return len(self.Ns)
00050
00051
00052 p0 = [0,0,0]; p1 = [1,0,0]; p2 = [0,1,0];
00053
00054 triangle = Triangle(p0, p1, p2)
00055
00056 fe = CrouzeixRaviart(triangle)
00057 fe.compute_basis_functions()
00058 print fe.nbf()
00059
00060 for i in range(0,fe.nbf()):
00061 print "N(%d) = "%i, fe.N(i).eval().printc()
00062 print "grad(N(%d)) = "%i, grad(fe.N(i)).eval().printc()
00063 print "dof(%d) = "%i, fe.dof(i).eval().printc()
00064
00065
00066
00067
00068
00069