SyFi 0.3
SyFi::CrouzeixRaviart Class Reference

#include <CrouzeixRaviart.h>

Inheritance diagram for SyFi::CrouzeixRaviart:
SyFi::StandardFE SyFi::StandardFE SyFi::FE SyFi::FE SyFi::FE SyFi::FE SyFi::_object SyFi::_object SyFi::_object SyFi::_object

List of all members.

Public Member Functions

 CrouzeixRaviart ()
 CrouzeixRaviart (Polygon &p, unsigned int order=1)
virtual ~CrouzeixRaviart ()
void compute_basis_functions ()
def __init__
def compute_basis_functions

Public Attributes

 this

Static Private Attributes

dictionary __swig_setmethods__ = {}
tuple __setattr__ = lambdaself,name,value:_swig_setattr(self, CrouzeixRaviart, name, value)
dictionary __swig_getmethods__ = {}
tuple __getattr__ = lambdaself,name:_swig_getattr(self, CrouzeixRaviart, name)
 __repr__ = _swig_repr
 __swig_destroy__ = _SyFi.delete_CrouzeixRaviart
 __del__ = lambdaself:None;

Detailed Description

Proxy of C++ SyFi::CrouzeixRaviart class

Definition at line 12 of file CrouzeixRaviart.h.


Constructor & Destructor Documentation

SyFi::CrouzeixRaviart::CrouzeixRaviart ( )

Definition at line 15 of file CrouzeixRaviart.cpp.

References SyFi::StandardFE::order.

                                          : StandardFE()
        {
                order = 1;
        }
SyFi::CrouzeixRaviart::CrouzeixRaviart ( Polygon p,
unsigned int  order = 1 
)

Definition at line 20 of file CrouzeixRaviart.cpp.

References compute_basis_functions().

virtual SyFi::CrouzeixRaviart::~CrouzeixRaviart ( ) [inline, virtual]

Definition at line 17 of file CrouzeixRaviart.h.

{}
def SyFi::CrouzeixRaviart::__init__ (   self,
  args 
)
__init__(self) -> CrouzeixRaviart
__init__(self, Polygon p, unsigned int order = 1) -> CrouzeixRaviart
__init__(self, Polygon p) -> CrouzeixRaviart

Reimplemented from SyFi::StandardFE.

Definition at line 2275 of file SyFi.py.

02275                              : 
02276         """
02277         __init__(self) -> CrouzeixRaviart
02278         __init__(self, Polygon p, unsigned int order = 1) -> CrouzeixRaviart
02279         __init__(self, Polygon p) -> CrouzeixRaviart
02280         """
02281         this = _SyFi.new_CrouzeixRaviart(*args)
02282         try: self.this.append(this)
02283         except: self.this = this

Member Function Documentation

void SyFi::CrouzeixRaviart::compute_basis_functions ( ) [virtual]

Reimplemented from SyFi::StandardFE.

Definition at line 25 of file CrouzeixRaviart.cpp.

References SyFi::bernstein(), SyFi::StandardFE::description, SyFi::dirac(), SyFi::StandardFE::dofs, SyFi::Line::integrate(), SyFi::Triangle::line(), SyFi::StandardFE::Ns, SyFi::StandardFE::order, SyFi::StandardFE::p, SyFi::Polygon::str(), SyFi::sub(), SyFi::t, demos::geom_test::tetrahedron, SyFi::Tetrahedron::triangle(), demos::crouzeixraviart::triangle, SyFi::Polygon::vertex(), SyFi::x, SyFi::y, and SyFi::z.

Referenced by check_CrouzeixRaviart(), SyFi::VectorCrouzeixRaviart::compute_basis_functions(), and CrouzeixRaviart().

        {

                // remove previously computed basis functions and dofs
                Ns.clear();
                dofs.clear();

                if ( p == NULL )
                {
                        throw(std::logic_error("You need to set a polygon before the basisfunctions can be computed"));
                }

                if (order != 1)
                {
                        throw(std::logic_error("Only Crouziex-Raviart elements of order 1 is possible"));
                }

                // see e.g. Brezzi and Fortin book page 116 for the definition

                if ( p->str().find("ReferenceLine") != string::npos )
                {
                        cout <<"Can not define the Raviart-Thomas element on a line"<<endl;
                }
                else if ( p->str().find("Triangle") != string::npos )
                {

                        description = "CrouzeixRaviart_2D";

                        Triangle& triangle = (Triangle&)(*p);

                        // create the polynomial space
                        GiNaC::ex polynom_space = bernstein(1, triangle, "a");
                        GiNaC::ex polynom = polynom_space.op(0);
                        GiNaC::lst variables = GiNaC::ex_to<GiNaC::lst>(polynom_space.op(1));
                        GiNaC::ex basis = polynom_space.op(2);

                        // create the dofs
                        GiNaC::symbol t("t");
                        for (int j=0; j< 3; j++)
                        {

                                // solve the linear system to compute
                                // each of the basis functions
                                GiNaC::lst equations;
                                for (int i=0; i< 3; i++)
                                {

                                        Line line = triangle.line(i);
                                        //                          GiNaC::ex dofi = line.integrate(polynom);
                                        GiNaC::lst midpoint = GiNaC::lst(
                                                (line.vertex(0).op(0) +  line.vertex(1).op(0))/2,
                                                (line.vertex(0).op(1) +  line.vertex(1).op(1))/2);
                                        dofs.insert(dofs.end(), midpoint);

                                        //                          GiNaC::ex dofi = polynom.subs( x == midpoint.op(0)).subs( y == midpoint.op(1));
                                        GiNaC::ex dofi = line.integrate(polynom);
                                        equations.append( dofi == dirac(i,j));

                                        if (j == 1)
                                        {
                                                //                                  GiNaC::lst d = GiNaC::lst(line.vertex(0) ,  line.vertex(1));
                                                dofs.insert(dofs.end(), midpoint);
                                        }

                                }
                                GiNaC::ex sub = lsolve(equations, variables);
                                GiNaC::ex Ni = polynom.subs(sub);
                                Ns.insert(Ns.end(),Ni);
                        }

                }
                else if ( p->str().find("Tetrahedron") != string::npos )
                {

                        description = "CrouzeixRaviart_3D";

                        Tetrahedron& tetrahedron = (Tetrahedron&)(*p);
                        GiNaC::ex polynom_space = bernstein(1, tetrahedron, "a");
                        GiNaC::ex polynom = polynom_space.op(0);
                        GiNaC::lst variables = GiNaC::ex_to<GiNaC::lst>(polynom_space.op(1));
                        GiNaC::ex basis = polynom_space.op(2);

                        GiNaC::ex bernstein_pol;

                        GiNaC::symbol t("t");
                        // dofs related to edges
                        for (int j=0; j< 4; j++)
                        {

                                GiNaC::lst equations;
                                for (int i=0; i< 4; i++)
                                {
                                        Triangle triangle = tetrahedron.triangle(i);
                                        GiNaC::lst midpoint = GiNaC::lst(
                                                (triangle.vertex(0).op(0) + triangle.vertex(1).op(0) + triangle.vertex(2).op(0))/3,
                                                (triangle.vertex(0).op(1) + triangle.vertex(1).op(1) + triangle.vertex(2).op(1))/3,
                                                (triangle.vertex(0).op(2) + triangle.vertex(1).op(2) + triangle.vertex(2).op(2))/3
                                                );

                                        GiNaC::ex dofi = polynom.subs(x == midpoint.op(0)).subs(y == midpoint.op(1)).subs(z == midpoint.op(2));
                                        equations.append( dofi == dirac(i,j));

                                        if ( j == 1 )
                                        {
                                                //                                  GiNaC::lst d = GiNaC::lst(triangle.vertex(0), triangle.vertex(1), triangle.vertex(2));
                                                dofs.insert(dofs.end(), midpoint);
                                        }
                                }
                                GiNaC::ex sub = lsolve(equations, variables);
                                GiNaC::ex Ni = polynom.subs(sub);
                                Ns.insert(Ns.end(),Ni);

                        }
                }
        }
def SyFi::CrouzeixRaviart::compute_basis_functions (   self)
compute_basis_functions(self)

Reimplemented from SyFi::StandardFE.

Definition at line 2286 of file SyFi.py.

02287                                      :
02288         """compute_basis_functions(self)"""
02289         return _SyFi.CrouzeixRaviart_compute_basis_functions(self)


Member Data Documentation

SyFi::CrouzeixRaviart::__del__ = lambdaself:None; [static, private]

Reimplemented from SyFi::StandardFE.

Definition at line 2285 of file SyFi.py.

tuple SyFi::CrouzeixRaviart::__getattr__ = lambdaself,name:_swig_getattr(self, CrouzeixRaviart, name) [static, private]

Reimplemented from SyFi::StandardFE.

Definition at line 2273 of file SyFi.py.

SyFi::CrouzeixRaviart::__repr__ = _swig_repr [static, private]

Reimplemented from SyFi::StandardFE.

Definition at line 2274 of file SyFi.py.

tuple SyFi::CrouzeixRaviart::__setattr__ = lambdaself,name,value:_swig_setattr(self, CrouzeixRaviart, name, value) [static, private]

Reimplemented from SyFi::StandardFE.

Definition at line 2270 of file SyFi.py.

SyFi::CrouzeixRaviart::__swig_destroy__ = _SyFi.delete_CrouzeixRaviart [static, private]

Reimplemented from SyFi::StandardFE.

Definition at line 2284 of file SyFi.py.

dictionary SyFi::CrouzeixRaviart::__swig_getmethods__ = {} [static, private]

Reimplemented from SyFi::StandardFE.

Definition at line 2271 of file SyFi.py.

dictionary SyFi::CrouzeixRaviart::__swig_setmethods__ = {} [static, private]

Reimplemented from SyFi::StandardFE.

Definition at line 2268 of file SyFi.py.

Reimplemented from SyFi::StandardFE.

Definition at line 2279 of file SyFi.py.


The documentation for this class was generated from the following files:
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Defines