parameter setting for full 3000 second run

Author: Daniel Ruprecht

.travis.yml

@@ -0,0 +1,7 @@
+language: python
+python:
+  - "2.7"
+# command to install dependencies
+install: "pip install -r requirements.txt"
+#
+script: py.test -v tests/

examples/boussinesq_2d_imex/playground.py

@@ -28,22 +28,22 @@
 
     # This comes as read-in for the level class
     lparams = {}
-    lparams['restol'] = 1E-14
+    lparams['restol'] = 1E-8
 
     swparams = {}
     swparams['collocation_class'] = collclass.CollGaussLobatto
     swparams['num_nodes'] = 3
     swparams['do_LU'] = True
 
     sparams = {}
-    sparams['maxiter'] = 4
+    sparams['maxiter'] = 12
 
     # setup parameters "in time"
     t0     = 0
-    #Tend   = 3000
-    #Nsteps =  500
-    Tend   = 30
-    Nsteps =  5
+    Tend   = 3000
+    Nsteps =  500
+    #Tend   = 30
+    #Nsteps =  5
     dt = Tend/float(Nsteps)
 
     # This comes as read-in for the problem class
@@ -58,7 +58,7 @@
     pparams['order_upw'] = [5, 1]
     pparams['gmres_maxiter'] = [50, 50]
     pparams['gmres_restart'] = [20, 20]
-    pparams['gmres_tol']     = [1e-6, 1e-6]
+    pparams['gmres_tol']     = [1e-8, 1e-8]
 
     # This comes as read-in for the transfer operations
     tparams = {}

requirements.txt

@@ -0,0 +1,3 @@
+py.test
+numpy
+scipy

tests/test_collocation.py

@@ -0,0 +1,84 @@
+import pySDC.Collocation
+from pySDC.CollocationClasses import *
+import numpy as np
+import pytest
+
+# py.test excludes classes with a constructor by default, so define parameter here
+t_start = 0.0
+t_end   = 1.0
+classes = [ ["CollGaussLegendre", 2, 2]]
+#classes = [ ["CollGaussLegendre", 2, 12], ["CollGaussLobatto", 2, 12], ["CollGaussRadau_Right", 2, 12] ]
+
+class TestCollocation:
+
+  # TEST 1:
+  # Check that the quadrature rule integrates polynomials up to order p-1 exactly
+  # -----------------------------------------------------------------------------
+  def test_1(self):
+    for type in classes:
+      for M in range(type[1],type[2]+1):
+        coll = getattr(pySDC.CollocationClasses, type[0])(M, t_start, t_end)
+        
+        # some basic consistency tests
+        assert np.size(coll.nodes)==np.size(coll.weights), "For node type " + type[0] + ", number of entries in nodes and weights is different"
+        assert np.size(coll.nodes)==M, "For node type " + type[0] + ", requesting M nodes did not produce M entries in nodes and weights"
+
+
+        # generate random set of polynomial coefficients
+        poly_coeff = np.random.rand(coll.order-1)
+        # evaluate polynomial at collocation nodes
+        poly_vals  = np.polyval(poly_coeff, coll.nodes)
+        # use python's polyint function to compute anti-derivative of polynomial
+        poly_int_coeff = np.polyint(poly_coeff)
+        # Compute integral from 0.0 to 1.0
+        int_ex = np.polyval(poly_int_coeff, t_end) - np.polyval(poly_int_coeff, t_start)
+        # use quadrature rule to compute integral
+        int_coll = coll.evaluate(coll.weights, poly_vals)
+        # For large values of M, substantial differences from different round of error have to be considered
+        assert abs(int_ex - int_coll) < 1e-10, "For node type " + type[0] + ", failed to integrate polynomial of degree " + str(coll.order-1) + " exactly. Error: %5.3e" % abs(int_ex - int_coll)
+
+
+  # TEST 2:
+  # Check that the Qmat entries are equal to the sum of Smat entries
+  # ----------------------------------------------------------------
+  def test_2(self):
+    for type in classes:
+      for M in range(type[1],type[2]+1):
+        coll = getattr(pySDC.CollocationClasses, type[0])(M, t_start, t_end)
+        Q = coll.Qmat
+        Q = Q[1:,1:]
+        S = coll.Smat    
+        S = S[1:,1:]
+        #print Q
+        #print S
+        assert np.shape(Q) == np.shape(S), "For node type " + type[0] + ", Qmat and Smat have different shape"
+        shape = np.shape(Q)
+        assert shape[0] == shape[1], "For node type " + type[0] + ", Qmat / Smat are not quadratic"
+        for i in range(0,M):
+            Ssum = np.sum(S[0:i,:], axis=0)
+            #print Ssum
+            #print Q[i,:]
+            # ...the matrices do not have size MxM, but rather (M+1)x(M+1)... how does the summation property look like here
+
+  # TEST 3:
+  # Check that the partial quadrature rules from Qmat entries have order equal to number of nodes M
+  # -----------------------------------------------------------------------------------------------
+  def test_3(self):
+    for type in classes:
+      for M in range(type[1],type[2]+1):
+        coll = getattr(pySDC.CollocationClasses, type[0])(M, t_start, t_end)
+        Q = coll.Qmat
+        Q = Q[1:,1:]
+        # as in TEST 1, create and integrate a polynomial with random coefficients, but now of degree M-1
+        poly_coeff = np.random.rand(M-1)
+        poly_vals  = np.polyval(poly_coeff, coll.nodes)
+        poly_int_coeff = np.polyint(poly_coeff)
+        #print Q
+        for i in range(0,M):
+            int_ex = np.polyval(poly_int_coeff, coll.nodes[i]) - np.polyval(poly_int_coeff, t_start)
+            print np.shape(Q)
+            #int_coll = np.dot(poly_vals, Q[i,:])
+
+    # TEST 4:
+    # Check that the quadrature rules for [tau_m, tau_m+1] defined by Smat have order equal to number of nodes M
+    # ----------------------------------------------------------------------------------------------------------