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

[luis at geodynamics.org]

Author: luis
Date: 2008-10-15 02:08:18 -0700 (Wed, 15 Oct 2008)
New Revision: 13069

Added:
   cs/cigma/trunk/src/py_Quadrature.cpp
Log:
Python wrapper for the C++ Quadrature class

Added: cs/cigma/trunk/src/py_Quadrature.cpp
===================================================================
--- cs/cigma/trunk/src/py_Quadrature.cpp	                        (rev 0)
+++ cs/cigma/trunk/src/py_Quadrature.cpp	2008-10-15 09:08:18 UTC (rev 13069)
@@ -0,0 +1,162 @@
+#include <iostream>
+#include "Quadrature.h"
+#include "numpy_util.h"
+
+using namespace cigma;
+using namespace boost::python;
+
+struct pyQuadrature : Quadrature
+{
+    pyQuadrature() : na_points(0), na_weights(0)
+    {
+        PyErr_SetString(PyExc_RuntimeError, "Constructor called without arguments");
+        throw_error_already_set();
+    }
+
+    pyQuadrature(int npts, int ndim) : Quadrature(), na_points(0), na_weights(0)
+    {
+        std::vector<int> point_dims;
+        point_dims.push_back(npts);
+        point_dims.push_back(ndim);
+        na_points = numpy_util::makeArray(point_dims, PyArray_DOUBLE);
+        na_weights = numpy_util::makeArray(npts, PyArray_DOUBLE);
+        this->npts = npts;
+        this->ndim = ndim;
+        this->points = (double *)numpy_util::data(na_points);
+        this->weights = (double *)numpy_util::data(na_weights);
+        for (int q = 0; q < npts; q++)
+        {   
+            this->weights[q] = 0;
+            for (int j = 0; j < ndim; j++)
+            {   
+                this->points[ndim*q + j] = 0;
+            }
+        }
+    }
+
+    ~pyQuadrature()
+    {
+        //std::cout << "Calling ~pyQuadrature()" << std::endl;
+        npts = 0;
+        ndim = 0;
+        points = 0;
+        weights = 0;
+    }
+
+    double getPoint(int q, int j) const
+    {
+        if ((q >= n_points()) || (q < 0))
+        {
+            PyErr_SetString(PyExc_IndexError, "First index is out of bounds");
+            throw_error_already_set();
+        }
+        if ((j >= n_dim()) || (j < 0))
+        {
+            PyErr_SetString(PyExc_IndexError, "Second index is out of bounds");
+            throw_error_already_set();
+        }
+        return Quadrature::getPoint(q,j);
+    }
+
+    double getWeight(int q) const
+    {
+        if ((q >= n_points()) || (q < 0))
+        {
+            PyErr_SetString(PyExc_IndexError, "Index is out of bounds");
+            throw_error_already_set();
+        }
+        return Quadrature::getWeight(q);
+    }
+
+    void setPoint(int q, int j, double val)
+    {
+        if ((q >= n_points()) || (q < 0))
+        {
+            PyErr_SetString(PyExc_IndexError, "First index is out of bounds");
+            throw_error_already_set();
+        }
+        if ((j >= n_dim()) || (j < 0))
+        {
+            PyErr_SetString(PyExc_IndexError, "Second index is out of bounds");
+            throw_error_already_set();
+        }
+        Quadrature::setPoint(q, j, val);
+    }
+
+    void setWeight(int q, double val)
+    {
+        if ((q >= n_points()) || (q < 0))
+        {
+            PyErr_SetString(PyExc_IndexError, "Index is out of bounds");
+            throw_error_already_set();
+        }
+        Quadrature::setWeight(q, val);
+    }
+
+    numeric::array getPoints()
+    {
+        //return numpy_util::makeArray<double>(weights, n_points());
+        return na_points;
+    }
+
+    numeric::array getWeights()
+    {
+        //return numpy_util::makeArray<double>(points, n_points() * n_dim());
+        return na_weights;
+    }
+
+    void setPoints(numeric::array x)
+    {
+        numpy_util::check_type(x, PyArray_DOUBLE);
+        numpy_util::check_rank(x, 2);
+        numpy_util::check_size(x, npts * ndim);
+        double *data = (double *)numpy_util::data(x);
+        for (int q = 0; q < npts; q++)
+        {
+            for (int j = 0; j < ndim; j++)
+            {
+                int n = ndim*q + j;
+                points[n] = data[n];
+            }
+        }
+    }
+
+    void setWeights(numeric::array w)
+    {
+        numpy_util::check_type(w, PyArray_DOUBLE);
+        numpy_util::check_rank(w, 1);
+        numpy_util::check_size(w, npts);
+        double *data = (double *)numpy_util::data(w);
+        for (int q = 0; q < npts; q++)
+        {
+            weights[q] = data[q];
+        }
+    }
+
+private:
+    numeric::array na_points;
+    numeric::array na_weights;
+};
+
+
+void export_Quadrature()
+{
+    using namespace boost::python;
+
+    typedef pyQuadrature pq;
+    class_<pyQuadrature, boost::noncopyable>("Quadrature", init<int,int>())
+        
+        .def("n_points", &pq::n_points)
+        .def("n_dim", &pq::n_dim)
+        
+        .def("getPoint", &pq::getPoint)
+        .def("getWeight", &pq::getWeight)
+
+        .def("setPoint", &pq::setPoint)
+        .def("setWeight", &pq::setWeight)
+
+        .add_property("points", &pq::getPoints, &pq::setPoints)
+        .add_property("weights", &pq::getWeights, &pq::setWeights)
+        ;
+}
+