[cig-commits] r11722 - cs/benchmark/cigma/trunk/src

[luis at geodynamics.org]

Author: luis
Date: 2008-04-02 11:01:39 -0700 (Wed, 02 Apr 2008)
New Revision: 11722

Modified:
   cs/benchmark/cigma/trunk/src/QuadratureRule.cpp
   cs/benchmark/cigma/trunk/src/QuadratureRule.h
Log:
Merged code from child FE class back into QuadratureRule


Modified: cs/benchmark/cigma/trunk/src/QuadratureRule.cpp
===================================================================
--- cs/benchmark/cigma/trunk/src/QuadratureRule.cpp	2008-04-02 18:01:38 UTC (rev 11721)
+++ cs/benchmark/cigma/trunk/src/QuadratureRule.cpp	2008-04-02 18:01:39 UTC (rev 11722)
@@ -11,11 +11,16 @@
     meshPart = 0;
     points = 0;
     jxw = 0;
+    basis_tab = 0;
+    basis_jet = 0;
 }
 
 QuadratureRule::~QuadratureRule()
 {
+    if (points != 0) delete points;
     if (jxw != 0) delete [] jxw;
+    if (basis_tab != 0) delete [] basis_tab;
+    if (basis_jet != 0) delete [] basis_jet;
 }
 
 
@@ -24,11 +29,13 @@
 void QuadratureRule::set_mesh(MeshPart *mesh)
 {
     assert(mesh != 0);
+    assert(mesh->cell != 0);
     this->meshPart = mesh;
 }
 
-void QuadratureRule::set_quadrature_points(QuadraturePoints *points)
+void QuadratureRule::set_quadrature(double *wts, double *pts, int npts, int ndim)
 {
+    /*
     assert(meshPart != 0);
     assert(meshPart->cell != 0);
 
@@ -36,10 +43,45 @@
     assert(points != 0);
     assert(points->n_points() > 0);
     assert(points->n_refdim() == meshPart->cell->n_celldim());
+    */
 
-    jxw = new double[points->n_points()];
+    points = new QuadraturePoints();
+    points->set_quadrature(wts, pts, npts, ndim);
+
+    jxw = new double[npts];
 }
 
+void QuadratureRule::initialize_basis_tab()
+{
+    assert(meshPart != 0);
+    assert(meshPart->cell != 0);
+    assert(points != 0);
+    assert(points->n_refdim() == meshPart->cell->n_celldim());
+    
+    const int nq = points->n_points();
+    const int ndofs = meshPart->cell->n_nodes();
+
+    // get shape function values at known quadrature points
+    basis_tab = new double[nq * ndofs];
+    meshPart->cell->shape(nq, points->qpts, basis_tab);
+}
+
+void QuadratureRule::initialize_basis_jet()
+{
+    assert(meshPart != 0);
+    assert(meshPart->cell != 0);
+    assert(points != 0);
+    assert(points->n_refdim() == meshPart->cell->n_celldim());
+    
+    const int nq = points->n_points();
+    const int ndofs = meshPart->cell->n_nodes();
+    const int ndim = points->n_dim();
+
+    // get shape function derivatives at known quadrature points
+    basis_jet = new double[nq * ndofs * ndim];
+    meshPart->cell->grad_shape(nq, points->qpts, basis_jet);
+}
+
 // ---------------------------------------------------------------------------
 
 void QuadratureRule::update_jxw()

Modified: cs/benchmark/cigma/trunk/src/QuadratureRule.h
===================================================================
--- cs/benchmark/cigma/trunk/src/QuadratureRule.h	2008-04-02 18:01:38 UTC (rev 11721)
+++ cs/benchmark/cigma/trunk/src/QuadratureRule.h	2008-04-02 18:01:39 UTC (rev 11722)
@@ -20,7 +20,9 @@
 
 public:
     void set_mesh(MeshPart *mesh);
-    virtual void set_quadrature_points(QuadraturePoints *pts);
+    void set_quadrature(double *wts, double *pts, int npts, int ndim);
+    void initialize_basis_tab();
+    void initialize_basis_jet();
 
 public:
     void select_cell(int e);
@@ -33,7 +35,9 @@
 public:
     MeshPart *meshPart;
     QuadraturePoints *points;
-    double *jxw;
+    double *jxw;        // [nq x 1]
+    double *basis_tab;  // [nq x ndofs]
+    double *basis_jet;  // [nq x ndofs x celldim]
 };