Ports odespy to python 3 (fix #11)

odespy/PyDSTool.py

@@ -1,7 +1,9 @@
 # Author: Liwei Wang
 """
 """
-from solvers import *
+from __future__ import print_function
+from __future__ import absolute_import
+from .solvers import *
 import numpy as np
 
 class Pyds(Solver):
@@ -24,10 +26,10 @@ class Pyds(Solver):
 
     def initialize(self):
         try:
-            import PyDSTool
+            from . import PyDSTool
         except ImportError:
-            raise ImportError,'''
-        PyDSTool is not installed - required for solvers from PyDSTool'''
+            raise ImportError('''
+        PyDSTool is not installed - required for solvers from PyDSTool''')
 
     def solve(self, time_points, terminate=None):
         # Common parts as superclass
@@ -45,7 +47,7 @@ def solve(self, time_points, terminate=None):
         # through Python dictionaries with string keys.
 
         # Start setting for PyDSTool
-        import PyDSTool
+        from . import PyDSTool
         neq, f, u0 = self.neq, self.f, self.U0
 
         # Initialize variables as trajectories in PyDSTOOL
@@ -148,7 +150,7 @@ def initialize_for_solve(self):
     method = Vode_pyds(f)
     method.set_initial_condition([0.,1.])
     u,t = method.solve(np.linspace(0.,10.,50))
-    print u
+    print(u)
     import scitools.std as st
     st.plot(t,u[:,0])
-    print max(u[:,0]-np.sin(t))
+    print(max(u[:,0]-np.sin(t)))

odespy/RungeKutta.py

@@ -1,4 +1,6 @@
-from solvers import Solver, Adaptive
+from __future__ import print_function
+from __future__ import absolute_import
+from .solvers import Solver, Adaptive
 import numpy as np
 
 def _calculate_order_1_level(coefficients):
@@ -134,7 +136,7 @@ def middle(x,y,z):    # Auxilary function
         k = np.zeros((k_len, self.neq), self.dtype)  # intern stages
 
         if self.verbose > 0:
-            print 'advance solution in [%s, %s], h=%g' % (t_n, t_next, h)
+            print('advance solution in [%s, %s], h=%g' % (t_n, t_next, h))
 
         # Loop until next time point is reached
         while (abs(t - t_n) < abs(t_next - t_n)):
@@ -150,7 +152,7 @@ def middle(x,y,z):    # Auxilary function
 
             self.info['rejected'] += 1  # reduced below if accepted
             if self.verbose > 0:
-                print '  u(t=%g)=%g: ' % (t+h, u_new),
+                print('  u(t=%g)=%g: ' % (t+h, u_new), end=' ')
 
             # local error between 2 levels
             error = h*np.abs(np.dot(factors_error, k))
@@ -171,18 +173,18 @@ def middle(x,y,z):    # Auxilary function
                 self.info['rejected'] -= 1
 
                 if self.verbose > 0:
-                    print 'accepted, ',
+                    print('accepted, ', end=' ')
             else:
                 if self.verbose > 0:
-                    print 'rejected, ',
+                    print('rejected, ', end=' ')
 
             if self.verbose > 0:
-                print 'err=%s, ' % str(error),
+                print('err=%s, ' % str(error), end=' ')
                 if hasattr(self, 'u_exact') and callable(self.u_exact):
-                    print 'exact-err=%s, ' % \
-                          (np.asarray(self.u_exact(t+h))-u_new),
+                    print('exact-err=%s, ' % \
+                          (np.asarray(self.u_exact(t+h))-u_new), end=' ')
                 if h <= self.min_step:
-                    print 'h=min_step!! ',
+                    print('h=min_step!! ', end=' ')
 
 
            # Replace 0 values by 1e-16 since we will divide by error
@@ -209,7 +211,7 @@ def middle(x,y,z):    # Auxilary function
             h = min(h, t_next - t_intermediate[-1])
 
             if self.verbose > 0:
-                print 'new h=%g' % h
+                print('new h=%g' % h)
 
             if h == 0:
                 break
@@ -367,16 +369,16 @@ def validate_data(self):
         # Check for dimension of user-defined butcher table.
         array_shape = self.butcher_tableau.shape
         if len(array_shape) is not 2:
-            raise ValueError,'''
-        Illegal input! Your input butcher_tableau should be a 2d-array!'''
+            raise ValueError('''
+        Illegal input! Your input butcher_tableau should be a 2d-array!''')
         else:
             m,n = array_shape
             if m not in (n, n + 1):
-                raise ValueError, '''\
+                raise ValueError('''\
         The dimension of 2d-array <method_yours_array> should be:
         1. Either (n, n), --> For 1-level RungeKutta methods
         2. Or (n+1, n),   --> For 2-levels RungeKutta methods
-        The shape of your input array is (%d, %d).''' % (m,n)
+        The shape of your input array is (%d, %d).''' % (m,n))
         self._butcher_tableau = self.butcher_tableau
 
         # Check for user-defined order,
@@ -393,19 +395,19 @@ def validate_data(self):
             if array_shape[0] == array_shape[1] + 1:
                 # 2-level RungeKutta methods
                 if type(self.method_order) is int:
-                    raise ValueError, error_2level
+                    raise ValueError(error_2level)
                 try:
                     order1, order2 = self.method_order
                     if abs(order1-order2) != 1 or \
                             order1 < 1 or order2 < 1:
-                        raise ValueError, error_2level
+                        raise ValueError(error_2level)
                 except:
-                    raise ValueError,error_2level
+                    raise ValueError(error_2level)
             else:
                 # 1-level RungeKutta methods
                 if type(self.method_order) is not int or \
                         self.method_order < 1:
-                    raise ValueError,error_1level
+                    raise ValueError(error_1level)
             self._method_order = self.method_order
 
         else:   # method_order is not specified
@@ -418,14 +420,14 @@ def validate_data(self):
         for i in range(1,array_shape[1] - 1):
             if not np.allclose(self.butcher_tableau[i][0],\
                                sum(self.butcher_tableau[i][1:])):
-                raise ValueError, '''
+                raise ValueError('''
         Inconsistent data in Butcher_Tableau!
         In each lines of stage-coefficients, first number should be
         equal to the sum of other numbers.
         That is, for a butcher_table with K columns,
             a[i][0] == a[i][1] + a[i][2] + ... + a[i][K - 1]
             where 1 <= i <= K - 1
         Your input for line %d is :%s
-        ''' % (i,str(self.butcher_tableau[i]))
+        ''' % (i,str(self.butcher_tableau[i])))
 
         return True

odespy/__init__.py

@@ -5,6 +5,7 @@
 supported.  A wide range of numerical methods for ODEs are offered:
 
 """
+from __future__ import absolute_import
 
 # Insert tutorial from ../doc/src/odespy/odespy.rst
 
@@ -1255,13 +1256,13 @@ def f(u, t):
 with :math:`f(u,t)` implemented in Fortran.
 '''
 
-from solvers import *
-from RungeKutta import *
-from rkc import *
-from rkf45 import *
-from odepack import *
-from radau5 import *
-import problems
+from .solvers import *
+from .RungeKutta import *
+from .rkc import *
+from .rkf45 import *
+from .odepack import *
+from .radau5 import *
+from . import problems
 
 # Update doc strings with common info
 class_, doc_str, classname = None, None, None
@@ -1283,7 +1284,7 @@ def f(u, t):
 
 # Do not pollute namespace
 del class_, doc_str, classname, classnames, toc, typeset_toc, \
-    table_of_parameters, name, obj, inspect
+    table_of_parameters, inspect
 
 if __name__ == '__main__':
     from os.path import join

odespy/demos/demo_Lsodi_1.py

@@ -9,6 +9,7 @@
 This example is the typical usage of Lsodi with 
 user-supplied functions composed in Python.
 """
+from __future__ import print_function
 from odespy import *
 import scitools.std as st
 import numpy as np
@@ -52,7 +53,7 @@ def jac(u, t, s):
 u,t = m.solve(time_points)
 st.plot(t, u[:,0], 'b-', title="Lsodi with Python functions",
         legend="with res, adda, ydoti & jac", hold="on")
-print 'Max error for test case 1 is %g' % max(u[-1] - exact_final)
+print('Max error for test case 1 is %g' % max(u[-1] - exact_final))
 
 # Test case 2: Lsodi, with res, ydoti & adda
 m = method(res=res, rtol=rtol, atol=atol, ydoti=ydoti,
@@ -61,6 +62,6 @@ def jac(u, t, s):
 u,t = m.solve(time_points)
 st.plot(t, u[:,0], 'r*', title="Lsodi with Python functions",
         legend="with res, adda & ydoti", hold="on")
-print 'Max error for test case 1 is %g' % max(u[-1] - exact_final)
+print('Max error for test case 1 is %g' % max(u[-1] - exact_final))
 
 

odespy/demos/demo_Lsodi_1_fortran.py

@@ -8,6 +8,7 @@
 user-supplied functions composed in Fortran code.
 
 """
+from __future__ import print_function
 from odespy import *
 import scitools.std as st
 import numpy as np
@@ -83,7 +84,7 @@
 u,t = m.solve(time_points)
 st.plot(t, u[:,0], 'r-', title="Lsodi with Fortran subroutines",
         legend="with res, adda, ydoti & jac", hold="on")
-print 'Max error for test case 1 is %g' % max(u[-1] - exact_final)
+print('Max error for test case 1 is %g' % max(u[-1] - exact_final))
 
 # Test case 2: Lsodi, with res, ydoti & adda
 m = method(res_f77=res_f77, rtol=rtol, atol=atol, ydoti=ydoti,
@@ -92,6 +93,6 @@
 u,t = m.solve(time_points)
 st.plot(t, u[:,0], 'g*', title="Lsodi with Fortran subroutines",
         legend="with res, adda & ydoti", hold="on")
-print 'Max error for test case 2 is %g' % max(u[-1] - exact_final)
+print('Max error for test case 2 is %g' % max(u[-1] - exact_final))
 
 os.remove('callback.so')

odespy/demos/demo_Lsodi_2.py

@@ -27,6 +27,7 @@
         0.5000E+00  0.0000E+00  0.0000E+00  0.0000E+00  0.0000E+00
 
 """
+from __future__ import print_function
 from odespy import *
 #import scitools.basics,easyviz as st
 import scitools.std as st
@@ -113,7 +114,7 @@ def jac_banded(u, t, s, ml, mu):   # Banded jacobian
 u,t = m.solve(time_points)
 st.plot(t, u[:,0], 'r-', title="Lsodi with Python functions",
         legend="with res, adda_full & jac_full", hold="on")
-print 'Max error with test case 1 is %g' % max(u[-1] - exact_final)
+print('Max error with test case 1 is %g' % max(u[-1] - exact_final))
 
 
 # Test case 2: Lsodi, with res & adda_full
@@ -122,7 +123,7 @@ def jac_banded(u, t, s, ml, mu):   # Banded jacobian
 u,t = m.solve(time_points)
 st.plot(t, u[:,0], 'b*', title="Lsodi with Python functions",
         legend="with res & adda_full", hold="on")
-print 'Max error with test case 2 is %g' % max(u[-1] - exact_final)
+print('Max error with test case 2 is %g' % max(u[-1] - exact_final))
 
 # Test case 3: Lsodi, with res, adda_banded, ml, mu, jac_banded
 m = method(res=res, rtol=rtol, atol=atol, 
@@ -133,7 +134,7 @@ def jac_banded(u, t, s, ml, mu):   # Banded jacobian
 u,t = m.solve(time_points)
 st.plot(t, u[:,0], 'go', title="Lsodi with Python functions",
         legend="with res, adda_banded, jac_banded, ml, mu", hold="on")
-print 'Max error with test case 3 is %g' % max(u[-1] - exact_final)
+print('Max error with test case 3 is %g' % max(u[-1] - exact_final))
 
 # Test case 4: Lsodi, with res, adda_banded, ml, mu
 m = method(res=res, rtol=rtol, atol=atol, 
@@ -143,7 +144,7 @@ def jac_banded(u, t, s, ml, mu):   # Banded jacobian
 u,t = m.solve(time_points)
 st.plot(t, u[:,0], 'y-', title="Lsodi with Python functions",
         legend="with res, adda_banded, ml, mu", hold="on")
-print 'Max error with test case 4 is %g' % max(u[-1] - exact_final)
+print('Max error with test case 4 is %g' % max(u[-1] - exact_final))
 
 
 

odespy/demos/demo_Lsodis_1.py

@@ -26,6 +26,7 @@
      = r4d*(y(nm1)**2-y(1)**2)+eodsq*(y(1)-2*y(NEQ)+y(nm1))
      where r4d = 1/(4*delx), eodsq = eta/delx**2 and nm1 = NEQ-1.
 """
+from __future__ import print_function
 from odespy import *
 import numpy as np
 
@@ -73,11 +74,11 @@ def jac(y, t, s, j, ia, ja):
 m = method(res=res, adda_lsodis=adda, atol=atol, rtol=rtol, jac_lsodis=jac)
 m.set_initial_condition(u0)
 y, t = m.solve(time_points)
-print 'Max error for test case 1 is %g' % max(y[-1] - exact_final)
+print('Max error for test case 1 is %g' % max(y[-1] - exact_final))
 
 # Test case 2:  With res & adda
 m = method(res=res, adda_lsodis=adda, atol=atol, rtol=rtol,lrw=4000,liw=100)
 m.set_initial_condition(u0)
 y, t = m.solve(time_points)
-print 'Max error for test case 2 is %g' % max(y[-1] - exact_final)
+print('Max error for test case 2 is %g' % max(y[-1] - exact_final))
 

odespy/demos/demo_Lsodis_2.py

@@ -19,6 +19,7 @@
     An ODE system is generated by a simplified Galerkin treatment of the spatial
     variable x.
 """
+from __future__ import print_function
 from odespy import *
 import numpy as np
 
@@ -69,25 +70,25 @@ def jac(y, t, s, j, ia, ja):
            ia=ia, ja=ja, ic=ic, jc=jc)
 m.set_initial_condition(u0)
 y, t = m.solve(time_points)
-print 'Max error for test case 1 is %g' % max(y[-1] - exact_final)
+print('Max error for test case 1 is %g' % max(y[-1] - exact_final))
 
 # Test case 2: with res, adda, ia, ja & jac
 m = method(res=res, adda_lsodis=adda, atol=atol, rtol=rtol, jac_lsodis=jac,
            ia=ia, ja=ja)
 m.set_initial_condition(u0)
 y, t = m.solve(time_points)
-print 'Max error for test case 2 is %g' % max(y[-1] - exact_final)
+print('Max error for test case 2 is %g' % max(y[-1] - exact_final))
 
 # Test case 3: with res, adda & jac
 m = method(res=res, adda_lsodis=adda, atol=atol, rtol=rtol, jac_lsodis=jac)
 m.set_initial_condition(u0)
 y, t = m.solve(time_points)
-print 'Max error for test case 3 is %g' % max(y[-1] - exact_final)
+print('Max error for test case 3 is %g' % max(y[-1] - exact_final))
 
 # Test case 4:  With res & adda
 m = method(res=res, adda_lsodis=adda, atol=atol, rtol=rtol,lrw=4000,liw=100)
 m.set_initial_condition(u0)
 y, t = m.solve(time_points)
-print 'Max error for test case 4 is %g' % max(y[-1] - exact_final)
+print('Max error for test case 4 is %g' % max(y[-1] - exact_final))
 
 

odespy/demos/demo_Lsoibt_2.py

@@ -32,6 +32,7 @@
 Terms involving boundary values (subscripts 0 or 100) are dropped
 from the equations for k = 1 and k = 99 above.
 """
+from __future__ import print_function
 from odespy import *
 import scitools.std as st
 import numpy as np
@@ -209,7 +210,7 @@ def jac(y, t, s):
 u_final = u[-1].reshape(99,3)
 u1, u2, u3 = u_final[:, 0], u_final[:, 1], u_final[:, 2]
 max_error = max(max(u1 - u1_exact), max(u2 - u2_exact), max(u3 - u3_exact))
-print 'Max error with Test case 1 is %g' % max_error
+print('Max error with Test case 1 is %g' % max_error)
 
 # Test case 2: Lsoibt, with res, adda, mb, nb
 m = method(rtol=rtol, atol=atol, res=res, adda_lsoibt=adda, 
@@ -222,4 +223,4 @@ def jac(y, t, s):
 u_final = u[-1].reshape(99,3)
 u1, u2, u3 = u_final[:, 0], u_final[:, 1], u_final[:, 2]
 max_error = max(max(u1 - u1_exact), max(u2 - u2_exact), max(u3 - u3_exact))
-print 'Max error with Test case 2 is %g' % max_error
+print('Max error with Test case 2 is %g' % max_error)

odespy/demos/demo_MyRungeKutta.py

@@ -4,6 +4,7 @@
 This example intends to show users how to apply MyRungeKutta to define 
 own RungKutta solvers.
 """
+from __future__ import print_function
 from odespy import *
 #import scitools.basics,easyviz as st
 import scitools.std as st
@@ -45,7 +46,7 @@ def f(u,t):
     method.set_initial_condition(u0)
     u,t = method.solve(time_points)
     error = abs((u[-1] - np.exp(-1.))/np.exp(-1.))
-    print 'Error is %g with solver %s' % (error, m)
+    print('Error is %g with solver %s' % (error, m))