#include <BrezziDouglasMarini.h>
Public Member Functions | |
BrezziDouglasMarini () | |
BrezziDouglasMarini (Polygon &p, int order=1, bool pointwise=true) | |
virtual | ~BrezziDouglasMarini () |
virtual void | compute_basis_functions () |
Public Attributes | |
bool | pointwise |
GiNaC::lst | dof_repr |
Definition at line 12 of file BrezziDouglasMarini.h.
SyFi::BrezziDouglasMarini::BrezziDouglasMarini | ( | ) |
Definition at line 16 of file BrezziDouglasMarini.cpp.
References SyFi::StandardFE::description.
00016 : StandardFE() 00017 { 00018 description = "BrezziDouglasMarini"; 00019 }
SyFi::BrezziDouglasMarini::BrezziDouglasMarini | ( | Polygon & | p, | |
int | order = 1 , |
|||
bool | pointwise = true | |||
) |
Definition at line 21 of file BrezziDouglasMarini.cpp.
References compute_basis_functions(), and pointwise.
00021 : StandardFE(p, order) 00022 { 00023 pointwise = pointwise_; 00024 compute_basis_functions(); 00025 }
virtual SyFi::BrezziDouglasMarini::~BrezziDouglasMarini | ( | ) | [inline, virtual] |
void SyFi::BrezziDouglasMarini::compute_basis_functions | ( | ) | [virtual] |
Reimplemented from SyFi::StandardFE.
Definition at line 27 of file BrezziDouglasMarini.cpp.
References SyFi::bernsteinv(), SyFi::collapse(), SyFi::StandardFE::description, SyFi::StandardFE::dofs, SyFi::inner(), SyFi::interior_coordinates(), SyFi::istr(), SyFi::Triangle::line(), SyFi::matrix_from_equations(), SyFi::normal(), SyFi::StandardFE::Ns, SyFi::StandardFE::order, SyFi::StandardFE::p, pointwise, SyFi::Polygon::str(), SyFi::t, demos::crouzeixraviart::triangle, SyFi::x, and SyFi::y.
Referenced by BrezziDouglasMarini().
00028 { 00029 00030 if ( order < 1 ) 00031 { 00032 throw(std::logic_error("Brezzi-Douglas-Marini elements must be of order 1 or higher.")); 00033 } 00034 00035 if ( p == NULL ) 00036 { 00037 throw(std::logic_error("You need to set a polygon before the basisfunctions can be computed")); 00038 } 00039 00040 GiNaC::ex nsymb = GiNaC::symbol("n"); 00041 if (pointwise) 00042 { 00043 00044 if ( p->str().find("ReferenceLine") != string::npos ) 00045 { 00046 00047 cout <<"Can not define the Brezzi-Douglas-Marini element on a line"<<endl; 00048 00049 } 00050 00051 if ( p->str().find("ReferenceTetrahedron") != string::npos ) 00052 { 00053 00054 cout <<"Can not define the Brezzi-Douglas-Marini element on a tetrahedron, see Nedelec2Hdiv"<<endl; 00055 00056 } 00057 else if ( p->str().find("Triangle") != string::npos ) 00058 { 00059 00060 description = istr("BrezziDouglasMarini_", order) + "_2D"; 00061 00062 Triangle& triangle = (Triangle&)(*p); 00063 GiNaC::lst equations; 00064 GiNaC::lst variables; 00065 GiNaC::lst polynom_space = bernsteinv(2,order, triangle, "b"); 00066 GiNaC::ex pspace = polynom_space.op(0); 00067 00068 variables = collapse(GiNaC::ex_to<GiNaC::lst>(polynom_space.op(1))); 00069 00070 GiNaC::symbol t("t"); 00071 GiNaC::ex dofi; 00072 // dofs related to edges 00073 for (int i=0; i< 3; i++) 00074 { 00075 Line line = triangle.line(i); 00076 GiNaC::lst normal_vec = normal(triangle, i); 00077 GiNaC::ex Vn = inner(pspace, normal_vec); 00078 GiNaC::lst points = interior_coordinates(line, order); 00079 00080 GiNaC::ex point; 00081 for (unsigned int j=0; j< points.nops(); j++) 00082 { 00083 point = points.op(j); 00084 dofi = Vn.subs(x == point.op(0)).subs(y == point.op(1)); 00085 GiNaC::ex eq = dofi == GiNaC::numeric(0); 00086 equations.append(eq); 00087 dofs.insert(dofs.end(), GiNaC::lst(point,nsymb)); 00088 } 00089 } 00090 00091 // std::cout <<"no variables "<<variables.nops()<<std::endl; 00092 // std::cout <<"no equations "<<equations.nops()<<std::endl; 00093 00094 // invert the matrix: 00095 // GiNaC has a bit strange way to invert a matrix. 00096 // It solves the system AA^{-1} = Id. 00097 // It seems that this way is the only way to do 00098 // properly with the solve_algo::gauss flag. 00099 // 00100 GiNaC::matrix b; GiNaC::matrix A; 00101 matrix_from_equations(equations, variables, A, b); 00102 00103 unsigned int ncols = A.cols(); 00104 GiNaC::matrix vars_sq(ncols, ncols); 00105 00106 // matrix of symbols 00107 for (unsigned r=0; r<ncols; ++r) 00108 for (unsigned c=0; c<ncols; ++c) 00109 vars_sq(r, c) = GiNaC::symbol(); 00110 00111 GiNaC::matrix id(ncols, ncols); 00112 00113 // identity 00114 const GiNaC::ex _ex1(1); 00115 for (unsigned i=0; i<ncols; ++i) 00116 id(i, i) = _ex1; 00117 00118 // invert the matrix 00119 GiNaC::matrix m_inv(ncols, ncols); 00120 m_inv = A.solve(vars_sq, id, GiNaC::solve_algo::gauss); 00121 00122 for (unsigned int i=0; i<dofs.size(); i++) 00123 { 00124 b.let_op(i) = GiNaC::numeric(1); 00125 GiNaC::ex xx = m_inv.mul(GiNaC::ex_to<GiNaC::matrix>(b)); 00126 00127 GiNaC::lst subs; 00128 for (unsigned int ii=0; ii<xx.nops(); ii++) 00129 { 00130 subs.append(variables.op(ii) == xx.op(ii)); 00131 } 00132 GiNaC::ex Nj1 = pspace.op(0).subs(subs); 00133 GiNaC::ex Nj2 = pspace.op(1).subs(subs); 00134 Ns.insert(Ns.end(), GiNaC::matrix(2,1,GiNaC::lst(Nj1,Nj2))); 00135 b.let_op(i) = GiNaC::numeric(0); 00136 } 00137 } 00138 } 00139 00140 }
GiNaC::lst SyFi::BrezziDouglasMarini::dof_repr |
Definition at line 16 of file BrezziDouglasMarini.h.
Definition at line 15 of file BrezziDouglasMarini.h.
Referenced by BrezziDouglasMarini(), and compute_basis_functions().