SyFi
0.3
|
00001 #!/usr/bin/python 00002 00003 from swiginac import * 00004 from SyFi import * 00005 00006 initSyFi(3) 00007 00008 x = cvar.x; y = cvar.y; z = cvar.z # fetch some global variables 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