SyFi
0.3
|
Classes | |
class | ReferenceSimplex |
class | Lagrange |
Functions | |
def | bernstein_space |
def | create_point_set |
def | create_matrix |
Variables | |
tuple | x = Symbol('x') |
tuple | y = Symbol('y') |
tuple | z = Symbol('z') |
tuple | t = ReferenceSimplex(2) |
f = x+y | |
tuple | fe = Lagrange(2,2) |
int | u = 0 |
list | us = [] |
tuple | ui = Symbol("u_%d" % i) |
tuple | J = zeronm(fe.nbf(), fe.nbf()) |
list | Fi = u*fe.N[i] |
list | uj = us[j] |
tuple | integrands = diff(Fi, uj) |
def fem_sympy_core.bernstein_space | ( | order, | |
nsd | |||
) |
Definition at line 39 of file fem_sympy_core.py.
00039 00040 def bernstein_space(order, nsd): 00041 if nsd > 3: 00042 raise RuntimeError("Bernstein only implemented in 1D, 2D, and 3D") 00043 sum = 0 00044 basis = [] 00045 coeff = [] 00046 00047 if nsd == 1: 00048 b1, b2 = x, 1-x 00049 for o1 in range(0,order+1): 00050 for o2 in range(0,order+1): 00051 if o1 + o2 == order: 00052 aij = Symbol("a_%d_%d" % (o1,o2)) 00053 sum += aij*binomial(order,o1)*pow(b1, o1)*pow(b2, o2) 00054 basis.append(binomial(order,o1)*pow(b1, o1)*pow(b2, o2)) 00055 coeff.append(aij) 00056 00057 00058 if nsd == 2: 00059 b1, b2, b3 = x, y, 1-x-y 00060 for o1 in range(0,order+1): 00061 for o2 in range(0,order+1): 00062 for o3 in range(0,order+1): 00063 if o1 + o2 + o3 == order: 00064 aij = Symbol("a_%d_%d_%d" % (o1,o2,o3)) 00065 fac = factorial(order)/(factorial(o1)*factorial(o2)*factorial(o3)) 00066 sum += aij*fac*pow(b1, o1)*pow(b2, o2)*pow(b3, o3) 00067 basis.append(fac*pow(b1, o1)*pow(b2, o2)*pow(b3, o3)) 00068 coeff.append(aij) 00069 00070 if nsd == 3: 00071 b1, b2, b3, b4 = x, y, z, 1-x-y-z 00072 for o1 in range(0,order+1): 00073 for o2 in range(0,order+1): 00074 for o3 in range(0,order+1): 00075 for o4 in range(0,order+1): 00076 if o1 + o2 + o3 + o4 == order: 00077 aij = Symbol("a_%d_%d_%d_%d" % (o1,o2,o3,o4)) 00078 fac = factorial(order)/(factorial(o1)*factorial(o2)*factorial(o3)*factorial(o4)) 00079 sum += aij*fac*pow(b1, o1)*pow(b2, o2)*pow(b3, o3)*pow(b4, o4) 00080 basis.append(fac*pow(b1, o1)*pow(b2, o2)*pow(b3, o3)*pow(b4, o4)) 00081 coeff.append(aij) 00082 00083 00084 return sum, coeff, basis
def fem_sympy_core.create_matrix | ( | equations, | |
coeffs | |||
) |
Definition at line 117 of file fem_sympy_core.py.
00117 00118 def create_matrix(equations, coeffs): 00119 A = Matrix(len(equations), len(equations)) 00120 i = 0; j = 0 00121 for j in range(0, len(coeffs)): 00122 c = coeffs[j] 00123 for i in range(0, len(equations)): 00124 e = equations[i] 00125 d = diff(e, c) 00126 A[i,j] = d 00127 return A 00128 00129
def fem_sympy_core.create_point_set | ( | order, | |
nsd | |||
) |
Definition at line 85 of file fem_sympy_core.py.
00085 00086 def create_point_set(order, nsd): 00087 h = Rational(1,order) 00088 set = [] 00089 00090 if nsd == 1: 00091 for i in range(0, order+1): 00092 x = i*h 00093 if x <= 1: 00094 set.append((x,y)) 00095 00096 if nsd == 2: 00097 for i in range(0, order+1): 00098 x = i*h 00099 for j in range(0, order+1): 00100 y = j*h 00101 if x + y <= 1: 00102 set.append((x,y)) 00103 00104 if nsd == 3: 00105 for i in range(0, order+1): 00106 x = i*h 00107 for j in range(0, order+1): 00108 y = j*h 00109 for k in range(0, order+1): 00110 z = j*h 00111 if x + y + z <= 1: 00112 set.append((x,y,z)) 00113 00114 return set 00115 00116
Definition at line 177 of file fem_sympy_core.py.
tuple fem_sympy_core::fe = Lagrange(2,2) |
Definition at line 188 of file fem_sympy_core.py.
list fem_sympy_core::Fi = u*fe.N[i] |
Definition at line 206 of file fem_sympy_core.py.
tuple fem_sympy_core::integrands = diff(Fi, uj) |
Definition at line 209 of file fem_sympy_core.py.
tuple fem_sympy_core::J = zeronm(fe.nbf(), fe.nbf()) |
Definition at line 204 of file fem_sympy_core.py.
tuple fem_sympy_core::t = ReferenceSimplex(2) |
Definition at line 176 of file fem_sympy_core.py.
int fem_sympy_core::u = 0 |
Definition at line 195 of file fem_sympy_core.py.
tuple fem_sympy_core::ui = Symbol("u_%d" % i) |
Definition at line 199 of file fem_sympy_core.py.
list fem_sympy_core::uj = us[j] |
Definition at line 208 of file fem_sympy_core.py.
list fem_sympy_core::us = [] |
Definition at line 197 of file fem_sympy_core.py.
tuple fem_sympy_core::x = Symbol('x') |
Definition at line 8 of file fem_sympy_core.py.
Referenced by _p_SyFi__ArnoldFalkWintherWeakSymPTo_p_SyFi__StandardFE(), _p_SyFi__ArnoldFalkWintherWeakSymSigmaTo_p_SyFi__StandardFE(), _p_SyFi__ArnoldFalkWintherWeakSymUTo_p_SyFi__StandardFE(), _p_SyFi__BoxTo_p_SyFi__Polygon(), _p_SyFi__BubbleTo_p_SyFi__StandardFE(), _p_SyFi__CrouzeixRaviartTo_p_SyFi__StandardFE(), _p_SyFi__DiscontinuousLagrangeTo_p_SyFi__FE(), _p_SyFi__DiscontinuousLagrangeTo_p_SyFi__Lagrange(), _p_SyFi__HermiteTo_p_SyFi__StandardFE(), _p_SyFi__LagrangeTo_p_SyFi__StandardFE(), _p_SyFi__LineTo_p_SyFi__Polygon(), _p_SyFi__MixedFETo_p_SyFi__FE(), _p_SyFi__Nedelec2HdivTo_p_SyFi__StandardFE(), _p_SyFi__NedelecTo_p_SyFi__StandardFE(), _p_SyFi__P0To_p_SyFi__StandardFE(), _p_SyFi__RaviartThomasTo_p_SyFi__StandardFE(), _p_SyFi__RectangleTo_p_SyFi__Polygon(), _p_SyFi__ReferenceBoxTo_p_SyFi__Box(), _p_SyFi__ReferenceLineTo_p_SyFi__Line(), _p_SyFi__ReferenceRectangleTo_p_SyFi__Rectangle(), _p_SyFi__ReferenceTetrahedronTo_p_SyFi__Tetrahedron(), _p_SyFi__ReferenceTriangleTo_p_SyFi__Triangle(), _p_SyFi__RobustTo_p_SyFi__StandardFE(), _p_SyFi__SimplexTo_p_SyFi__Polygon(), _p_SyFi__SpaceTimeDomainTo_p_SyFi__Polygon(), _p_SyFi__SpaceTimeElementTo_p_SyFi__StandardFE(), _p_SyFi__StandardFETo_p_SyFi__FE(), _p_SyFi__TensorLagrangeTo_p_SyFi__StandardFE(), _p_SyFi__TensorP0To_p_SyFi__StandardFE(), _p_SyFi__TetrahedronTo_p_SyFi__Polygon(), _p_SyFi__TriangleTo_p_SyFi__Polygon(), _p_SyFi__VectorCrouzeixRaviartTo_p_SyFi__StandardFE(), _p_SyFi__VectorDiscontinuousLagrangeTo_p_SyFi__FE(), _p_SyFi__VectorDiscontinuousLagrangeTo_p_SyFi__VectorLagrange(), _p_SyFi__VectorLagrangeTo_p_SyFi__StandardFE(), _p_SyFi__VectorP0To_p_SyFi__StandardFE(), SyFi.barycenter_line(), SyFi.barycenter_tetrahedron(), SyFi.barycenter_triangle(), code_gen2D(), SyFi.coeff(), SyFi.coeffs(), SyFi::Hermite.compute_basis_functions(), SyFi::Lagrange.compute_basis_functions(), SyFi::CrouzeixRaviart.compute_basis_functions(), SyFi::Nedelec.compute_basis_functions(), SyFi::Bubble.compute_basis_functions(), SyFi::Nedelec2Hdiv.compute_basis_functions(), SyFi::RaviartThomas.compute_basis_functions(), SyFi::Robust.compute_basis_functions(), SyFi::BrezziDouglasMarini.compute_basis_functions(), SyFi::SpaceTimeElement.compute_basis_functions(), SyFi::Robust.compute_basis_functions_old(), SyFi.div(), SyFi.ex2equations(), SyFi.grad(), SyFi.homogenous_pol(), SyFi.initSyFi(), SyFi::SpaceTimeDomain.integrate(), SyFi.legendre(), SyFi.legendre1D(), main(), SyFi.matvec(), pickExpression(), SyFi.pol(), SyFi.pol2basisandcoeff(), SyFi.polb(), Ptv.Ptv(), SyFi::Triangle.repr(), SyFi::Rectangle.repr(), SyFi::Tetrahedron.repr(), SyFi::Box.repr(), SyFi::Simplex.repr(), std_list_Sl_GiNaC_ex_Sg__pop(), std_list_Sl_std_pair_Sl_GiNaC_symbol_Sc_GiNaC_ex_Sg__Sg__pop(), std_map_Sl_GiNaC_ex_Sc_GiNaC_ex_Sc_GiNaC_ex_is_less_Sg____setitem____SWIG_1(), std_map_Sl_GiNaC_ex_Sc_int_Sc_GiNaC_ex_is_less_Sg____setitem____SWIG_1(), std_vector_Sl_GiNaC_ex_Sg__pop(), SWIG_CanCastAsInteger(), Swig_var_x_get(), Swig_var_x_set(), and variants().
tuple fem_sympy_core::y = Symbol('y') |
Definition at line 9 of file fem_sympy_core.py.
Referenced by code_gen2D(), SyFi.coeff(), SyFi.coeffs(), SyFi::CrouzeixRaviart.compute_basis_functions(), SyFi::Hermite.compute_basis_functions(), SyFi::Lagrange.compute_basis_functions(), SyFi::Nedelec.compute_basis_functions(), SyFi::Nedelec2Hdiv.compute_basis_functions(), SyFi::Robust.compute_basis_functions(), SyFi::RaviartThomas.compute_basis_functions(), SyFi::BrezziDouglasMarini.compute_basis_functions(), SyFi::Robust.compute_basis_functions_old(), SyFi.div(), SyFi.ex2equations(), SyFi.grad(), SyFi.homogenous_pol(), SyFi.initSyFi(), SyFi.legendre(), main(), pickExpression(), SyFi.pol(), SyFi.pol2basisandcoeff(), SyFi.polb(), Ptv.Ptv(), Swig_var_y_get(), Swig_var_y_set(), and variants().
tuple fem_sympy_core::z = Symbol('z') |
Definition at line 10 of file fem_sympy_core.py.
Referenced by SyFi.coeff(), SyFi.coeffs(), SyFi::CrouzeixRaviart.compute_basis_functions(), SyFi::Lagrange.compute_basis_functions(), SyFi::Nedelec.compute_basis_functions(), SyFi::Hermite.compute_basis_functions(), SyFi::Nedelec2Hdiv.compute_basis_functions(), SyFi::RaviartThomas.compute_basis_functions(), SyFi.div(), SyFi.ex2equations(), SyFi.grad(), SyFi.homogenous_pol(), SyFi.initSyFi(), SyFi.legendre(), main(), pickExpression(), SyFi.pol(), SyFi.pol2basisandcoeff(), SyFi.polb(), Ptv.Ptv(), Swig_var_z_get(), Swig_var_z_set(), and variants().