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

[luis at geodynamics.org]

Author: luis
Date: 2008-12-17 02:33:15 -0800 (Wed, 17 Dec 2008)
New Revision: 13756

Modified:
   cs/cigma/trunk/src/Quadrature.cpp
   cs/cigma/trunk/src/Quadrature.h
Log:
Added two methods to cigma::Quadrature class

 * Quadrature::setFromCombinedData()
 * Quadrature::check_points()

Modified: cs/cigma/trunk/src/Quadrature.cpp
===================================================================
--- cs/cigma/trunk/src/Quadrature.cpp	2008-12-17 10:33:14 UTC (rev 13755)
+++ cs/cigma/trunk/src/Quadrature.cpp	2008-12-17 10:33:15 UTC (rev 13756)
@@ -17,9 +17,12 @@
 shared_ptr<Quadrature> Quadrature::default_rule(Cell::type cell_type)
 {
     TRI_LOG_STR("Quadrature::default_rule()");
-    TRI_LOG(Cell::type2string(cell_type));
 
+    string cellname = Cell::type2string(cell_type);
+    TRI_LOG(cellname);
+
     shared_ptr<Quadrature> Q;
+
     if (cell_type == Cell::TET4)
     {
         // tet4_qr(3)
@@ -271,6 +274,23 @@
     }
 }
 
+void Quadrature::setFromCombinedData(const cigma::array<double>& wx)
+{
+    assert(wx.n_points() > 0);
+    assert(wx.n_dim() > 1);
+
+    this->reinit(wx.n_points(), wx.n_dim() - 1);
+
+    for (int q = 0; q < npts; q++)
+    {
+        this->setWeight(q, wx(q,0));
+        for (int j = 0; j < ndim; j++)
+        {
+            this->setPoint(q, j, wx(q,1+j));
+        }
+    }
+}
+
 // ----------------------------------------------------------------------------
 
 
@@ -311,7 +331,7 @@
             if (cell_type == Cell::NONE)
             {
                 string cell_name;
-                int status = h5->readAttrString(loc.c_str(), "CellType", cell_name);
+                int status = h5->getStringAttr(loc.c_str(), "CellType", cell_name);
                 TRI_LOG(status);
                 if (status == 0)
                 {
@@ -364,23 +384,12 @@
                    << p_quad << "'" << std::ends;
             throw cigma::Exception("Quadrature::setPath", stream.str());
         }
-
-        const int q_npts = wx->n_points();
-        const int q_ndim = -1 + wx->n_dim();
-
-        this->reinit(q_npts, q_ndim);
-
-        for (int q = 0; q < npts; q++)
+        else
         {
-            this->setWeight(q, wx->data(q, 0));
-            for (int j = 0; j < ndim; j++)
-            {
-                this->setPoint(q, j, wx->data(q, 1+j));
-            }
+            this->setFromCombinedData(*wx);
+            delete wx;
         }
-        delete wx;
 
-        // check points
         bool all_points_in_cell = points_in_cell();
         TRI_LOG(all_points_in_cell);
         if (!all_points_in_cell)
@@ -467,15 +476,7 @@
         shared_ptr<Quadrature> Q = Quadrature::default_rule(cell_type);
         this->reinit(Q->npts, Q->ndim);
         this->setData(Q->points, Q->weights);
-
-        // check points
-        if (!points_in_cell())
-        {
-            std::ostringstream stream;
-            stream << "Quadrature points not compatible with cell ";
-            stream << "'" << cell_name << "'" << std::ends;
-            throw cigma::Exception("Quadrature::setCell", stream.str());
-        }
+        this->check_points();
     }
     else
     {
@@ -519,4 +520,16 @@
     return true;
 }
 
+void Quadrature::check_points() const
+{
+    if (!points_in_cell())
+    {
+        std::ostringstream stream;
+        stream << "Quadrature points not combatible with cell '"
+               << Cell::type2string(cell_type) << "'"
+               << std::ends;
+        throw cigma::Exception("Quadrature::check_points", stream.str());
+    }
+}
+
 // ----------------------------------------------------------------------------

Modified: cs/cigma/trunk/src/Quadrature.h
===================================================================
--- cs/cigma/trunk/src/Quadrature.h	2008-12-17 10:33:14 UTC (rev 13755)
+++ cs/cigma/trunk/src/Quadrature.h	2008-12-17 10:33:15 UTC (rev 13756)
@@ -5,8 +5,9 @@
 #include <boost/shared_ptr.hpp>
 
 #include "Cell.h"
+#include "core_array.h"
+#include "core_args.h"
 #include "DataPath.h"
-#include "core_args.h"
 
 namespace cigma
 {
@@ -30,10 +31,13 @@
     double getPoint(int q, int j) const;
     double getWeight(int q) const;
     Cell::type getCellType() const;
+
     bool points_in_cell() const;
+    void check_points() const;
 
     void reinit(int npts, int ndim);
-    void setData(double *points, double *weights); // XXX: rename to init(double*,double*)
+    void setData(double *points, double *weights);
+    void setFromCombinedData(const cigma::array<double>& wx);
 
     void setPoint(int q, int j, double val);
     void setWeight(int q, double val);