[cig-commits] r13493 - cs/cigma/trunk/src

[luis at geodynamics.org]

Author: luis
Date: 2008-12-09 18:12:31 -0800 (Tue, 09 Dec 2008)
New Revision: 13493

Modified:
   cs/cigma/trunk/src/FE.cpp
   cs/cigma/trunk/src/FE.h
Log:
Updates to cigma::FE

Modified: cs/cigma/trunk/src/FE.cpp
===================================================================
--- cs/cigma/trunk/src/FE.cpp	2008-12-10 02:12:29 UTC (rev 13492)
+++ cs/cigma/trunk/src/FE.cpp	2008-12-10 02:12:31 UTC (rev 13493)
@@ -1,5 +1,6 @@
 #include "FE.h"
 using namespace cigma;
+using boost::shared_ptr;
 
 FE::FE()
 {
@@ -10,9 +11,58 @@
 
 FE::FE(const FE& other)
 {
+    this->init(other.quadrature);
 }
 
 FE::~FE()
 {
+    if (jxw != 0) { delete [] jxw; jxw = 0; }
+    if (basis_tab != 0) { delete [] basis_tab; basis_tab = 0; }
+    if (basis_jet != 0) { delete [] basis_jet; basis_jet = 0; }
 }
 
+void FE::init(shared_ptr<Quadrature> Q)
+{
+    this->quadrature = shared_ptr<Quadrature>(new Quadrature(*Q));
+    this->cell = Cell::New(Q->getCellType());
+    this->init_basis();
+    //this->init_basis_jet();
+}
+
+void FE::init_basis()
+{
+    assert(cell);
+    assert(quadrature);
+    const int nq = quadrature->n_points();
+    const int ndofs = cell->n_nodes();
+    if (basis_tab != 0) { delete [] basis_tab; }
+    basis_tab = new double[nq * ndofs];
+    cell->shape(nq, quadrature->points, basis_tab);
+}
+
+void FE::init_basis_jet()
+{
+    assert(cell);
+    assert(quadrature);
+    const int nq = quadrature->n_points();
+    const int ndofs = cell->n_nodes();
+    const int ndim = quadrature->n_dim();
+    if (basis_jet != 0) { delete [] basis_jet; }
+    basis_jet = new double[nq * ndofs * ndim];
+    cell->grad_shape(nq, quadrature->points, basis_jet);
+}
+
+void FE::update_jxw()
+{
+    const int nq = quadrature->n_points();
+    const int ndim = quadrature->n_dim();
+    for (int q = 0; q < nq; q++)
+    {
+        double jac[3][3];
+        const double wq = quadrature->getWeight(q);
+        const double *xq = &(quadrature->points[ndim*q]); // XXX: better way?
+        this->jxw[q] = wq * cell->jacobian(xq, jac);
+    }
+}
+
+

Modified: cs/cigma/trunk/src/FE.h
===================================================================
--- cs/cigma/trunk/src/FE.h	2008-12-10 02:12:29 UTC (rev 13492)
+++ cs/cigma/trunk/src/FE.h	2008-12-10 02:12:31 UTC (rev 13493)
@@ -19,7 +19,11 @@
     FE(const FE& other);
     ~FE();
 
-    void init(boost::shared_ptr<Quadrature> Q, boost::shared_ptr<Cell> cell);
+    //void init(boost::shared_ptr<Quadrature> Q, boost::shared_ptr<Cell> cell);
+    void init(boost::shared_ptr<Quadrature> Q);
+    void init_basis();
+    void init_basis_jet();
+    void update_jxw();
 
 public:
     boost::shared_ptr<Cell> cell;