[Cython/1D] 1D solver can be interrupted by ctrl-c

By calling a pure Python function each time OneDim::eval is called, we can
catch KeyboardInterrupt exceptions and abort the 1D solver loop, returning
control to Python.

Partially addresses Issue 93.

include/cantera/oneD/OneDim.h

@@ -12,6 +12,7 @@ namespace Cantera
 
 class MultiJac;
 class MultiNewton;
+class Func1;
 
 /**
  * Container class for multiple-domain 1D problems. Each domain is
@@ -215,6 +216,13 @@ public:
     }
     void saveStats();
 
+    //! Set a function that will be called every time #eval is called.
+    //! Can be used to provide keyboard interrupt support in the high-level
+    //! language interfaces.
+    void setInterrupt(Func1* interrupt) {
+        m_interrupt = interrupt;
+    }
+
 protected:
 
     void evalSSJacobian(doublereal* x, doublereal* xnew);
@@ -247,6 +255,9 @@ protected:
     // options
     int m_ss_jac_age, m_ts_jac_age;
 
+    //! Function called at the start of every call to #eval.
+    Func1* m_interrupt;
+
 private:
 
     // statistics

interfaces/cython/cantera/_cantera.pxd

@@ -436,7 +436,7 @@ cdef extern from "cantera/oneD/Sim1D.h":
         void showSolution() except +
         void setTimeStep(double, size_t, int*) except +
         void getInitialSoln() except +
-        void solve(int, cbool) except +
+        void solve(int, cbool) except +translate_exception
         void refine(int) except +
         void setRefineCriteria(size_t, double, double, double, double)
         void save(string, string, string, int) except +
@@ -452,6 +452,7 @@ cdef extern from "cantera/oneD/Sim1D.h":
         void setMaxTimeStep(double)
         void setGridMin(int, double) except +
         void setFixedTemperature(double)
+        void setInterrupt(CxxFunc1*) except +
 
 
 cdef extern from "wrappers.h":

interfaces/cython/cantera/func1.pyx

@@ -9,7 +9,7 @@ cdef double func_callback(double t, void* obj, void** err):
     """
     try:
         return (<Func1>obj).callable(t)
-    except Exception as e:
+    except BaseException as e:
         exc_type, exc_value = sys.exc_info()[:2]
 
         # Stash the exception info to prevent it from being garbage collected

interfaces/cython/cantera/interrupts.py

@@ -0,0 +1,7 @@
+def no_op(t):
+    """
+    This function does nothing. It is used as an interrupt in the 1D solver
+    C++ loop where a pure Python function is needed in order for
+    KeyboardInterrupt events to be captured.
+    """
+    return 0.0

interfaces/cython/cantera/onedim.pyx

@@ -1,4 +1,5 @@
 import csv
+import interrupts
 
 cdef class Domain1D:
     cdef CxxDomain1D* domain
@@ -445,6 +446,7 @@ cdef class Sim1D:
     cdef CxxSim1D* sim
     cdef readonly object domains
     cdef object _initialized
+    cdef Func1 interrupt
 
     def __cinit__(self, *args, **kwargs):
         self.sim = NULL
@@ -457,9 +459,21 @@ cdef class Sim1D:
 
         self.sim = new CxxSim1D(D)
         self.domains = tuple(domains)
-
+        self.set_interrupt(interrupts.no_op)
         self._initialized = False
 
+    def set_interrupt(self, f):
+        """
+        Set an interrupt function to be called each time that OneDim::eval is
+        called. The signature of *f* is `float f(float)`. The default
+        interrupt function is used to trap KeyboardInterrupt exceptions so
+        that `ctrl-c` can be used to break out of the C++ solver loop.
+        """
+        if not isinstance(f, Func1):
+            f = Func1(f)
+        self.interrupt = f
+        self.sim.setInterrupt(self.interrupt.func)
+
     def domain_index(self, dom):
         """
         Get the index of a domain, specified either by name or as a Domain1D

src/oneD/OneDim.cpp

@@ -2,6 +2,7 @@
 #include "cantera/oneD/MultiNewton.h"
 #include "cantera/oneD/OneDim.h"
 
+#include "cantera/numerics/Func1.h"
 #include "cantera/base/ctml.h"
 
 #include <fstream>
@@ -22,7 +23,7 @@ OneDim::OneDim()
       m_nd(0), m_bw(0), m_size(0),
       m_init(false),
       m_ss_jac_age(10), m_ts_jac_age(20),
-      m_nevals(0), m_evaltime(0.0)
+      m_nevals(0), m_interrupt(0), m_evaltime(0.0)
 {
     //writelog("OneDim default constructor\n");
     m_newt = new MultiNewton(1);
@@ -41,7 +42,7 @@ OneDim::OneDim(vector<Domain1D*> domains) :
     m_nd(0), m_bw(0), m_size(0),
     m_init(false),
     m_ss_jac_age(10), m_ts_jac_age(20),
-    m_nevals(0), m_evaltime(0.0)
+    m_nevals(0), m_interrupt(0), m_evaltime(0.0)
 {
     //writelog("OneDim constructor\n");
 
@@ -282,6 +283,9 @@ Domain1D* OneDim::pointDomain(size_t i)
 void OneDim::eval(size_t j, double* x, double* r, doublereal rdt, int count)
 {
     clock_t t0 = clock();
+    if (m_interrupt) {
+        m_interrupt->eval(m_nevals);
+    }
     fill(r, r + m_size, 0.0);
     fill(m_mask.begin(), m_mask.end(), 0);
     if (rdt < 0.0) {