/[escript]/trunk/escriptcore/test/python/test_linearPDEs.py
ViewVC logotype

Diff of /trunk/escriptcore/test/python/test_linearPDEs.py

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

trunk/escript/test/python/test_linearPDEs.py revision 2548 by jfenwick, Mon Jul 20 06:20:06 2009 UTC trunk/escriptcore/test/python/test_linearPDEs.py revision 5148 by caltinay, Mon Sep 15 01:25:23 2014 UTC
# Line 1  Line 1 
1    # -*- coding: utf-8 -*-
2    
3  ########################################################  ##############################################################################
4  #  #
5  # Copyright (c) 2003-2009 by University of Queensland  # Copyright (c) 2003-2014 by University of Queensland
6  # Earth Systems Science Computational Center (ESSCC)  # http://www.uq.edu.au
 # http://www.uq.edu.au/esscc  
7  #  #
8  # Primary Business: Queensland, Australia  # Primary Business: Queensland, Australia
9  # Licensed under the Open Software License version 3.0  # Licensed under the Open Software License version 3.0
10  # http://www.opensource.org/licenses/osl-3.0.php  # http://www.opensource.org/licenses/osl-3.0.php
11  #  #
12  ########################################################  # Development until 2012 by Earth Systems Science Computational Center (ESSCC)
13    # Development 2012-2013 by School of Earth Sciences
14    # Development from 2014 by Centre for Geoscience Computing (GeoComp)
15    #
16    ##############################################################################
17    
18  __copyright__="""Copyright (c) 2003-2008 by University of Queensland  __copyright__="""Copyright (c) 2003-2014 by University of Queensland
19  Earth Systems Science Computational Center (ESSCC)  http://www.uq.edu.au
 http://www.uq.edu.au/esscc  
20  Primary Business: Queensland, Australia"""  Primary Business: Queensland, Australia"""
21  __license__="""Licensed under the Open Software License version 3.0  __license__="""Licensed under the Open Software License version 3.0
22  http://www.opensource.org/licenses/osl-3.0.php"""  http://www.opensource.org/licenses/osl-3.0.php"""
# Line 22  __url__="https://launchpad.net/escript-f Line 25  __url__="https://launchpad.net/escript-f
25  """  """
26  Test suite for linearPDEs class  Test suite for linearPDEs class
27    
 The tests must be linked with a Domain class object in the setUp method:  
   
    from esys.finley import Rectangle  
    class Test_LinearPDEOnFinley(Test_LinearPDE):  
        def setUp(self):  
            self.domain = Rectangle(10,10,2)  
        def tearDown(self):  
            del self.domain  
    suite = unittest.TestSuite()  
    suite.addTest(unittest.makeSuite(Test_LinearPDEOnFinley))  
    unittest.TextTestRunner(verbosity=2).run(suite)  
   
 @var __author__: name of author  
 @var __copyright__: copyrights  
 @var __license__: licence agreement  
 @var __url__: url entry point on documentation  
 @var __version__: version  
 @var __date__: date of the version  
28  """  """
29    
30  __author__="Lutz Gross, l.gross@uq.edu.au"  __author__="Lutz Gross, l.gross@uq.edu.au"
31    
32  from esys.escript.util import Lsup,kronecker,interpolate,whereZero, outer, swap_axes  from esys.escript.util import Lsup,kronecker,interpolate,whereZero, outer, swap_axes
33  from esys.escript import Function,FunctionOnBoundary,FunctionOnContactZero,Solution,ReducedSolution,Vector,ContinuousFunction,Scalar, ReducedFunction,ReducedFunctionOnBoundary,ReducedFunctionOnContactZero,Data, Tensor4, Tensor  from esys.escript import Function,FunctionOnBoundary,FunctionOnContactZero,Solution,ReducedSolution,Vector,ContinuousFunction,Scalar, ReducedFunction,ReducedFunctionOnBoundary,ReducedFunctionOnContactZero,Data, Tensor4, Tensor, getEscriptParamInt, canInterpolate
34  from esys.escript.linearPDEs import LinearPDE,IllegalCoefficientValue,Poisson, IllegalCoefficientFunctionSpace, TransportPDE, IllegalCoefficient, Helmholtz, LameEquation, SolverOptions  from esys.escript.linearPDEs import SolverBuddy, LinearPDE,IllegalCoefficientValue,Poisson, IllegalCoefficientFunctionSpace, TransportPDE, IllegalCoefficient, Helmholtz, LameEquation, SolverOptions
35  import numpy  import numpy
36  import unittest  import esys.escriptcore.utestselect as unittest
37    
38  class Test_linearPDEs(unittest.TestCase):  class Test_linearPDEs(unittest.TestCase):
39      TOL=1.e-6      TOL=1.e-6
40      SOLVER_TOL=1.e-10      SOLVER_TOL=1.e-10
41      DEBUG=False      DEBUG=False
42      VERBOSE=False      VERBOSE=False
43        _domainCanInterpolateAdvanced=None
44    
45        # Can the domain interpolate from ReducedFunction to Function?
46        def specialInterpolationSupported(self):
47            if self._domainCanInterpolateAdvanced is None:
48                self._domainCanInterpolateAdvanced=canInterpolate(ReducedFunction(self.domain), Function(self.domain))
49            return self._domainCanInterpolateAdvanced
50    
51      def check(self,arg,ref_arg,tol=None):      def check(self,arg,ref_arg,tol=None):
52          """          """
53          checks if arg and ref_arg are nearly identical using the L{Lsup<esys.escript.util.Lsup>}          checks if arg and ref_arg are nearly identical using the `Lsup`
54          """          """
55          if tol==None: tol=self.TOL          if tol==None: tol=self.TOL
56          return Lsup(arg-ref_arg)<=tol*Lsup(ref_arg)          return Lsup(arg-ref_arg)<=tol*Lsup(ref_arg)
57    
58        def checkIfNotEmpty(self, pde, coefs):
59            for coef in coefs:
60                self.assertTrue(pde.getCoefficient(coef).isEmpty(),
61                        "%s is not empty"%coef)
62    
63        def checkIfReducedNotEmpty(self, pde, coefs):
64            for coef in coefs:
65                coef += "_reduced"
66                self.assertTrue(pde.getCoefficient(coef).isEmpty(),
67                        "%s is not empty"%coef)
68        
69        def checkContactsNotEmpty(self, pde):
70            if self.domain.supportsContactElements():
71                self.checkIfNotEmpty(pde, ["d_contact", "y_contact"])
72                self.checkIfReducedNotEmpty(pde, ["d_contact", "y_contact"])
73            
74  class Test_LameEquation(Test_linearPDEs):  class Test_LameEquation(Test_linearPDEs):
   
75      def test_config(self):      def test_config(self):
76          mypde=LameEquation(self.domain,debug=self.DEBUG)          mypde=LameEquation(self.domain,debug=self.DEBUG)
77          d=self.domain.getDim()          d=self.domain.getDim()
78          self.failUnlessEqual((mypde.getNumEquations(), mypde.getNumSolutions(), mypde.getSolverOptions().isSymmetric()),(d,d,True),"set up incorrect")          self.assertEqual((mypde.getNumEquations(), mypde.getNumSolutions(), mypde.getSolverOptions().isSymmetric()),(d,d,True),"set up incorrect")
79    
80      def test_setCoefficient_q(self):      def test_setCoefficient_q(self):
81          mypde=LameEquation(self.domain,debug=self.DEBUG)          mypde=LameEquation(self.domain,debug=self.DEBUG)
82          x=self.domain.getX()          x=self.domain.getX()
83          mypde.setValue(q=x)          mypde.setValue(q=x)
   
84          q_ref=interpolate(x,Solution(self.domain))          q_ref=interpolate(x,Solution(self.domain))
85          self.failUnless(self.check(mypde.getCoefficient("A"),0),"A is not 0")  
86          self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")          self.checkIfNotEmpty(mypde, ["B", "C", "D", "X", "Y", "y", "d", "r"])
87          self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")          self.checkIfReducedNotEmpty(mypde, ["A", "B", "C", "D", "X", "Y", "y", "d"])
88          self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")          self.checkContactsNotEmpty(mypde)
89          self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")  
90          self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")          self.assertTrue(self.check(mypde.getCoefficient("A"),0),"A is not 0")
91          self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")          self.assertTrue(self.check(mypde.getCoefficient("q"),q_ref),"q is not empty")
92          self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")  
         self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")  
         self.failUnless(self.check(mypde.getCoefficient("q"),q_ref),"q is not empty")  
         self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")  
93    
94      def test_setCoefficient_r(self):      def test_setCoefficient_r(self):
95          mypde=LameEquation(self.domain,debug=self.DEBUG)          mypde=LameEquation(self.domain,debug=self.DEBUG)
96          x=self.domain.getX()          x=self.domain.getX()
97          mypde.setValue(r=x)          mypde.setValue(r=x)
   
98          r_ref=interpolate(x,Solution(self.domain))          r_ref=interpolate(x,Solution(self.domain))
99          self.failUnless(self.check(mypde.getCoefficient("A"),0),"A is not 0")  
100          self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")          self.checkIfNotEmpty(mypde, ["B", "C", "D", "X", "Y", "y", "d", "q"])
101          self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")          self.checkIfReducedNotEmpty(mypde, ["A", "B", "C", "D", "X", "Y", "y", "d"])
102          self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")          self.checkContactsNotEmpty(mypde)
103          self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")          
104          self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")          self.assertTrue(self.check(mypde.getCoefficient("A"),0),"A is not 0")
105          self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")          self.assertTrue(self.check(mypde.getCoefficient("r"),r_ref),"r is not x")
         self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")  
         self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")  
         self.failUnless(self.check(mypde.getCoefficient("r"),r_ref),"r is nor x")  
         self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")  
106    
107    
108      def test_setCoefficient_F(self):      def test_setCoefficient_F(self):
# Line 134  class Test_LameEquation(Test_linearPDEs) Line 111  class Test_LameEquation(Test_linearPDEs)
111          mypde.setValue(F=x)          mypde.setValue(F=x)
112    
113          Y_ref=interpolate(x,Function(self.domain))          Y_ref=interpolate(x,Function(self.domain))
114          self.failUnless(self.check(mypde.getCoefficient("A"),0),"A is not 0")  
115          self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")          self.checkIfNotEmpty(mypde, ["B", "C", "D", "X", "y", "d", "q", "r"])
116          self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")          self.checkIfReducedNotEmpty(mypde, ["A", "B", "C", "D", "X", "Y", "y", "d"])
117          self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")          self.checkContactsNotEmpty(mypde)
118          self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")  
119          self.failUnless(self.check(mypde.getCoefficient("Y"),Y_ref),"Y is not x")          self.assertTrue(self.check(mypde.getCoefficient("A"),0),"A is not 0")
120          self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")          self.assertTrue(self.check(mypde.getCoefficient("Y"),Y_ref),"Y is not x")
         self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")  
         self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")  
         self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")  
121    
122      def test_setCoefficient_f(self):      def test_setCoefficient_f(self):
123          mypde=LameEquation(self.domain,debug=self.DEBUG)          mypde=LameEquation(self.domain,debug=self.DEBUG)
# Line 163  class Test_LameEquation(Test_linearPDEs) Line 125  class Test_LameEquation(Test_linearPDEs)
125          mypde.setValue(f=x)          mypde.setValue(f=x)
126    
127          y_ref=interpolate(x,FunctionOnBoundary(self.domain))          y_ref=interpolate(x,FunctionOnBoundary(self.domain))
128          self.failUnless(self.check(mypde.getCoefficient("A"),0),"A is not 0")          
129          self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")          self.checkIfNotEmpty(mypde, ["B", "C", "D", "X", "Y", "d", "q", "r"])
130          self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")          self.checkIfReducedNotEmpty(mypde, ["A", "B", "C", "D", "X", "Y", "y", "d"])
131          self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")          self.checkContactsNotEmpty(mypde)
132          self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")  
133          self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")          self.assertTrue(self.check(mypde.getCoefficient("A"),0),"A is not 0")
134          self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")          self.assertTrue(self.check(mypde.getCoefficient("y"),y_ref),"y is not x[0]")
         self.failUnless(self.check(mypde.getCoefficient("y"),y_ref),"d is not x[0]")  
         self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")  
         self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"X_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")  
         self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")  
135    
136      def test_setCoefficient_sigma(self):      def test_setCoefficient_sigma(self):
137          mypde=LameEquation(self.domain,debug=self.DEBUG)          mypde=LameEquation(self.domain,debug=self.DEBUG)
138          x=self.domain.getX()          x=self.domain.getX()
139          mypde.setValue(sigma=outer(x,x))          mypde.setValue(sigma=outer(x,x))
140    
141            self.checkIfNotEmpty(mypde, ["B", "C", "D", "Y", "y", "d", "q", "r"])
142            self.checkIfReducedNotEmpty(mypde, ["A", "B", "C", "D", "X", "Y", "y", "d"])
143            self.checkContactsNotEmpty(mypde)
144    
145          X_ref=interpolate(outer(x,x),Function(self.domain))          X_ref=interpolate(outer(x,x),Function(self.domain))
146          self.failUnless(self.check(mypde.getCoefficient("A"),0),"A is not 0")          self.assertTrue(self.check(mypde.getCoefficient("A"),0),"A is not 0")
147          self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")          self.assertTrue(self.check(mypde.getCoefficient("X"),X_ref),"X is not x X x")
         self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")  
         self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")  
         self.failUnless(self.check(mypde.getCoefficient("X"),X_ref),"X is not x X x")  
         self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")  
         self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")  
         self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")  
         self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"X_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")  
         self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")  
148    
149      def test_setCoefficient_lambda(self):      def test_setCoefficient_lambda(self):
150          mypde=LameEquation(self.domain,debug=self.DEBUG)          mypde=LameEquation(self.domain,debug=self.DEBUG)
# Line 225  class Test_LameEquation(Test_linearPDEs) Line 156  class Test_LameEquation(Test_linearPDEs)
156          k3Xk3=outer(k3,k3)          k3Xk3=outer(k3,k3)
157          A_ref=x[0]*k3Xk3          A_ref=x[0]*k3Xk3
158    
159          self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")          self.checkIfNotEmpty(mypde, ["B", "C", "D", "X", "Y", "y", "d", "q", "r"])
160          self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")          self.checkIfReducedNotEmpty(mypde, ["A", "B", "C", "D", "X", "Y", "y", "d"])
161          self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")          self.checkContactsNotEmpty(mypde)
162          self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")  
163          self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")          self.assertTrue(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
         self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")  
         self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")  
         self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")  
         self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")  
         self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")  
164    
165      def test_setCoefficient_mu(self):      def test_setCoefficient_mu(self):
166          mypde=LameEquation(self.domain,debug=self.DEBUG)          mypde=LameEquation(self.domain,debug=self.DEBUG)
# Line 258  class Test_LameEquation(Test_linearPDEs) Line 172  class Test_LameEquation(Test_linearPDEs)
172          k3Xk3=outer(k3,k3)          k3Xk3=outer(k3,k3)
173          A_ref=x[0]*(swap_axes(k3Xk3,0,3)+swap_axes(k3Xk3,1,3))          A_ref=x[0]*(swap_axes(k3Xk3,0,3)+swap_axes(k3Xk3,1,3))
174    
175          self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")          self.checkIfNotEmpty(mypde, ["B", "C", "D", "X", "Y", "y", "d", "q", "r"])
176          self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")          self.checkIfReducedNotEmpty(mypde, ["A", "B", "C", "D", "X", "Y", "y", "d"])
177          self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")          self.checkContactsNotEmpty(mypde)
178          self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")  
179          self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")          self.assertTrue(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
         self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")  
         self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")  
         self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")  
         self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")  
         self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")  
180    
181      def test_setCoefficient_lambdamu(self):      def test_setCoefficient_lambdamu(self):
182          mypde=LameEquation(self.domain,debug=self.DEBUG)          mypde=LameEquation(self.domain,debug=self.DEBUG)
# Line 290  class Test_LameEquation(Test_linearPDEs) Line 187  class Test_LameEquation(Test_linearPDEs)
187          k3Xk3=outer(k3,k3)          k3Xk3=outer(k3,k3)
188          A_ref=x[0]*k3Xk3+x[1]*(swap_axes(k3Xk3,0,3)+swap_axes(k3Xk3,1,3))          A_ref=x[0]*k3Xk3+x[1]*(swap_axes(k3Xk3,0,3)+swap_axes(k3Xk3,1,3))
189    
190          self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")          self.checkIfNotEmpty(mypde, ["B", "C", "D", "X", "Y", "y", "d", "q", "r"])
191          self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")          self.checkIfReducedNotEmpty(mypde, ["A", "B", "C", "D", "X", "Y", "y", "d"])
192          self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")          self.checkContactsNotEmpty(mypde)
193          self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")  
194          self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")          self.assertTrue(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
         self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")  
         self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")  
         self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")  
         self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")  
         self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")  
195    
196      def test_solve(self):      def test_solve(self):
197         d=self.domain.getDim()         d=self.domain.getDim()
# Line 324  class Test_LameEquation(Test_linearPDEs) Line 204  class Test_LameEquation(Test_linearPDEs)
204         mypde.setValue(q=msk,r=u_ex,lame_mu=3,lame_lambda=50,f=(2*3+50*d)*FunctionOnBoundary(self.domain).getNormal())         mypde.setValue(q=msk,r=u_ex,lame_mu=3,lame_lambda=50,f=(2*3+50*d)*FunctionOnBoundary(self.domain).getNormal())
205    
206         u=mypde.getSolution()         u=mypde.getSolution()
207         self.failUnless(self.check(u,u_ex,10*self.TOL),"incorrect solution")         self.assertTrue(self.check(u,u_ex,10*self.TOL),"incorrect solution")
208    
209  class Test_Helmholtz(Test_linearPDEs):  class Test_Helmholtz(Test_linearPDEs):
210    
211      def test_config(self):      def test_config(self):
212          mypde=Helmholtz(self.domain,debug=self.DEBUG)          mypde=Helmholtz(self.domain,debug=self.DEBUG)
213          self.failUnlessEqual((mypde.getNumEquations(), mypde.getNumSolutions(), mypde.getSolverOptions().isSymmetric()),(1,1,True),"set up incorrect")          self.assertEqual((mypde.getNumEquations(), mypde.getNumSolutions(), mypde.getSolverOptions().isSymmetric()),(1,1,True),"set up incorrect")
214      def test_setCoefficient_q(self):      def test_setCoefficient_q(self):
215          mypde=Helmholtz(self.domain,debug=self.DEBUG)          mypde=Helmholtz(self.domain,debug=self.DEBUG)
216          x=self.domain.getX()          x=self.domain.getX()
# Line 339  class Test_Helmholtz(Test_linearPDEs): Line 219  class Test_Helmholtz(Test_linearPDEs):
219          q_ref=interpolate(whereZero(x[0]),Solution(self.domain))          q_ref=interpolate(whereZero(x[0]),Solution(self.domain))
220          A_ref=kronecker(self.domain)          A_ref=kronecker(self.domain)
221    
222          self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")          self.checkIfNotEmpty(mypde, ["B", "C", "D", "X", "Y", "y", "d", "r"])
223          self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")          self.checkIfReducedNotEmpty(mypde, ["A", "B", "C", "D", "X", "Y", "y", "d"])
224          self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")          self.checkContactsNotEmpty(mypde)
225          self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")  
226          self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")          self.assertTrue(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
227          self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")          self.assertTrue(self.check(mypde.getCoefficient("q"),q_ref),"q is not empty")
         self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")  
         self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")  
         self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")  
         self.failUnless(self.check(mypde.getCoefficient("q"),q_ref),"q is not empty")  
         self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")  
228    
229      def test_setCoefficient_r(self):      def test_setCoefficient_r(self):
230          mypde=Helmholtz(self.domain,debug=self.DEBUG)          mypde=Helmholtz(self.domain,debug=self.DEBUG)
# Line 369  class Test_Helmholtz(Test_linearPDEs): Line 233  class Test_Helmholtz(Test_linearPDEs):
233    
234          r_ref=interpolate(x[0],Solution(self.domain))          r_ref=interpolate(x[0],Solution(self.domain))
235          A_ref=kronecker(self.domain)          A_ref=kronecker(self.domain)
236          self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")  
237          self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")          self.checkIfNotEmpty(mypde, ["B", "C", "D", "X", "Y", "y", "d", "q"])
238          self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")          self.checkIfReducedNotEmpty(mypde, ["A", "B", "C", "D", "X", "Y", "y", "d"])
239          self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")          self.checkContactsNotEmpty(mypde)
240          self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")  
241          self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")          self.assertTrue(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
242          self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")          self.assertTrue(self.check(mypde.getCoefficient("r"),r_ref),"r is nor x[0]")
         self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")  
         self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")  
         self.failUnless(self.check(mypde.getCoefficient("r"),r_ref),"r is nor x[0]")  
         self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")  
243    
244    
245      def test_setCoefficient_f(self):      def test_setCoefficient_f(self):
# Line 400  class Test_Helmholtz(Test_linearPDEs): Line 249  class Test_Helmholtz(Test_linearPDEs):
249    
250          Y_ref=interpolate(x[0],Function(self.domain))          Y_ref=interpolate(x[0],Function(self.domain))
251          A_ref=kronecker(self.domain)          A_ref=kronecker(self.domain)
252          self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")  
253          self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")          self.checkIfNotEmpty(mypde, ["B", "C", "D", "X", "y", "d", "q", "r"])
254          self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")          self.checkIfReducedNotEmpty(mypde, ["A", "B", "C", "D", "X", "Y", "y", "d"])
255          self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")          self.checkContactsNotEmpty(mypde)
256          self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")  
257          self.failUnless(self.check(mypde.getCoefficient("Y"),Y_ref),"Y is not x[0]")          self.assertTrue(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
258          self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")          self.assertTrue(self.check(mypde.getCoefficient("Y"),Y_ref),"Y is not x[0]")
         self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")  
         self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")  
         self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")  
259    
260      def test_setCoefficient_alpha(self):      def test_setCoefficient_alpha(self):
261          mypde=Helmholtz(self.domain,debug=self.DEBUG)          mypde=Helmholtz(self.domain,debug=self.DEBUG)
# Line 430  class Test_Helmholtz(Test_linearPDEs): Line 264  class Test_Helmholtz(Test_linearPDEs):
264    
265          d_ref=interpolate(x[0],FunctionOnBoundary(self.domain))          d_ref=interpolate(x[0],FunctionOnBoundary(self.domain))
266          A_ref=kronecker(self.domain)          A_ref=kronecker(self.domain)
267          self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")  
268          self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")          self.checkIfNotEmpty(mypde, ["B", "C", "D", "X", "Y", "y", "q", "r"])
269          self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")          self.checkIfReducedNotEmpty(mypde, ["A", "B", "C", "D", "X", "Y", "y", "d"])
270          self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")          self.checkContactsNotEmpty(mypde)
271          self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")  
272          self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")          self.assertTrue(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
273          self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")          self.assertTrue(self.check(mypde.getCoefficient("d"),d_ref),"d is not x[0]")
         self.failUnless(self.check(mypde.getCoefficient("d"),d_ref),"d is not x[0]")  
         self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")  
         self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"X_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")  
         self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")  
274    
275      def test_setCoefficient_g(self):      def test_setCoefficient_g(self):
276          mypde=Helmholtz(self.domain,debug=self.DEBUG)          mypde=Helmholtz(self.domain,debug=self.DEBUG)
# Line 460  class Test_Helmholtz(Test_linearPDEs): Line 279  class Test_Helmholtz(Test_linearPDEs):
279    
280          y_ref=interpolate(x[0],FunctionOnBoundary(self.domain))          y_ref=interpolate(x[0],FunctionOnBoundary(self.domain))
281          A_ref=kronecker(self.domain)          A_ref=kronecker(self.domain)
282          self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")  
283          self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")          self.checkIfNotEmpty(mypde, ["B", "C", "D", "X", "Y", "d", "q", "r"])
284          self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")          self.checkIfReducedNotEmpty(mypde, ["A", "B", "C", "D", "X", "Y", "y", "d"])
285          self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")          self.checkContactsNotEmpty(mypde)
286          self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")  
287          self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")          self.assertTrue(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
288          self.failUnless(self.check(mypde.getCoefficient("y"),y_ref),"y is not x[0]")          self.assertTrue(self.check(mypde.getCoefficient("y"),y_ref),"y is not x[0]")
         self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")  
         self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")  
         self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")  
289    
290      def test_setCoefficient_omega(self):      def test_setCoefficient_omega(self):
291          mypde=Helmholtz(self.domain,debug=self.DEBUG)          mypde=Helmholtz(self.domain,debug=self.DEBUG)
# Line 490  class Test_Helmholtz(Test_linearPDEs): Line 294  class Test_Helmholtz(Test_linearPDEs):
294    
295          D_ref=interpolate(x[0],Function(self.domain))          D_ref=interpolate(x[0],Function(self.domain))
296          A_ref=kronecker(self.domain)          A_ref=kronecker(self.domain)
297          self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")  
298          self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")          self.checkIfNotEmpty(mypde, ["B", "C", "X", "Y", "y", "d", "q", "r"])
299          self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")          self.checkIfReducedNotEmpty(mypde, ["A", "B", "C", "D", "X", "Y", "y", "d"])
300          self.failUnless(self.check(mypde.getCoefficient("D"),D_ref),"D is not x[0]")          self.checkContactsNotEmpty(mypde)
301          self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")  
302          self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")          self.assertTrue(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
303          self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")          self.assertTrue(self.check(mypde.getCoefficient("D"),D_ref),"D is not x[0]")
         self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")  
         self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")  
         self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")  
304    
305      def test_solve(self):      def test_solve(self):
306         d=self.domain.getDim()         d=self.domain.getDim()
# Line 520  class Test_Helmholtz(Test_linearPDEs): Line 309  class Test_Helmholtz(Test_linearPDEs):
309         mypde=Helmholtz(self.domain)         mypde=Helmholtz(self.domain)
310         mypde.setValue(f=3,omega=3,alpha=2,g=2)         mypde.setValue(f=3,omega=3,alpha=2,g=2)
311         u=mypde.getSolution()         u=mypde.getSolution()
312         self.failUnless(self.check(u,u_ex,10*self.TOL),"incorrect solution")         self.assertTrue(self.check(u,u_ex,10*self.TOL),"incorrect solution")
313    
314  class Test_Poisson(Test_linearPDEs):  class Test_Poisson(Test_linearPDEs):
315    
316      def test_config(self):      def test_config(self):
317          mypde=Poisson(self.domain,debug=self.DEBUG)          mypde=Poisson(self.domain,debug=self.DEBUG)
318          self.failUnlessEqual((mypde.getNumEquations(), mypde.getNumSolutions(), mypde.getSolverOptions().isSymmetric()),(1,1,True),"set up incorrect")          self.assertEqual((mypde.getNumEquations(), mypde.getNumSolutions(), mypde.getSolverOptions().isSymmetric()),(1,1,True),"set up incorrect")
319    
320      def test_setCoefficient_q(self):      def test_setCoefficient_q(self):
321          mypde=Poisson(self.domain,debug=self.DEBUG)          mypde=Poisson(self.domain,debug=self.DEBUG)
322          x=self.domain.getX()          x=self.domain.getX()
323          q_ref=interpolate(whereZero(x[0]),Solution(self.domain))          q_ref=interpolate(whereZero(x[0]),Solution(self.domain))
324          A_ref=kronecker(self.domain)          A_ref=kronecker(self.domain)
325          mypde.setValue(q=whereZero(x[0]))          mypde.setValue(q=whereZero(x[0]))
326          self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")  
327          self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")          self.checkIfNotEmpty(mypde, ["B", "C", "D", "X", "Y", "y", "d", "r"])
328          self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")          self.checkIfReducedNotEmpty(mypde, ["A", "B", "C", "D", "X", "Y", "y", "d"])
329          self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")          self.checkContactsNotEmpty(mypde)
330          self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")  
331          self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")          self.assertTrue(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
332          self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")          self.assertTrue(self.check(mypde.getCoefficient("q"),q_ref),"q is not empty")
333          self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")  
         self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")  
         self.failUnless(self.check(mypde.getCoefficient("q"),q_ref),"q is not empty")  
         self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")  
334      def test_setCoefficient_f(self):      def test_setCoefficient_f(self):
335          mypde=Poisson(self.domain,debug=self.DEBUG)          mypde=Poisson(self.domain,debug=self.DEBUG)
336          x=self.domain.getX()          x=self.domain.getX()
337          Y_ref=interpolate(x[0],Function(self.domain))          Y_ref=interpolate(x[0],Function(self.domain))
338          A_ref=kronecker(self.domain)          A_ref=kronecker(self.domain)
339          mypde.setValue(f=x[0])          mypde.setValue(f=x[0])
340          self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")  
341          self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")          self.checkIfNotEmpty(mypde, ["B", "C", "D", "X", "y", "d", "q", "r"])
342          self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")          self.checkIfReducedNotEmpty(mypde, ["A", "B", "C", "D", "X", "Y", "y", "d"])
343          self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")          self.checkContactsNotEmpty(mypde)
344          self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")  
345          self.failUnless(self.check(mypde.getCoefficient("Y"),Y_ref),"Y is not x[0]")          self.assertTrue(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
346          self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")          self.assertTrue(self.check(mypde.getCoefficient("Y"),Y_ref),"Y is not x[0]")
347          self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")  
         self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")  
         self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")  
348      def test_setCoefficient_f_reduced(self):      def test_setCoefficient_f_reduced(self):
349          mypde=Poisson(self.domain,debug=self.DEBUG)          mypde=Poisson(self.domain,debug=self.DEBUG)
350          x=self.domain.getX()          x=self.domain.getX()
351          Y_ref=interpolate(x[0],ReducedFunction(self.domain))          Y_ref=interpolate(x[0],ReducedFunction(self.domain))
352          A_ref=kronecker(self.domain)          A_ref=kronecker(self.domain)
353          mypde.setValue(f_reduced=x[0])          mypde.setValue(f_reduced=x[0])
354          self.failUnless(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")  
355          self.failUnless(mypde.getCoefficient("B").isEmpty(),"B is not empty")          self.checkIfNotEmpty(mypde, ["B", "C", "D", "X", "Y", "y", "d", "q", "r"])
356          self.failUnless(mypde.getCoefficient("C").isEmpty(),"C is not empty")          self.checkIfReducedNotEmpty(mypde, ["A", "B", "C", "D", "X", "y", "d"])
357          self.failUnless(mypde.getCoefficient("D").isEmpty(),"D is not empty")          self.checkContactsNotEmpty(mypde)
358          self.failUnless(mypde.getCoefficient("X").isEmpty(),"X is not empty")  
359          self.failUnless(mypde.getCoefficient("Y").isEmpty(),"Y is not empty")          self.assertTrue(self.check(mypde.getCoefficient("A"),A_ref),"A is not kronecker")
360          self.failUnless(mypde.getCoefficient("y").isEmpty(),"y is not empty")          self.assertTrue(self.check(mypde.getCoefficient("Y_reduced"),Y_ref),"Y_reduced is not x[0]")
361          self.failUnless(mypde.getCoefficient("d").isEmpty(),"d is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty")  
         self.failUnless(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty")  
         self.failUnless(self.check(mypde.getCoefficient("Y_reduced"),Y_ref),"Y_reduced is not x[0]")  
         self.failUnless(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty")  
         self.failUnless(mypde.getCoefficient("q").isEmpty(),"q is not empty")  
         self.failUnless(mypde.getCoefficient("r").isEmpty(),"r is not empty")  
362      def test_solve(self):      def test_solve(self):
363         d=self.domain.getDim()         d=self.domain.getDim()
364         cf=ContinuousFunction(self.domain)         cf=ContinuousFunction(self.domain)
# Line 636  class Test_Poisson(Test_linearPDEs): Line 384  class Test_Poisson(Test_linearPDEs):
384         mypde=Poisson(self.domain)         mypde=Poisson(self.domain)
385         mypde.setValue(f=f,q=msk)         mypde.setValue(f=f,q=msk)
386         u=mypde.getSolution()         u=mypde.getSolution()
387         self.failUnless(self.check(u,u_ex,10*self.TOL),"incorrect solution")         self.assertTrue(self.check(u,u_ex,10*self.TOL),"incorrect solution")
388    
389  class Test_LinearPDE_noLumping(Test_linearPDEs):  class Test_LinearPDE_noLumping(Test_linearPDEs):
390      N=4      N=4
     def test_SolverOptions(self):  
         so=SolverOptions()  
391    
392          self.failUnless(so.getLevelMax() == 10, "initial  LevelMax is wrong.")      def test_SolverOptions(self):
393          self.failUnlessRaises(ValueError,so.setLevelMax,-1)          sb=SolverBuddy()
394          so.setLevelMax(3)          so=SolverOptions
395          self.failUnless(so.getLevelMax() == 3, "LevelMax is wrong.")          self.assertTrue(sb.getSmoother() == so.GAUSS_SEIDEL, "initial Smoother is wrong.")
396            self.assertRaises(ValueError,sb.setSmoother,-1)
397          self.failUnless(so.getCoarseningThreshold() == 0.05, "initial  CoarseningThreshold is wrong.")          sb.setSmoother(so.GAUSS_SEIDEL)
398          self.failUnlessRaises(ValueError,so.setCoarseningThreshold,-1)          self.assertTrue(sb.getSmoother() == so.GAUSS_SEIDEL, "Gauss-Seidel smoother is not set.")
399          so.setCoarseningThreshold(0.1)          sb.setSmoother(so.JACOBI)
400          self.failUnless(so.getCoarseningThreshold() == 0.1, "CoarseningThreshold is wrong.")          self.assertTrue(sb.getSmoother() == so.JACOBI, "Jacobi smoother is not set.")
401    
402            self.assertTrue(sb.getLevelMax() == 100, "initial LevelMax is wrong.")
403            self.assertRaises(ValueError,sb.setLevelMax,-1)
404            sb.setLevelMax(20)
405            self.assertTrue(sb.getLevelMax() == 20, "LevelMax is wrong.")
406    
407            self.assertTrue(sb.getCoarseningThreshold() == 0.25, "initial CoarseningThreshold is wrong.")
408            self.assertRaises(ValueError,sb.setCoarseningThreshold,-1)
409            sb.setCoarseningThreshold(0.1)
410            self.assertTrue(sb.getCoarseningThreshold() == 0.1, "CoarseningThreshold is wrong.")
411            
412            self.assertTrue(sb.getMinCoarseMatrixSize() == 500, "initial Minimum Coarse Matrix Size is wrong.")
413            self.assertRaises(ValueError,sb.setMinCoarseMatrixSize,-1)
414            sb.setMinCoarseMatrixSize(1000)
415            self.assertTrue(sb.getMinCoarseMatrixSize() == 1000, "Minimum Coarse Matrix Size is wrong.")
416    
417            self.assertTrue(sb.getNumSweeps() == 1, "initial number of sweeps is wrong.")
418            self.assertRaises(ValueError,sb.setNumSweeps,-1)
419            sb.setNumSweeps(3)
420            self.assertTrue(sb.getNumSweeps() == 3, "Sweeps is wrong.")
421    
422            self.assertTrue(sb.getNumPreSweeps() == 1, "initial PreSweeps is wrong.")
423            self.assertRaises(ValueError,sb.setNumPreSweeps,-1)
424            sb.setNumPreSweeps(4)
425            self.assertTrue(sb.getNumPreSweeps() == 4, "PreSweeps is wrong.")
426    
427            self.assertTrue(sb.getNumPostSweeps() == 1, "initial PostSweeps is wrong.")
428            self.assertRaises(ValueError,sb.setNumPostSweeps,-1)
429            sb.setNumPostSweeps(5)
430            self.assertTrue(sb.getNumPostSweeps() == 5, "PostSweeps is wrong.")
431    
432            self.assertTrue(sb.getTolerance() == 1.e-8, "initial Tolerance is wrong.")
433            self.assertRaises(ValueError,sb.setTolerance,-1)
434            sb.setTolerance(0.2)
435            self.assertTrue(sb.getTolerance() == 0.2, "Tolerance is wrong.")
436    
437            self.assertTrue(sb.getAbsoluteTolerance() == 0., "initial AbsoluteTolerance is wrong.")
438            self.assertRaises(ValueError,sb.setAbsoluteTolerance,-1)
439            sb.setAbsoluteTolerance(0.3)
440            self.assertTrue(sb.getAbsoluteTolerance() == 0.3, "AbsoluteTolerance is wrong.")
441    
442            self.assertTrue(sb.getInnerTolerance() == 0.9, "initial InnerTolerance is wrong.")
443            self.assertRaises(ValueError,sb.setInnerTolerance,-1)
444            sb.setInnerTolerance(0.4)
445            self.assertTrue(sb.getInnerTolerance() == 0.4, "InnerTolerance is wrong.")
446    
447            self.assertTrue(sb.getDropTolerance() == 0.01, "initial DropTolerance is wrong.")
448            self.assertRaises(ValueError,sb.setDropTolerance,-1)
449            sb.setDropTolerance(0.5)
450            self.assertTrue(sb.getDropTolerance() == 0.5, "DropDropTolerance is wrong.")
451    
452            self.assertTrue(sb.getDropStorage() == 2., "initial DropStorage is wrong.")
453            self.assertRaises(ValueError,sb.setDropStorage,-1)
454            sb.setDropStorage(10)
455            self.assertTrue(sb.getDropStorage() == 10, "DropStorage is wrong.")
456                    
457          self.failUnless(so.getMinCoarseMatrixSize() == 500, "initial  Minimum Coarse Matrix Size is wrong.")          self.assertTrue(sb.getRelaxationFactor() == 0.3, "initial RelaxationFactor is wrong.")
458          self.failUnlessRaises(ValueError,so.setMinCoarseMatrixSize,-1)          self.assertRaises(ValueError,sb.setRelaxationFactor,-1)
459          so.setMinCoarseMatrixSize(1000)          sb.setRelaxationFactor(0.1)
460          self.failUnless(so.getMinCoarseMatrixSize() == 1000, "Minimum Coarse Matrix Size is wrong.")          self.assertTrue(sb.getRelaxationFactor() == 0.1, "Relaxation is wrong.")
461    
462          self.failUnless(so.getNumSweeps() == 2, "initial  Sweeps is wrong.")          self.assertTrue(sb.getIterMax() == 100000, "initial IterMax is wrong.")
463          self.failUnlessRaises(ValueError,so.setNumSweeps,-1)          self.assertRaises(ValueError,sb.setIterMax,0)
464          so.setNumSweeps(3)          sb.setIterMax(11)
465          self.failUnless(so.getNumSweeps() == 3, "Sweeps is wrong.")          self.assertTrue(sb.getIterMax() == 11, "IterMax is wrong.")
466    
467          self.failUnless(so.getNumPreSweeps() == 2, "initial  PreSweeps is wrong.")          self.assertTrue(sb.getInnerIterMax() == 10, "initial InnerIterMax is wrong.")
468          self.failUnlessRaises(ValueError,so.setNumPreSweeps,-1)          self.assertRaises(ValueError,sb.setInnerIterMax,0)
469          so.setNumPreSweeps(4)          sb.setInnerIterMax(12)
470          self.failUnless(so.getNumPreSweeps() == 4, "PreSweeps is wrong.")          self.assertTrue(sb.getInnerIterMax() == 12, "InnerIterMax is wrong.")
471    
472          self.failUnless(so.getNumPostSweeps() == 2, "initial  PreSweeps is wrong.")          self.assertTrue(sb.getTruncation() == 20, "initial Truncation is wrong.")
473          self.failUnlessRaises(ValueError,so.setNumPostSweeps,-1)          self.assertRaises(ValueError,sb.setTruncation,0)
474          so.setNumPostSweeps(5)          sb.setTruncation(13)
475          self.failUnless(so.getNumPostSweeps() == 5, "PreSweeps is wrong.")          self.assertTrue(sb.getTruncation() == 13, "Truncation is wrong.")
476    
477          self.failUnless(so.getTolerance() == 1.e-8, "initial Tolerance is wrong.")          self.assertTrue(sb.getRestart() == 0, "initial Truncation is wrong.")
478          self.failUnlessRaises(ValueError,so.setTolerance,-1)          self.assertRaises(ValueError,sb.setTruncation,0)
479          so.setTolerance(0.2)          sb.setRestart(14)
480          self.failUnless(so.getTolerance() == 0.2, "Tolerance is wrong.")          self.assertTrue(sb.getRestart() == 14, "Truncation is wrong.")
481            sb.setRestart(0)
482          self.failUnless(so.getAbsoluteTolerance() == 0., "initial  AbsoluteTolerance is wrong.")          self.assertTrue(sb.getRestart() == 0, "Truncation is wrong.")
         self.failUnlessRaises(ValueError,so.setAbsoluteTolerance,-1)  
         so.setAbsoluteTolerance(0.3)  
         self.failUnless(so.getAbsoluteTolerance() == 0.3, "AbsoluteTolerance is wrong.")  
   
         self.failUnless(so.getInnerTolerance() == 0.9, "initial InnerTolerance is wrong.")  
         self.failUnlessRaises(ValueError,so.setInnerTolerance,-1)  
         so.setInnerTolerance(0.4)  
         self.failUnless(so.getInnerTolerance() == 0.4, "InnerTolerance is wrong.")  
   
         self.failUnless(so.getDropTolerance() == 0.01, "initial DropTolerance is wrong.")  
         self.failUnlessRaises(ValueError,so.setDropTolerance,-1)  
         so.setDropTolerance(0.5)  
         self.failUnless(so.getDropTolerance() == 0.5, "DropDropTolerance is wrong.")  
   
         self.failUnless(so.getDropStorage() == 2., "initial DropStorage is wrong.")  
         self.failUnlessRaises(ValueError,so.setDropStorage,-1)  
         so.setDropStorage(10)  
         self.failUnless(so.getDropStorage() == 10, "DropStorage is wrong.")  
483                    
484          self.failUnless(so.getRelaxationFactor() == 0.3, "initial RelaxationFactor is wrong.")          self.assertTrue(not sb.isVerbose(), "initial verbosity flag is wrong.")
485          self.failUnlessRaises(ValueError,so.setRelaxationFactor,-1)          sb.setVerbosityOn()
486          so.setRelaxationFactor(0.1)          self.assertTrue(sb.isVerbose(), "verbosity (1) flag is wrong.")
487          self.failUnless(so.getRelaxationFactor() == 0.1, "Relaxation is wrong.")          sb.setVerbosityOff()
488            self.assertTrue(not sb.isVerbose(), "verbosity (2) flag is wrong.")
489            sb.setVerbosity(True)
490          self.failUnless(so.getIterMax() == 100000, "initial IterMax is wrong.")          self.assertTrue(sb.isVerbose(), "verbosity (3) flag is wrong.")
491          self.failUnlessRaises(ValueError,so.setIterMax,0)          sb.setVerbosity(False)
492          so.setIterMax(11)          self.assertTrue(not sb.isVerbose(), "verbosity (4) flag is wrong.")
493          self.failUnless(so.getIterMax() == 11, "IterMax is wrong.")  
494            self.assertTrue(not sb.isSymmetric(), "initial symmetry flag is wrong.")
495          self.failUnless(so.getInnerIterMax() == 10, "initial InnerIterMax is wrong.")          sb.setSymmetryOn()
496          self.failUnlessRaises(ValueError,so.setInnerIterMax,0)          self.assertTrue(sb.isSymmetric(), "symmetry (1) flag is wrong.")
497          so.setInnerIterMax(12)          sb.setSymmetryOff()
498          self.failUnless(so.getInnerIterMax() == 12, "InnerIterMax is wrong.")          self.assertTrue(not sb.isSymmetric(), "symmetry (2) flag is wrong.")
499            sb.setSymmetry(True)
500          self.failUnless(so.getTruncation() == 20, "initial Truncation is wrong.")          self.assertTrue(sb.isSymmetric(), "symmetry (3) flag is wrong.")
501          self.failUnlessRaises(ValueError,so.setTruncation,0)          sb.setSymmetry(False)
502          so.setTruncation(13)          self.assertTrue(not sb.isSymmetric(), "symmetry (4) flag is wrong.")
503          self.failUnless(so.getTruncation() == 13, "Truncation is wrong.")  
504            self.assertTrue(sb.adaptInnerTolerance(), "initial InnerToleranceAdaption flag is wrong.")
505          self.failUnless(so.getRestart() == None, "initial Truncation is wrong.")          sb.setInnerToleranceAdaptionOn()
506          self.failUnlessRaises(ValueError,so.setTruncation,0)          self.assertTrue(sb.adaptInnerTolerance(), "InnerToleranceAdaption (1) flag is wrong.")
507          so.setRestart(14)          sb.setInnerToleranceAdaptionOff()
508          self.failUnless(so.getRestart() == 14, "Truncation is wrong.")          self.assertTrue(not sb.adaptInnerTolerance(), "InnerToleranceAdaption (2) flag is wrong.")
509          so.setRestart(None)          sb.setInnerToleranceAdaption(adapt=True)
510          self.failUnless(so.getRestart() == None, "Truncation is wrong.")          self.assertTrue(sb.adaptInnerTolerance(), "InnerToleranceAdaption (3) flag is wrong.")
511                sb.setInnerToleranceAdaption(adapt=False)
512          self.failUnless(not so.isVerbose(), "initial verbosity flag is wrong.")          self.assertTrue(not sb.adaptInnerTolerance(), "InnerToleranceAdaption (4) flag is wrong.")
         so.setVerbosityOn()  
         self.failUnless(so.isVerbose(), "verbosity (1) flag is wrong.")  
         so.setVerbosityOff()  
         self.failUnless(not so.isVerbose(), "verbosity (2) flag is wrong.")  
         so.setVerbosity(verbose=True)  
         self.failUnless(so.isVerbose(), "verbosity (3) flag is wrong.")  
         so.setVerbosity(verbose=False)  
         self.failUnless(not so.isVerbose(), "verbosity (4) flag is wrong.")  
   
         self.failUnless(not so.isSymmetric(), "initial symmetry flag is wrong.")  
         so.setSymmetryOn()  
         self.failUnless(so.isSymmetric(), "symmetry (1) flag is wrong.")  
         so.setSymmetryOff()  
         self.failUnless(not so.isSymmetric(), "symmetry (2) flag is wrong.")  
         so.setSymmetry(flag=True)  
         self.failUnless(so.isSymmetric(), "symmetry (3) flag is wrong.")  
         so.setSymmetry(flag=False)  
         self.failUnless(not so.isSymmetric(), "symmetry (4) flag is wrong.")  
   
         self.failUnless(so.adaptInnerTolerance(), "initial InnerToleranceAdaption flag is wrong.")  
         so.setInnerToleranceAdaptionOn()  
         self.failUnless(so.adaptInnerTolerance(), "InnerToleranceAdaption (1) flag is wrong.")  
         so.setInnerToleranceAdaptionOff()  
         self.failUnless(not so.adaptInnerTolerance(), "InnerToleranceAdaption (2) flag is wrong.")  
         so.setInnerToleranceAdaption(adapt=True)  
         self.failUnless(so.adaptInnerTolerance(), "InnerToleranceAdaption (3) flag is wrong.")  
         so.setInnerToleranceAdaption(adapt=False)  
         self.failUnless(not so.adaptInnerTolerance(), "InnerToleranceAdaption (4) flag is wrong.")  
513            
514          self.failUnless(not so.acceptConvergenceFailure(), "initial acceptConvergenceFailure flag is wrong.")          self.assertTrue(not sb.acceptConvergenceFailure(), "initial acceptConvergenceFailure flag is wrong.")
515          so.setAcceptanceConvergenceFailureOn()          sb.setAcceptanceConvergenceFailureOn()
516          self.failUnless(so.acceptConvergenceFailure(), "acceptConvergenceFailure (1) flag is wrong.")          self.assertTrue(sb.acceptConvergenceFailure(), "acceptConvergenceFailure (1) flag is wrong.")
517          so.setAcceptanceConvergenceFailureOff()          sb.setAcceptanceConvergenceFailureOff()
518          self.failUnless(not so.acceptConvergenceFailure(), "acceptConvergenceFailure (2) flag is wrong.")          self.assertTrue(not sb.acceptConvergenceFailure(), "acceptConvergenceFailure (2) flag is wrong.")
519          so.setAcceptanceConvergenceFailure(accept=True)          sb.setAcceptanceConvergenceFailure(accept=True)
520          self.failUnless(so.acceptConvergenceFailure(), "acceptConvergenceFailure (3) flag is wrong.")          self.assertTrue(sb.acceptConvergenceFailure(), "acceptConvergenceFailure (3) flag is wrong.")
521          so.setAcceptanceConvergenceFailure(accept=False)          sb.setAcceptanceConvergenceFailure(accept=False)
522          self.failUnless(not so.acceptConvergenceFailure(), "acceptConvergenceFailure (4) flag is wrong.")            self.assertTrue(not sb.acceptConvergenceFailure(), "acceptConvergenceFailure (4) flag is wrong.")  
523                    
524          self.failUnless(so.getReordering() == 30, "initial Reordering is wrong.")          self.assertTrue(sb.getReordering() == so.DEFAULT_REORDERING, "initial Reordering is wrong.")
525          self.failUnlessRaises(ValueError,so.setReordering,-1)          self.assertRaises(ValueError,sb.setReordering,-1)
526          so.setReordering(so.NO_REORDERING)          sb.setReordering(so.NO_REORDERING)
527          self.failUnless(so.getReordering() == 17, "NO_REORDERING is not set.")          self.assertTrue(sb.getReordering() == so.NO_REORDERING, "NO_REORDERING is not set.")
528          so.setReordering(so.MINIMUM_FILL_IN)          sb.setReordering(so.MINIMUM_FILL_IN)
529          self.failUnless(so.getReordering() == 18, "MINIMUM_FILL_IN is not set.")          self.assertTrue(sb.getReordering() == so.MINIMUM_FILL_IN, "MINIMUM_FILL_IN is not set.")
530          so.setReordering(so.NESTED_DISSECTION)          sb.setReordering(so.NESTED_DISSECTION)
531          self.failUnless(so.getReordering() == 19, "NESTED_DISSECTION is not set.")          self.assertTrue(sb.getReordering() == so.NESTED_DISSECTION, "NESTED_DISSECTION is not set.")
532          so.setReordering(so.DEFAULT_REORDERING)          sb.setReordering(so.DEFAULT_REORDERING)
533          self.failUnless(so.getReordering() == 30, "DEFAULT_REORDERING is not set.")          self.assertTrue(sb.getReordering() == so.DEFAULT_REORDERING, "DEFAULT_REORDERING is not set.")
534                    
535          self.failUnless(so.getPackage() == 0, "initial solver package is wrong.")          self.assertTrue(sb.getPackage() == so.DEFAULT, "initial solver package is wrong.")
536          self.failUnlessRaises(ValueError,so.setPackage,-1)          self.assertRaises(ValueError,sb.setPackage,-1)
537          so.setPackage(so.PASO)          sb.setPackage(so.PASO)
538          self.failUnless(so.getPackage() == 21, "PASO is not set.")          self.assertTrue(sb.getPackage() == so.PASO, "PASO is not set.")
539          so.setPackage(so.SUPER_LU)          sb.setPackage(so.CUSP)
540          self.failUnless(so.getPackage() == 31, "SUPER_LU is not set.")          self.assertTrue(sb.getPackage() == so.CUSP, "CUSP is not set.")
541          so.setPackage(so.PASTIX)          sb.setPackage(so.SUPER_LU)
542          self.failUnless(so.getPackage() == 32, "PASTIX is not set.")          self.assertTrue(sb.getPackage() == so.SUPER_LU, "SUPER_LU is not set.")
543          so.setPackage(so.MKL)          sb.setPackage(so.PASTIX)
544          self.failUnless(so.getPackage() == 15, "MKL is not set.")          self.assertTrue(sb.getPackage() == so.PASTIX, "PASTIX is not set.")
545          so.setPackage(so.UMFPACK)          sb.setPackage(so.MKL)
546          self.failUnless(so.getPackage() == 16, "UMFPACK is not set.")          self.assertTrue(sb.getPackage() == so.MKL, "MKL is not set.")
547          so.setPackage(so.TRILINOS)          sb.setPackage(so.UMFPACK)
548          self.failUnless(so.getPackage() == 24, "TRILINOS is not set.")          self.assertTrue(sb.getPackage() == so.UMFPACK, "UMFPACK is not set.")
549            sb.setPackage(so.TRILINOS)
550          self.failUnless(so.getSolverMethod() == 0, "initial SolverMethod is wrong.")          self.assertTrue(sb.getPackage() == so.TRILINOS, "TRILINOS is not set.")
551          self.failUnlessRaises(ValueError,so.setSolverMethod,-1)  
552          so.setSolverMethod(so.DIRECT)          self.assertTrue(sb.getSolverMethod() == so.DEFAULT, "initial SolverMethod is wrong.")
553          self.failUnless(so.getSolverMethod() == 1, "DIRECT is not set.")          self.assertRaises(ValueError,sb.setSolverMethod,-1)
554          so.setSolverMethod(so.CHOLEVSKY)          sb.setSolverMethod(so.DIRECT)
555          self.failUnless(so.getSolverMethod() == 2, "CHOLEVSKY is not set.")          self.assertTrue(sb.getSolverMethod() == so.DIRECT, "DIRECT is not set.")
556          so.setSolverMethod(so.PCG)          sb.setSolverMethod(so.CHOLEVSKY)
557          self.failUnless(so.getSolverMethod() == 3, "PCG is not set.")          self.assertTrue(sb.getSolverMethod() == so.CHOLEVSKY, "CHOLEVSKY is not set.")
558          so.setSolverMethod(so.CR)          sb.setSolverMethod(so.PCG)
559          self.failUnless(so.getSolverMethod() == 4, "CR is not set.")          self.assertTrue(sb.getSolverMethod() == so.PCG, "PCG is not set.")
560          so.setSolverMethod(so.CGS)          sb.setSolverMethod(so.CR)
561          self.failUnless(so.getSolverMethod() == 5, "CGS is not set.")          self.assertTrue(sb.getSolverMethod() == so.CR, "CR is not set.")
562          so.setSolverMethod(so.BICGSTAB)          sb.setSolverMethod(so.CGS)
563          self.failUnless(so.getSolverMethod() == 6, "BICGSTAB is not set.")          self.assertTrue(sb.getSolverMethod() == so.CGS, "CGS is not set.")
564          so.setSolverMethod(so.SSOR)          sb.setSolverMethod(so.BICGSTAB)
565          self.failUnless(so.getSolverMethod() == 7, "SSOR is not set.")          self.assertTrue(sb.getSolverMethod() == so.BICGSTAB, "BICGSTAB is not set.")
566          so.setSolverMethod(so.GMRES)          sb.setSolverMethod(so.GMRES)
567          self.failUnless(so.getSolverMethod() == 11, "GMRES is not set.")          self.assertTrue(sb.getSolverMethod() == so.GMRES, "GMRES is not set.")
568          so.setSolverMethod(so.PRES20)          sb.setSolverMethod(so.PRES20)
569          self.failUnless(so.getSolverMethod() == 12, "PRES20 is not set.")          self.assertTrue(sb.getSolverMethod() == so.PRES20, "PRES20 is not set.")
570          so.setSolverMethod(so.LUMPING)          sb.setSolverMethod(so.LUMPING)
571          self.failUnless(so.getSolverMethod() == 13, "LUMPING is not set.")          self.assertTrue(sb.getSolverMethod() == so.LUMPING, "LUMPING is not set.")
572          so.setSolverMethod(so.ITERATIVE)          sb.setSolverMethod(so.ITERATIVE)
573          self.failUnless(so.getSolverMethod() == 20, "ITERATIVE is not set.")          self.assertTrue(sb.getSolverMethod() == so.ITERATIVE, "ITERATIVE is not set.")
574          so.setSolverMethod(so.AMG)          sb.setSolverMethod(so.LSQR)
575          self.failUnless(so.getSolverMethod() == 22, "AMG is not set.")          self.assertTrue(sb.getSolverMethod() == so.LSQR, "LSQR is not set.")
576          so.setSolverMethod(so.NONLINEAR_GMRES)          sb.setSolverMethod(so.NONLINEAR_GMRES)
577          self.failUnless(so.getSolverMethod() == 25, "NONLINEAR_GMRES is not set.")          self.assertTrue(sb.getSolverMethod() == so.NONLINEAR_GMRES, "NONLINEAR_GMRES is not set.")
578          so.setSolverMethod(so.TFQMR)          sb.setSolverMethod(so.TFQMR)
579          self.failUnless(so.getSolverMethod() == 26, "TFQMR is not set.")          self.assertTrue(sb.getSolverMethod() == so.TFQMR, "TFQMR is not set.")
580          so.setSolverMethod(so.MINRES)          sb.setSolverMethod(so.MINRES)
581          self.failUnless(so.getSolverMethod() == 27, "MINRES is not set.")          self.assertTrue(sb.getSolverMethod() == so.MINRES, "MINRES is not set.")
582          so.setSolverMethod(so.GAUSS_SEIDEL)          sb.setSolverMethod(so.DEFAULT)
583          self.failUnless(so.getSolverMethod() == 28, "GAUSS_SEIDEL is not set.")          self.assertTrue(sb.getSolverMethod() == so.DEFAULT, "DEFAULT is not set.")
584          so.setSolverMethod(so.DEFAULT)  
585          self.failUnless(so.getSolverMethod() == 0, "DEFAULT is not set.")          self.assertTrue(sb.getPreconditioner() == so.JACOBI, "initial Preconditioner is wrong.")
586            self.assertRaises(ValueError,sb.setPreconditioner,-1)
587          self.failUnless(so.getPreconditioner() == 10, "initial Preconditioner is wrong.")          sb.setPreconditioner(so.ILU0)
588          self.failUnlessRaises(ValueError,so.setPreconditioner,-1)          self.assertTrue(sb.getPreconditioner() == so.ILU0, "ILU0 is not set.")
589          so.setPreconditioner(so.ILU0)          sb.setPreconditioner(so.ILUT)
590          self.failUnless(so.getPreconditioner() == 8, "ILU0 is not set.")          self.assertTrue(sb.getPreconditioner() == so.ILUT, "ILUT is not set.")
591          so.setPreconditioner(so.SSOR)          sb.setPreconditioner(so.JACOBI)
592          self.failUnless(so.getPreconditioner() == 7, "SSOR is not set.")          self.assertTrue(sb.getPreconditioner() == so.JACOBI, "JACOBI is not set.")
593          so.setPreconditioner(so.ILUT)          if getEscriptParamInt('DISABLE_AMG', 0):
594          self.failUnless(so.getPreconditioner() == 9, "ILUT is not set.")              print("AMG test disabled on MPI build")
595          so.setPreconditioner(so.JACOBI)          else:
596          self.failUnless(so.getPreconditioner() == 10, "JACOBI is not set.")              sb.setPreconditioner(so.AMG)
597          so.setPreconditioner(so.AMG)              self.assertTrue(sb.getPreconditioner() == so.AMG, "AMG is not set.")
598          self.failUnless(so.getPreconditioner() == 22, "AMG is not set.")          sb.setPreconditioner(so.REC_ILU)
599          so.setPreconditioner(so.REC_ILU)          self.assertTrue(sb.getPreconditioner() == so.REC_ILU, "REC_ILU is not set.")
600          self.failUnless(so.getPreconditioner() == 23, "REC_ILU is not set.")          sb.setPreconditioner(so.GAUSS_SEIDEL)
601          so.setPreconditioner(so.GAUSS_SEIDEL)          self.assertTrue(sb.getPreconditioner() == so.GAUSS_SEIDEL, "GAUSS_SEIDEL is not set.")
602          self.failUnless(so.getPreconditioner() == 28, "GAUSS_SEIDEL is not set.")          sb.setPreconditioner(so.RILU)
603          so.setPreconditioner(so.RILU)          self.assertTrue(sb.getPreconditioner() == so.RILU, "RILU is not set.")
604          self.failUnless(so.getPreconditioner() == 29, "RILU is not set.")          sb.setPreconditioner(so.AMLI)
605          so.setPreconditioner(so.NO_PRECONDITIONER)          self.assertTrue(sb.getPreconditioner() == so.AMLI, "AMLI is not set.")
606          self.failUnless(so.getPreconditioner() == 36, "NO_PRECONDITIONER is not set.")                  sb.setPreconditioner(so.NO_PRECONDITIONER)
607            self.assertTrue(sb.getPreconditioner() == so.NO_PRECONDITIONER, "NO_PRECONDITIONER is not set.")        
608          self.failUnless(so.getCoarsening() == 0, "initial Coarseningr is wrong.")  
609          self.failUnlessRaises(ValueError,so.setCoarsening,-1)          self.assertTrue(sb.getCoarsening() == so.DEFAULT, "initial Coarsening is wrong.")
610          so.setCoarsening(so.YAIR_SHAPIRA_COARSENING)          self.assertRaises(ValueError,sb.setCoarsening,-1)
611          self.failUnless(so.getCoarsening() == 33, "YAIR_SHAPIRA_COARSENING is not set.")          sb.setCoarsening(so.YAIR_SHAPIRA_COARSENING)
612          so.setCoarsening(so.RUGE_STUEBEN_COARSENING)          self.assertTrue(sb.getCoarsening() == so.YAIR_SHAPIRA_COARSENING, "YAIR_SHAPIRA_COARSENING is not set.")
613          self.failUnless(so.getCoarsening() == 34, "RUGE_STUEBEN_COARSENING is not set.")          sb.setCoarsening(so.RUGE_STUEBEN_COARSENING)
614          so.setCoarsening(so.AGGREGATION_COARSENING)          self.assertTrue(sb.getCoarsening() == so.RUGE_STUEBEN_COARSENING, "RUGE_STUEBEN_COARSENING is not set.")
615          self.failUnless(so.getCoarsening() == 35, "AGREGATION_COARSENING is not set.")          sb.setCoarsening(so.AGGREGATION_COARSENING)
616          so.setCoarsening(so.DEFAULT)          self.assertTrue(sb.getCoarsening() == so.AGGREGATION_COARSENING, "AGREGATION_COARSENING is not set.")
617          self.failUnless(so.getCoarsening() == 0, "DEFAULT is not set.")          sb.setCoarsening(so.STANDARD_COARSENING)
618            self.assertTrue(sb.getCoarsening() == so.STANDARD_COARSENING, "STANDARD_COARSENING is not set.")
619          self.failUnless(so.getDiagnostics("num_iter") == None, "initial num_iter is wrong.")          sb.setCoarsening(so.DEFAULT)
620          self.failUnless(so.getDiagnostics("num_inner_iter") == None, "initial num_inner_iter is wrong.")          self.assertTrue(sb.getCoarsening() == so.DEFAULT, "DEFAULT is not set.")
621          self.failUnless(so.getDiagnostics("time") == None, "initial time is wrong.")  
622          self.failUnless(so.getDiagnostics("set_up_time") == None, "initial set_up_time is wrong.")          self.assertTrue(sb.getDiagnostics("num_iter") == 0, "initial num_iter is wrong.")
623          self.failUnless(so.getDiagnostics("residual_norm") == None, "initial residual_norm is wrong.")          self.assertTrue(sb.getDiagnostics("num_inner_iter") == 0, "initial num_inner_iter is wrong.")
624          self.failUnless(so.getDiagnostics("converged") == None, "initial converged is wrong.")          self.assertTrue(sb.getDiagnostics("time") == 0, "initial time is wrong.")
625          self.failUnless(so.hasConverged() == None, "initial convergence flag is wrong.")          self.assertTrue(sb.getDiagnostics("set_up_time") == 0, "initial set_up_time is wrong.")
626          self.failUnless(so.getDiagnostics("cum_num_inner_iter") == 0, "initial cum_num_inner_iter is wrong.")          self.assertTrue(sb.getDiagnostics("residual_norm") == 0, "initial residual_norm is wrong.")
627          self.failUnless(so.getDiagnostics("cum_num_iter") == 0, "initial cum_num_iter is wrong.")          self.assertTrue(sb.getDiagnostics("converged") == 0, "initial converged is wrong.")
628          self.failUnless(so.getDiagnostics("cum_time") ==0, "initial cum_time is wrong.")          self.assertTrue(sb.hasConverged() == 0, "initial convergence flag is wrong.")
629          self.failUnless(so.getDiagnostics("cum_set_up_time") == 0, "initial cum_set_up_time is wrong.")          self.assertTrue(sb.getDiagnostics("cum_num_inner_iter") == 0, "initial cum_num_inner_iter is wrong.")
630            self.assertTrue(sb.getDiagnostics("cum_num_iter") == 0, "initial cum_num_iter is wrong.")
631          so._updateDiagnostics("num_iter",1)          self.assertTrue(sb.getDiagnostics("cum_time") ==0, "initial cum_time is wrong.")
632          so._updateDiagnostics("num_inner_iter",2)          self.assertTrue(sb.getDiagnostics("cum_set_up_time") == 0, "initial cum_set_up_time is wrong.")
633          so._updateDiagnostics("time",3)  
634          so._updateDiagnostics("set_up_time",4)          sb._updateDiagnostics("num_iter",1)
635          so._updateDiagnostics("residual_norm",5)          sb._updateDiagnostics("num_inner_iter",2)
636          so._updateDiagnostics("converged",True)          sb._updateDiagnostics("time",3)
637            sb._updateDiagnostics("set_up_time",4)
638          self.failUnless(so.getDiagnostics("num_iter") == 1, "num_iter is wrong.")          sb._updateDiagnostics("residual_norm",5)
639          self.failUnless(so.getDiagnostics("num_inner_iter") == 2, "num_inner_iter is wrong.")          sb._updateDiagnostics("converged",True)
640          self.failUnless(so.getDiagnostics("time") == 3, "time is wrong.")  
641          self.failUnless(so.getDiagnostics("set_up_time") == 4, "set_up_time is wrong.")          self.assertTrue(sb.getDiagnostics("num_iter") == 1, "num_iter is wrong.")
642          self.failUnless(so.getDiagnostics("residual_norm") == 5, "residual_norm is wrong.")          self.assertTrue(sb.getDiagnostics("num_inner_iter") == 2, "num_inner_iter is wrong.")
643          self.failUnless(so.getDiagnostics("converged"), "converged is wrong.")          self.assertTrue(sb.getDiagnostics("time") == 3, "time is wrong.")
644          self.failUnless(so.hasConverged(), "convergence flag is wrong.")          self.assertTrue(sb.getDiagnostics("set_up_time") == 4, "set_up_time is wrong.")
645          self.failUnless(so.getDiagnostics("cum_num_inner_iter") == 2, "cum_num_inner_iter is wrong.")          self.assertTrue(sb.getDiagnostics("residual_norm") == 5, "residual_norm is wrong.")
646          self.failUnless(so.getDiagnostics("cum_num_iter") == 1, "cum_num_iter is wrong.")          self.assertTrue(sb.getDiagnostics("converged"), "converged is wrong.")
647          self.failUnless(so.getDiagnostics("cum_time") ==3, "cum_time is wrong.")          self.assertTrue(sb.hasConverged(), "convergence flag is wrong.")
648          self.failUnless(so.getDiagnostics("cum_set_up_time") == 4, "cum_set_up_time is wrong.")            self.assertTrue(sb.getDiagnostics("cum_num_inner_iter") == 2, "cum_num_inner_iter is wrong.")
649            self.assertTrue(sb.getDiagnostics("cum_num_iter") == 1, "cum_num_iter is wrong.")
650            self.assertTrue(sb.getDiagnostics("cum_time") ==3, "cum_time is wrong.")
651            self.assertTrue(sb.getDiagnostics("cum_set_up_time") == 4, "cum_set_up_time is wrong.")  
652                    
653          so.resetDiagnostics()          sb.resetDiagnostics()
654          self.failUnless(so.getDiagnostics("num_iter") == None, "initial num_iter is wrong.")          self.assertTrue(sb.getDiagnostics("num_iter") == 0, "initial num_iter is wrong.")
655          self.failUnless(so.getDiagnostics("num_inner_iter") == None, "initial num_inner_iter is wrong.")          self.assertTrue(sb.getDiagnostics("num_inner_iter") == 0, "initial num_inner_iter is wrong.")
656          self.failUnless(so.getDiagnostics("time") == None, "initial time is wrong.")          self.assertTrue(sb.getDiagnostics("time") == 0, "initial time is wrong.")
657          self.failUnless(so.getDiagnostics("set_up_time") == None, "initial set_up_time is wrong.")          self.assertTrue(sb.getDiagnostics("set_up_time") == 0, "initial set_up_time is wrong.")
658          self.failUnless(so.getDiagnostics("residual_norm") == None, "initial residual_norm is wrong.")          self.assertTrue(sb.getDiagnostics("residual_norm") == 0, "initial residual_norm is wrong.")
659          self.failUnless(so.getDiagnostics("converged") == None, "initial converged is wrong.")          self.assertTrue(sb.getDiagnostics("converged") == 0, "initial converged is wrong.")
660          self.failUnless(so.hasConverged() == None, "initial convergence flag is wrong")                self.assertTrue(sb.hasConverged() == 0, "initial convergence flag is wrong")      
661          self.failUnless(so.getDiagnostics("cum_num_inner_iter") == 2, "cum_num_inner_iter is wrong.")          self.assertTrue(sb.getDiagnostics("cum_num_inner_iter") == 2, "cum_num_inner_iter is wrong.")
662          self.failUnless(so.getDiagnostics("cum_num_iter") == 1, "cum_num_iter is wrong.")          self.assertTrue(sb.getDiagnostics("cum_num_iter") == 1, "cum_num_iter is wrong.")
663          self.failUnless(so.getDiagnostics("cum_time") ==3, "cum_time is wrong.")          self.assertTrue(sb.getDiagnostics("cum_time") == 3, "cum_time is wrong.")
664          self.failUnless(so.getDiagnostics("cum_set_up_time") == 4, "cum_set_up_time is wrong.")          self.assertTrue(sb.getDiagnostics("cum_set_up_time") == 4, "cum_set_up_time is wrong.")
665    
666          so._updateDiagnostics("num_iter",10)          sb._updateDiagnostics("num_iter",10)
667          so._updateDiagnostics("num_inner_iter",20)          sb._updateDiagnostics("num_inner_iter",20)
668          so._updateDiagnostics("time",30)          sb._updateDiagnostics("time",30)
669          so._updateDiagnostics("set_up_time",40)          sb._updateDiagnostics("set_up_time",40)
670          so._updateDiagnostics("residual_norm",50)          sb._updateDiagnostics("residual_norm",50)
671          so._updateDiagnostics("converged",False)          sb._updateDiagnostics("converged",False)
672    
673          self.failUnless(so.getDiagnostics("num_iter") == 10, "num_iter is wrong.")          self.assertTrue(sb.getDiagnostics("num_iter") == 10, "num_iter is wrong.")
674          self.failUnless(so.getDiagnostics("num_inner_iter") == 20, "num_inner_iter is wrong.")          self.assertTrue(sb.getDiagnostics("num_inner_iter") == 20, "num_inner_iter is wrong.")
675          self.failUnless(so.getDiagnostics("time") == 30, "time is wrong.")          self.assertTrue(sb.getDiagnostics("time") == 30, "time is wrong.")
676          self.failUnless(so.getDiagnostics("set_up_time") == 40, "set_up_time is wrong.")          self.assertTrue(sb.getDiagnostics("set_up_time") == 40, "set_up_time is wrong.")
677          self.failUnless(so.getDiagnostics("residual_norm") == 50, "residual_norm is wrong.")          self.assertTrue(sb.getDiagnostics("residual_norm") == 50, "residual_norm is wrong.")
678          self.failUnless(not so.getDiagnostics("converged"), "converged is wrong.")          self.assertTrue(not sb.getDiagnostics("converged"), "converged is wrong.")
679          self.failUnless(not so.hasConverged(), "convergence flag is wrong.")          self.assertTrue(not sb.hasConverged(), "convergence flag is wrong.")
680          self.failUnless(so.getDiagnostics("cum_num_inner_iter") == 22, "cum_num_inner_iter is wrong.")          self.assertTrue(sb.getDiagnostics("cum_num_inner_iter") == 22, "cum_num_inner_iter is wrong.")
681          self.failUnless(so.getDiagnostics("cum_num_iter") == 11, "cum_num_iter is wrong.")          self.assertTrue(sb.getDiagnostics("cum_num_iter") == 11, "cum_num_iter is wrong.")
682          self.failUnless(so.getDiagnostics("cum_time") ==33, "cum_time is wrong.")          self.assertTrue(sb.getDiagnostics("cum_time") ==33, "cum_time is wrong.")
683          self.failUnless(so.getDiagnostics("cum_set_up_time") == 44, "cum_set_up_time is wrong.")            self.assertTrue(sb.getDiagnostics("cum_set_up_time") == 44, "cum_set_up_time is wrong.")  
684    
685          so.resetDiagnostics(all=True)          sb.resetDiagnostics(all=True)
686          self.failUnless(so.getDiagnostics("num_iter") == None, "initial num_iter is wrong.")          self.assertTrue(sb.getDiagnostics("num_iter") == 0, "initial num_iter is wrong.")
687          self.failUnless(so.getDiagnostics("num_inner_iter") == None, "initial num_inner_iter is wrong.")          self.assertTrue(sb.getDiagnostics("num_inner_iter") == 0, "initial num_inner_iter is wrong.")
688          self.failUnless(so.getDiagnostics("time") == None, "initial time is wrong.")          self.assertTrue(sb.getDiagnostics("time") == 0, "initial time is wrong.")
689          self.failUnless(so.getDiagnostics("set_up_time") == None, "initial set_up_time is wrong.")          self.assertTrue(sb.getDiagnostics("set_up_time") == 0, "initial set_up_time is wrong.")
690          self.failUnless(so.getDiagnostics("residual_norm") == None, "initial residual_norm is wrong.")          self.assertTrue(sb.getDiagnostics("residual_norm") == 0, "initial residual_norm is wrong.")
691          self.failUnless(so.getDiagnostics("converged") == None, "initial converged is wrong.")          self.assertTrue(sb.getDiagnostics("converged") == 0, "initial converged is wrong.")
692          self.failUnless(so.hasConverged() == None, "initial convergence flag is wrong.")          self.assertTrue(sb.hasConverged() == 0, "initial convergence flag is wrong.")
693          self.failUnless(so.getDiagnostics("cum_num_inner_iter") == 0, "initial cum_num_inner_iter is wrong.")          self.assertTrue(sb.getDiagnostics("cum_num_inner_iter") == 0, "initial cum_num_inner_iter is wrong.")
694          self.failUnless(so.getDiagnostics("cum_num_iter") == 0, "initial cum_num_iter is wrong.")          self.assertTrue(sb.getDiagnostics("cum_num_iter") == 0, "initial cum_num_iter is wrong.")
695          self.failUnless(so.getDiagnostics("cum_time") ==0, "initial cum_time is wrong.")          self.assertTrue(sb.getDiagnostics("cum_time") ==0, "initial cum_time is wrong.")
696          self.failUnless(so.getDiagnostics("cum_set_up_time") == 0, "initial cum_set_up_time is wrong.")          self.assertTrue(sb.getDiagnostics("cum_set_up_time") == 0, "initial cum_set_up_time is wrong.")
697                    
698      def test_setCoefficient_WithIllegalFunctionSpace(self):      def test_setCoefficient_WithIllegalFunctionSpace(self):
699          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
700          self.failUnlessRaises(IllegalCoefficientFunctionSpace, mypde.setValue, C=Vector(0.,FunctionOnBoundary(self.domain)))          self.assertRaises(IllegalCoefficientFunctionSpace, mypde.setValue, C=Vector(0.,FunctionOnBoundary(self.domain)))
701    
702      def test_setCoefficient_WithWrongName(self):      def test_setCoefficient_WithWrongName(self):
703          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
704          self.failUnlessRaises(IllegalCoefficient, mypde.setValue, ROMA=0.)          self.assertRaises(IllegalCoefficient, mypde.setValue, ROMA=0.)
705    
706      def test_resetCoefficient_WithWrongShape(self):      def test_resetCoefficient_WithWrongShape(self):
707          mypde=LinearPDE(self.domain,numEquations=2,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numEquations=2,debug=self.DEBUG)
708          self.failUnlessRaises(IllegalCoefficientValue, mypde.setValue, C=0.)          self.assertRaises(IllegalCoefficientValue, mypde.setValue, C=0.)
709    
710      def test_reducedOn(self):      def test_reducedOn(self):
711          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
# Line 955  class Test_LinearPDE_noLumping(Test_line Line 713  class Test_LinearPDE_noLumping(Test_line
713          mypde.setReducedOrderOn()          mypde.setReducedOrderOn()
714          mypde.setValue(A=kronecker(self.domain),D=x[0],Y=x[0])          mypde.setValue(A=kronecker(self.domain),D=x[0],Y=x[0])
715          u=mypde.getSolution()          u=mypde.getSolution()
716          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
717    
718      def test_attemptToChangeOrderAfterDefinedCoefficient(self):      def test_attemptToChangeOrderAfterDefinedCoefficient(self):
719          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
720          mypde.setValue(D=1.)          mypde.setValue(D=1.)
721          self.failUnlessRaises(RuntimeError, mypde.setReducedOrderOn)          self.assertRaises(RuntimeError, mypde.setReducedOrderOn)
722    
723      def test_reducedOnConfig(self):      def test_reducedOnConfig(self):
724          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
725          mypde.setReducedOrderOn()          mypde.setReducedOrderOn()
726          self.failUnlessEqual((mypde.getFunctionSpaceForSolution(), mypde.getFunctionSpaceForEquation()),(ReducedSolution(self.domain),ReducedSolution(self.domain)),"reduced function spaces expected.")          self.assertEqual((mypde.getFunctionSpaceForSolution(), mypde.getFunctionSpaceForEquation()),(ReducedSolution(self.domain),ReducedSolution(self.domain)),"reduced function spaces expected.")
727      #      #
728      #  set coefficients for scalars:      #  set coefficients for scalars:
729      #      #
# Line 974  class Test_LinearPDE_noLumping(Test_line Line 732  class Test_LinearPDE_noLumping(Test_line
732          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
733          mypde.setValue(A=numpy.ones((d,d)))          mypde.setValue(A=numpy.ones((d,d)))
734          coeff=mypde.getCoefficient("A")          coeff=mypde.getCoefficient("A")
735          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,d),Function(self.domain),1,1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,d),Function(self.domain),1,1))
736    
737            mypde.resetRightHandSideCoefficients()
738            self.assertFalse(mypde.getCoefficient("A").isEmpty(),"A is empty after reset of right hand side coefficients")
739    
740      def test_setCoefficient_B_Scalar(self):      def test_setCoefficient_B_Scalar(self):
741          d=self.domain.getDim()          d=self.domain.getDim()
742          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
743          mypde.setValue(B=numpy.ones((d,)))          mypde.setValue(B=numpy.ones((d,)))
744          coeff=mypde.getCoefficient("B")          coeff=mypde.getCoefficient("B")
745          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,),Function(self.domain),1,1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,),Function(self.domain),1,1))
746    
747            mypde.resetRightHandSideCoefficients()
748            self.assertFalse(mypde.getCoefficient("B").isEmpty(),"B is empty after reset of right hand side coefficients")
749    
750      def test_setCoefficient_C_Scalar(self):      def test_setCoefficient_C_Scalar(self):
751          d=self.domain.getDim()          d=self.domain.getDim()
752          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
753          mypde.setValue(C=numpy.ones((d,)))          mypde.setValue(C=numpy.ones((d,)))
754          coeff=mypde.getCoefficient("C")          coeff=mypde.getCoefficient("C")
755          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,),Function(self.domain),1,1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,),Function(self.domain),1,1))
756    
757            mypde.resetRightHandSideCoefficients()
758            self.assertFalse(mypde.getCoefficient("C").isEmpty(),"C is empty after reset of right hand side coefficients")
759    
760      def test_setCoefficient_D_Scalar(self):      def test_setCoefficient_D_Scalar(self):
761          d=self.domain.getDim()          d=self.domain.getDim()
762          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
763          mypde.setValue(D=1.)          mypde.setValue(D=1.)
764          coeff=mypde.getCoefficient("D")          coeff=mypde.getCoefficient("D")
765          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),Function(self.domain),1,1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),Function(self.domain),1,1))
766    
767            mypde.resetRightHandSideCoefficients()
768            self.assertFalse(mypde.getCoefficient("D").isEmpty(),"D is empty after reset of right hand side coefficients")
769    
770      def test_setCoefficient_X_Scalar(self):      def test_setCoefficient_X_Scalar(self):
771          d=self.domain.getDim()          d=self.domain.getDim()
772          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
773          mypde.setValue(X=numpy.ones((d,)))          mypde.setValue(X=numpy.ones((d,)))
774          coeff=mypde.getCoefficient("X")          coeff=mypde.getCoefficient("X")
775          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((d,),Function(self.domain),1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((d,),Function(self.domain),1))
776    
777            mypde.resetRightHandSideCoefficients()
778            self.assertTrue(mypde.getCoefficient("X").isEmpty(),"X is not empty after reset of right hand side coefficients")
779    
780      def test_setCoefficient_Y_Scalar(self):      def test_setCoefficient_Y_Scalar(self):
781          d=self.domain.getDim()          d=self.domain.getDim()
782          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
783          mypde.setValue(Y=1.)          mypde.setValue(Y=1.)
784          coeff=mypde.getCoefficient("Y")          coeff=mypde.getCoefficient("Y")
785          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),Function(self.domain),1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),Function(self.domain),1))
786    
787            mypde.resetRightHandSideCoefficients()
788            self.assertTrue(mypde.getCoefficient("Y").isEmpty(),"Y is not empty after reset of right hand side coefficients")
789    
790      def test_setCoefficient_y_Scalar(self):      def test_setCoefficient_y_Scalar(self):
791          d=self.domain.getDim()          d=self.domain.getDim()
792          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
793          mypde.setValue(y=1.)          mypde.setValue(y=1.)
794          coeff=mypde.getCoefficient("y")          coeff=mypde.getCoefficient("y")
795          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),FunctionOnBoundary(self.domain),1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),FunctionOnBoundary(self.domain),1))
796    
797            mypde.resetRightHandSideCoefficients()
798            self.assertTrue(mypde.getCoefficient("y").isEmpty(),"y is not empty after reset of right hand side coefficients")
799    
800      def test_setCoefficient_d_Scalar(self):      def test_setCoefficient_d_Scalar(self):
801          d=self.domain.getDim()          d=self.domain.getDim()
802          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
803          mypde.setValue(d=1.)          mypde.setValue(d=1.)
804          coeff=mypde.getCoefficient("d")          coeff=mypde.getCoefficient("d")
805          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),FunctionOnBoundary(self.domain),1,1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),FunctionOnBoundary(self.domain),1,1))
806    
807            mypde.resetRightHandSideCoefficients()
808            self.assertFalse(mypde.getCoefficient("d").isEmpty(),"d is empty after reset of right hand side coefficients")
809    
810      def test_setCoefficient_d_contact_Scalar(self):      def test_setCoefficient_d_contact_Scalar(self):
811          d=self.domain.getDim()          if self.domain.supportsContactElements():
812          mypde=LinearPDE(self.domain,debug=self.DEBUG)              d=self.domain.getDim()
813          mypde.setValue(d_contact=1.)              mypde=LinearPDE(self.domain,debug=self.DEBUG)
814          coeff=mypde.getCoefficient("d_contact")              mypde.setValue(d_contact=1.)
815          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),FunctionOnContactZero(self.domain),1,1))              coeff=mypde.getCoefficient("d_contact")
816                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),FunctionOnContactZero(self.domain),1,1))
817    
818                mypde.resetRightHandSideCoefficients()
819                self.assertFalse(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is empty after reset of right hand side coefficients")
820    
821      def test_setCoefficient_y_contact_Scalar(self):      def test_setCoefficient_y_contact_Scalar(self):
822          d=self.domain.getDim()          d=self.domain.getDim()
823          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          if self.domain.supportsContactElements():
824          mypde.setValue(y_contact=1.)              mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
825          coeff=mypde.getCoefficient("y_contact")              mypde.setValue(y_contact=1.)
826          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),FunctionOnContactZero(self.domain),1))              coeff=mypde.getCoefficient("y_contact")
827                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),FunctionOnContactZero(self.domain),1))
828    
829                mypde.resetRightHandSideCoefficients()
830                self.assertTrue(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty after reset of right hand side coefficients")
831    
832      def test_setCoefficient_A_reduced_Scalar(self):      def test_setCoefficient_A_reduced_Scalar(self):
833          d=self.domain.getDim()          d=self.domain.getDim()
834          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
835          mypde.setValue(A_reduced=numpy.ones((d,d)))          mypde.setValue(A_reduced=numpy.ones((d,d)))
836          coeff=mypde.getCoefficient("A_reduced")          coeff=mypde.getCoefficient("A_reduced")
837          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,d),ReducedFunction(self.domain),1,1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,d),ReducedFunction(self.domain),1,1))
838    
839            mypde.resetRightHandSideCoefficients()
840            self.assertFalse(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is empty after reset of right hand side coefficients")
841    
842      def test_setCoefficient_B_reduced_Scalar(self):      def test_setCoefficient_B_reduced_Scalar(self):
843          d=self.domain.getDim()          d=self.domain.getDim()
844          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
845          mypde.setValue(B_reduced=numpy.ones((d,)))          mypde.setValue(B_reduced=numpy.ones((d,)))
846          coeff=mypde.getCoefficient("B_reduced")          coeff=mypde.getCoefficient("B_reduced")
847          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1,1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1,1))
848    
849            mypde.resetRightHandSideCoefficients()
850            self.assertFalse(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is empty after reset of right hand side coefficients")
851    
852      def test_setCoefficient_C_reduced_Scalar(self):      def test_setCoefficient_C_reduced_Scalar(self):
853          d=self.domain.getDim()          d=self.domain.getDim()
854          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
855          mypde.setValue(C_reduced=numpy.ones((d,)))          mypde.setValue(C_reduced=numpy.ones((d,)))
856          coeff=mypde.getCoefficient("C_reduced")          coeff=mypde.getCoefficient("C_reduced")
857          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1,1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1,1))
858    
859            mypde.resetRightHandSideCoefficients()
860            self.assertFalse(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is empty after reset of right hand side coefficients")
861    
862      def test_setCoefficient_D_reduced_Scalar(self):      def test_setCoefficient_D_reduced_Scalar(self):
863          d=self.domain.getDim()          d=self.domain.getDim()
864          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
865          mypde.setValue(D_reduced=1.)          mypde.setValue(D_reduced=1.)
866          coeff=mypde.getCoefficient("D_reduced")          coeff=mypde.getCoefficient("D_reduced")
867          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),ReducedFunction(self.domain),1,1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),ReducedFunction(self.domain),1,1))
868    
869            mypde.resetRightHandSideCoefficients()
870            self.assertFalse(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is empty after reset of right hand side coefficients")
871    
872      def test_setCoefficient_X_reduced_Scalar(self):      def test_setCoefficient_X_reduced_Scalar(self):
873          d=self.domain.getDim()          d=self.domain.getDim()
874          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
875          mypde.setValue(X_reduced=numpy.ones((d,)))          mypde.setValue(X_reduced=numpy.ones((d,)))
876          coeff=mypde.getCoefficient("X_reduced")          coeff=mypde.getCoefficient("X_reduced")
877          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1))
878    
879            mypde.resetRightHandSideCoefficients()
880            self.assertTrue(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty after reset of right hand side coefficients")
881    
882      def test_setCoefficient_Y_reduced_Scalar(self):      def test_setCoefficient_Y_reduced_Scalar(self):
883          d=self.domain.getDim()          d=self.domain.getDim()
884          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
885          mypde.setValue(Y_reduced=1.)          mypde.setValue(Y_reduced=1.)
886          coeff=mypde.getCoefficient("Y_reduced")          coeff=mypde.getCoefficient("Y_reduced")
887          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),ReducedFunction(self.domain),1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),ReducedFunction(self.domain),1))
888    
889            mypde.resetRightHandSideCoefficients()
890            self.assertTrue(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty after reset of right hand side coefficients")
891    
892      def test_setCoefficient_y_reduced_Scalar(self):      def test_setCoefficient_y_reduced_Scalar(self):
893          d=self.domain.getDim()          d=self.domain.getDim()
894          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
895          mypde.setValue(y_reduced=1.)          mypde.setValue(y_reduced=1.)
896          coeff=mypde.getCoefficient("y_reduced")          coeff=mypde.getCoefficient("y_reduced")
897          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),ReducedFunctionOnBoundary(self.domain),1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),ReducedFunctionOnBoundary(self.domain),1))
898    
899            mypde.resetRightHandSideCoefficients()
900            self.assertTrue(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty after reset of right hand side coefficients")
901    
902      def test_setCoefficient_d_reduced_Scalar(self):      def test_setCoefficient_d_reduced_Scalar(self):
903          d=self.domain.getDim()          d=self.domain.getDim()
904          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
905          mypde.setValue(d_reduced=1.)          mypde.setValue(d_reduced=1.)
906          coeff=mypde.getCoefficient("d_reduced")          coeff=mypde.getCoefficient("d_reduced")
907          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),ReducedFunctionOnBoundary(self.domain),1,1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),ReducedFunctionOnBoundary(self.domain),1,1))
908    
909            mypde.resetRightHandSideCoefficients()
910            self.assertFalse(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty after reset of right hand side coefficients")
911    
912      def test_setCoefficient_d_contact_reduced_Scalar(self):      def test_setCoefficient_d_contact_reduced_Scalar(self):
913          d=self.domain.getDim()          if self.domain.supportsContactElements():
914          mypde=LinearPDE(self.domain,debug=self.DEBUG)              d=self.domain.getDim()
915          mypde.setValue(d_contact_reduced=1.)              mypde=LinearPDE(self.domain,debug=self.DEBUG)
916          coeff=mypde.getCoefficient("d_contact_reduced")              mypde.setValue(d_contact_reduced=1.)
917          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),ReducedFunctionOnContactZero(self.domain),1,1))              coeff=mypde.getCoefficient("d_contact_reduced")
918                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),ReducedFunctionOnContactZero(self.domain),1,1))
919    
920                mypde.resetRightHandSideCoefficients()
921                self.assertFalse(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty after reset of right hand side coefficients")
922      def test_setCoefficient_y_contact_reduced_Scalar(self):      def test_setCoefficient_y_contact_reduced_Scalar(self):
923          d=self.domain.getDim()          if self.domain.supportsContactElements():
924          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)              d=self.domain.getDim()
925          mypde.setValue(y_contact_reduced=1.)              mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
926          coeff=mypde.getCoefficient("y_contact_reduced")              mypde.setValue(y_contact_reduced=1.)
927          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),ReducedFunctionOnContactZero(self.domain),1))              coeff=mypde.getCoefficient("y_contact_reduced")
928                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),ReducedFunctionOnContactZero(self.domain),1))
929    
930                mypde.resetRightHandSideCoefficients()
931                self.assertTrue(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty after reset of right hand side coefficients")
932    
933      def test_setCoefficient_r_Scalar(self):      def test_setCoefficient_r_Scalar(self):
934          d=self.domain.getDim()          d=self.domain.getDim()
935          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
936          mypde.setValue(r=1.)          mypde.setValue(r=1.)
937          coeff=mypde.getCoefficient("r")          coeff=mypde.getCoefficient("r")
938          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((),Solution(self.domain),1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((),Solution(self.domain),1))
939    
940            mypde.resetRightHandSideCoefficients()
941            self.assertTrue(mypde.getCoefficient("r").isEmpty(),"r is not empty after reset of right hand side coefficients")
942    
943      def test_setCoefficient_q_Scalar(self):      def test_setCoefficient_q_Scalar(self):
944          d=self.domain.getDim()          d=self.domain.getDim()
945          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
946          mypde.setValue(q=1.)          mypde.setValue(q=1.)
947          coeff=mypde.getCoefficient("q")          coeff=mypde.getCoefficient("q")
948          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((),Solution(self.domain),1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((),Solution(self.domain),1))
949    
950            mypde.resetRightHandSideCoefficients()
951            self.assertFalse(mypde.getCoefficient("q").isEmpty(),"q is empty after reset of right hand side coefficients")
952    
953      def test_setCoefficient_r_Scalar_reducedOn(self):      def test_setCoefficient_r_Scalar_reducedOn(self):
954          d=self.domain.getDim()          d=self.domain.getDim()
955          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
956          mypde.setReducedOrderOn()          mypde.setReducedOrderOn()
957          mypde.setValue(r=1.)          mypde.setValue(r=1.)
958          coeff=mypde.getCoefficient("r")          coeff=mypde.getCoefficient("r")
959          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((),ReducedSolution(self.domain),1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((),ReducedSolution(self.domain),1))
960    
961            mypde.resetRightHandSideCoefficients()
962            self.assertTrue(mypde.getCoefficient("r").isEmpty(),"r is not empty after reset of right hand side coefficients")
963    
964      def test_setCoefficient_q_Scalar_reducedOn(self):      def test_setCoefficient_q_Scalar_reducedOn(self):
965          d=self.domain.getDim()          d=self.domain.getDim()
966          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
967          mypde.setReducedOrderOn()          mypde.setReducedOrderOn()
968          mypde.setValue(q=1.)          mypde.setValue(q=1.)
969          coeff=mypde.getCoefficient("q")          coeff=mypde.getCoefficient("q")
970          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((),ReducedSolution(self.domain),1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((),ReducedSolution(self.domain),1))
971    
972            mypde.resetRightHandSideCoefficients()
973            self.assertFalse(mypde.getCoefficient("q").isEmpty(),"q is empty after reset of right hand side coefficients")
974    
975      def test_setCoefficient_A_reduced_Scalar_usingA(self):      def test_setCoefficient_A_reduced_Scalar_usingA(self):
976          d=self.domain.getDim()          d=self.domain.getDim()
977          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
978          mypde.setValue(A=Data(numpy.ones((d,d)),ReducedFunction(self.domain)))          mypde.setValue(A=Data(numpy.ones((d,d)),ReducedFunction(self.domain)))
979          coeff=mypde.getCoefficient("A_reduced")          if self.specialInterpolationSupported():
980          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,d),ReducedFunction(self.domain),1,1))              coeff_name='A'
981                FS=Function
982            else:
983                coeff_name='A_reduced'
984                FS=ReducedFunction
985            coeff=mypde.getCoefficient(coeff_name)
986            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,d),FS(self.domain),1,1))
987            mypde.resetRightHandSideCoefficients()
988            self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
989    
990      def test_setCoefficient_B_reduced_Scalar_usingB(self):      def test_setCoefficient_B_reduced_Scalar_usingB(self):
991          d=self.domain.getDim()          d=self.domain.getDim()
992          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
993          mypde.setValue(B=Data(numpy.ones((d,)),ReducedFunction(self.domain)))          mypde.setValue(B=Data(numpy.ones((d,)),ReducedFunction(self.domain)))
994          coeff=mypde.getCoefficient("B_reduced")          if self.specialInterpolationSupported():
995          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1,1))              coeff_name='B'
996                FS=Function
997            else:
998                coeff_name='B_reduced'
999                FS=ReducedFunction
1000            coeff=mypde.getCoefficient(coeff_name)
1001            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,),FS(self.domain),1,1))
1002            mypde.resetRightHandSideCoefficients()
1003            self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
1004    
1005      def test_setCoefficient_C_reduced_Scalar_usingC(self):      def test_setCoefficient_C_reduced_Scalar_usingC(self):
1006          d=self.domain.getDim()          d=self.domain.getDim()
1007          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1008          mypde.setValue(C=Data(numpy.ones((d,)),ReducedFunction(self.domain)))          mypde.setValue(C=Data(numpy.ones((d,)),ReducedFunction(self.domain)))
1009          coeff=mypde.getCoefficient("C_reduced")          if self.specialInterpolationSupported():
1010          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1,1))              coeff_name='C'
1011                FS=Function
1012            else:
1013                coeff_name='C_reduced'
1014                FS=ReducedFunction
1015            coeff=mypde.getCoefficient(coeff_name)
1016            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,),FS(self.domain),1,1))
1017            mypde.resetRightHandSideCoefficients()
1018            self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
1019    
1020      def test_setCoefficient_D_reduced_Scalar_usingD(self):      def test_setCoefficient_D_reduced_Scalar_usingD(self):
1021          d=self.domain.getDim()          d=self.domain.getDim()
1022          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1023          mypde.setValue(D=Scalar(1.,ReducedFunction(self.domain)))          mypde.setValue(D=Scalar(1.,ReducedFunction(self.domain)))
1024          coeff=mypde.getCoefficient("D_reduced")          if self.specialInterpolationSupported():
1025          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),ReducedFunction(self.domain),1,1))              coeff_name='D'
1026                FS=Function
1027            else:
1028                coeff_name='D_reduced'
1029                FS=ReducedFunction
1030            coeff=mypde.getCoefficient(coeff_name)
1031            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),FS(self.domain),1,1))
1032            mypde.resetRightHandSideCoefficients()
1033            self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
1034    
1035      def test_setCoefficient_X_reduced_Scalar_usingX(self):      def test_setCoefficient_X_reduced_Scalar_usingX(self):
1036          d=self.domain.getDim()          d=self.domain.getDim()
1037          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1038          mypde.setValue(X_reduced=Data(numpy.ones((d,)),ReducedFunction(self.domain)))          mypde.setValue(X_reduced=Data(numpy.ones((d,)),ReducedFunction(self.domain)))
1039          coeff=mypde.getCoefficient("X_reduced")          coeff=mypde.getCoefficient('X_reduced')
1040          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1))
1041            mypde.resetRightHandSideCoefficients()
1042            self.assertTrue(mypde.getCoefficient('X_reduced').isEmpty(),"X_reduced is not empty after reset of right hand side coefficients")
1043    
1044      def test_setCoefficient_Y_reduced_Scalar_usingY(self):      def test_setCoefficient_Y_reduced_Scalar_usingY(self):
1045          d=self.domain.getDim()          d=self.domain.getDim()
1046          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1047          mypde.setValue(Y=Scalar(1.,ReducedFunction(self.domain)))          mypde.setValue(Y=Scalar(1.,ReducedFunction(self.domain)))
1048          coeff=mypde.getCoefficient("Y_reduced")          if self.specialInterpolationSupported():
1049          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),ReducedFunction(self.domain),1))              coeff_name='Y'
1050                FS=Function
1051            else:
1052                coeff_name='Y_reduced'
1053                FS=ReducedFunction
1054            coeff=mypde.getCoefficient(coeff_name)
1055            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),FS(self.domain),1))
1056            mypde.resetRightHandSideCoefficients()
1057            self.assertTrue(mypde.getCoefficient(coeff_name).isEmpty(),"%s is not empty after reset of right hand side coefficients"%coeff_name)
1058    
1059      def test_setCoefficient_y_reduced_Scalar_using_y(self):      def test_setCoefficient_y_reduced_Scalar_using_y(self):
1060          d=self.domain.getDim()          d=self.domain.getDim()
1061          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1062          mypde.setValue(y=Scalar(1.,ReducedFunctionOnBoundary(self.domain)))          mypde.setValue(y=Scalar(1.,ReducedFunctionOnBoundary(self.domain)))
1063          coeff=mypde.getCoefficient("y_reduced")          if self.specialInterpolationSupported():
1064          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),ReducedFunctionOnBoundary(self.domain),1))              coeff_name='y'
1065                FS=FunctionOnBoundary
1066            else:
1067                coeff_name='y_reduced'
1068                FS=ReducedFunctionOnBoundary
1069            coeff=mypde.getCoefficient(coeff_name)
1070            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),FS(self.domain),1))
1071            mypde.resetRightHandSideCoefficients()
1072            self.assertTrue(mypde.getCoefficient(coeff_name).isEmpty(),"%s is not empty after reset of right hand side coefficients"%coeff_name)
1073    
1074      def test_setCoefficient_d_reduced_Scalar_using_d(self):      def test_setCoefficient_d_reduced_Scalar_using_d(self):
1075          d=self.domain.getDim()          d=self.domain.getDim()
1076          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1077          mypde.setValue(d=Scalar(1.,ReducedFunctionOnBoundary(self.domain)))          mypde.setValue(d=Scalar(1.,ReducedFunctionOnBoundary(self.domain)))
1078          coeff=mypde.getCoefficient("d_reduced")          if self.specialInterpolationSupported():
1079          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),ReducedFunctionOnBoundary(self.domain),1,1))              coeff_name='d'
1080                FS=FunctionOnBoundary
1081            else:
1082                coeff_name='d_reduced'
1083                FS=ReducedFunctionOnBoundary
1084            coeff=mypde.getCoefficient(coeff_name)
1085            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),FS(self.domain),1,1))
1086            mypde.resetRightHandSideCoefficients()
1087            self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
1088    
1089      def test_setCoefficient_d_contact_reduced_Scalar_using_d_contact(self):      def test_setCoefficient_d_contact_reduced_Scalar_using_d_contact(self):
1090          d=self.domain.getDim()          if self.domain.supportsContactElements():
1091          mypde=LinearPDE(self.domain,debug=self.DEBUG)              d=self.domain.getDim()
1092          mypde.setValue(d_contact=Scalar(1.,ReducedFunctionOnContactZero(self.domain)))              mypde=LinearPDE(self.domain,debug=self.DEBUG)
1093          coeff=mypde.getCoefficient("d_contact_reduced")              mypde.setValue(d_contact=Scalar(1.,ReducedFunctionOnContactZero(self.domain)))
1094          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),ReducedFunctionOnContactZero(self.domain),1,1))              if self.specialInterpolationSupported():
1095                    coeff_name='d_contact'
1096                    FS=FunctionOnContactZero
1097                else:
1098                    coeff_name='d_contact_reduced'
1099                    FS=ReducedFunctionOnContactZero
1100                coeff=mypde.getCoefficient(coeff_name)
1101                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),FS(self.domain),1,1))
1102                mypde.resetRightHandSideCoefficients()
1103                self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
1104    
1105      def test_setCoefficient_y_contact_reduced_Scalar_using_y_contact(self):      def test_setCoefficient_y_contact_reduced_Scalar_using_y_contact(self):
1106          d=self.domain.getDim()          if self.domain.supportsContactElements():
1107          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)              d=self.domain.getDim()
1108          mypde.setValue(y_contact=Scalar(1.,ReducedFunctionOnContactZero(self.domain)))              mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1109          coeff=mypde.getCoefficient("y_contact_reduced")              mypde.setValue(y_contact=Scalar(1.,ReducedFunctionOnContactZero(self.domain)))
1110          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),ReducedFunctionOnContactZero(self.domain),1))              if self.specialInterpolationSupported():
1111                    coeff_name='y_contact'
1112                    FS=FunctionOnContactZero
1113                else:
1114                    coeff_name='y_contact_reduced'
1115                    FS=ReducedFunctionOnContactZero
1116                coeff=mypde.getCoefficient(coeff_name)
1117                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),FS(self.domain),1))
1118                mypde.resetRightHandSideCoefficients()
1119                self.assertTrue(mypde.getCoefficient(coeff_name).isEmpty(),"%s is not empty after reset of right hand side coefficients"%coeff_name)
1120    
1121      #      #
1122      #  set coefficients for systems:      #  set coefficients for systems:
1123      #      #
# Line 1184  class Test_LinearPDE_noLumping(Test_line Line 1126  class Test_LinearPDE_noLumping(Test_line
1126          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1127          mypde.setValue(A=numpy.ones((self.N,d,self.N,d)))          mypde.setValue(A=numpy.ones((self.N,d,self.N,d)))
1128          coeff=mypde.getCoefficient("A")          coeff=mypde.getCoefficient("A")
1129          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,d,self.N,d),Function(self.domain),self.N,self.N))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,d,self.N,d),Function(self.domain),self.N,self.N))
1130            mypde.resetRightHandSideCoefficients()
1131            self.assertFalse(mypde.getCoefficient("A").isEmpty(),"A is empty after reset of right hand side coefficients")
1132            
1133      def test_setCoefficient_B_System(self):      def test_setCoefficient_B_System(self):
1134          d=self.domain.getDim()          d=self.domain.getDim()
1135          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1136          mypde.setValue(B=numpy.ones((self.N,d,self.N)))          mypde.setValue(B=numpy.ones((self.N,d,self.N)))
1137          coeff=mypde.getCoefficient("B")          coeff=mypde.getCoefficient("B")
1138          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,d,self.N),Function(self.domain),self.N,self.N))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,d,self.N),Function(self.domain),self.N,self.N))
1139            mypde.resetRightHandSideCoefficients()
1140            self.assertFalse(mypde.getCoefficient("B").isEmpty(),"B is empty after reset of right hand side coefficients")
1141      def test_setCoefficient_C_System(self):      def test_setCoefficient_C_System(self):
1142          d=self.domain.getDim()          d=self.domain.getDim()
1143          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1144          mypde.setValue(C=numpy.ones((self.N,self.N,d)))          mypde.setValue(C=numpy.ones((self.N,self.N,d)))
1145          coeff=mypde.getCoefficient("C")          coeff=mypde.getCoefficient("C")
1146          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N,d),Function(self.domain),self.N,self.N))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N,d),Function(self.domain),self.N,self.N))
1147            mypde.resetRightHandSideCoefficients()
1148            self.assertFalse(mypde.getCoefficient("C").isEmpty(),"C is empty after reset of right hand side coefficients")
1149      def test_setCoefficient_D_System(self):      def test_setCoefficient_D_System(self):
1150          d=self.domain.getDim()          d=self.domain.getDim()
1151          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1152          mypde.setValue(D=numpy.ones((self.N,self.N)))          mypde.setValue(D=numpy.ones((self.N,self.N)))
1153          coeff=mypde.getCoefficient("D")          coeff=mypde.getCoefficient("D")
1154          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),Function(self.domain),self.N,self.N))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),Function(self.domain),self.N,self.N))
1155            mypde.resetRightHandSideCoefficients()
1156            self.assertFalse(mypde.getCoefficient("D").isEmpty(),"D is empty after reset of right hand side coefficients")
1157      def test_setCoefficient_X_System(self):      def test_setCoefficient_X_System(self):
1158          d=self.domain.getDim()          d=self.domain.getDim()
1159          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1160          mypde.setValue(X=numpy.ones((self.N,d)))          mypde.setValue(X=numpy.ones((self.N,d)))
1161          coeff=mypde.getCoefficient("X")          coeff=mypde.getCoefficient("X")
1162          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,d),Function(self.domain),self.N))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,d),Function(self.domain),self.N))
1163            mypde.resetRightHandSideCoefficients()
1164            self.assertTrue(mypde.getCoefficient("X").isEmpty(),"X is not empty after reset of right hand side coefficients")
1165      def test_setCoefficient_Y_System(self):      def test_setCoefficient_Y_System(self):
1166          d=self.domain.getDim()          d=self.domain.getDim()
1167          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1168          mypde.setValue(Y=numpy.ones((self.N,)))          mypde.setValue(Y=numpy.ones((self.N,)))
1169          coeff=mypde.getCoefficient("Y")          coeff=mypde.getCoefficient("Y")
1170          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),Function(self.domain),self.N))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),Function(self.domain),self.N))
1171            mypde.resetRightHandSideCoefficients()
1172            self.assertTrue(mypde.getCoefficient("Y").isEmpty(),"Y is not empty after reset of right hand side coefficients")
1173      def test_setCoefficient_y_System(self):      def test_setCoefficient_y_System(self):
1174          d=self.domain.getDim()          d=self.domain.getDim()
1175          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1176          mypde.setValue(y=numpy.ones((self.N,)))          mypde.setValue(y=numpy.ones((self.N,)))
1177          coeff=mypde.getCoefficient("y")          coeff=mypde.getCoefficient("y")
1178          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),FunctionOnBoundary(self.domain),self.N))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),FunctionOnBoundary(self.domain),self.N))
1179            mypde.resetRightHandSideCoefficients()
1180            self.assertTrue(mypde.getCoefficient("y").isEmpty(),"y is not empty after reset of right hand side coefficients")
1181      def test_setCoefficient_d_System(self):      def test_setCoefficient_d_System(self):
1182          d=self.domain.getDim()          d=self.domain.getDim()
1183          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1184          mypde.setValue(d=numpy.ones((self.N,self.N)))          mypde.setValue(d=numpy.ones((self.N,self.N)))
1185          coeff=mypde.getCoefficient("d")          coeff=mypde.getCoefficient("d")
1186          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),FunctionOnBoundary(self.domain),self.N,self.N))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),FunctionOnBoundary(self.domain),self.N,self.N))
1187            mypde.resetRightHandSideCoefficients()
1188            self.assertFalse(mypde.getCoefficient("d").isEmpty(),"d is empty after reset of right hand side coefficients")
1189      def test_setCoefficient_d_contact_System(self):      def test_setCoefficient_d_contact_System(self):
1190          d=self.domain.getDim()          if self.domain.supportsContactElements():
1191          mypde=LinearPDE(self.domain,debug=self.DEBUG)              d=self.domain.getDim()
1192          mypde.setValue(d_contact=numpy.ones((self.N,self.N)))              mypde=LinearPDE(self.domain,debug=self.DEBUG)
1193          coeff=mypde.getCoefficient("d_contact")              mypde.setValue(d_contact=numpy.ones((self.N,self.N)))
1194          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),FunctionOnContactZero(self.domain),self.N,self.N))              coeff=mypde.getCoefficient("d_contact")
1195                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),FunctionOnContactZero(self.domain),self.N,self.N))
1196                mypde.resetRightHandSideCoefficients()
1197                self.assertFalse(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is empty after reset of right hand side coefficients")
1198      def test_setCoefficient_y_contact_System(self):      def test_setCoefficient_y_contact_System(self):
1199          d=self.domain.getDim()          if self.domain.supportsContactElements():
1200          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)              d=self.domain.getDim()
1201          mypde.setValue(y_contact=numpy.ones((self.N,)))              mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1202          coeff=mypde.getCoefficient("y_contact")              mypde.setValue(y_contact=numpy.ones((self.N,)))
1203          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),FunctionOnContactZero(self.domain),self.N))              coeff=mypde.getCoefficient("y_contact")
1204                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),FunctionOnContactZero(self.domain),self.N))
1205                mypde.resetRightHandSideCoefficients()
1206                self.assertTrue(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty after reset of right hand side coefficients")
1207      def test_setCoefficient_A_reduced_System(self):      def test_setCoefficient_A_reduced_System(self):
1208          d=self.domain.getDim()          d=self.domain.getDim()
1209          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1210          mypde.setValue(A_reduced=numpy.ones((self.N,d,self.N,d)))          mypde.setValue(A_reduced=numpy.ones((self.N,d,self.N,d)))
1211          coeff=mypde.getCoefficient("A_reduced")          coeff=mypde.getCoefficient("A_reduced")
1212          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,d,self.N,d),ReducedFunction(self.domain),self.N,self.N))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,d,self.N,d),ReducedFunction(self.domain),self.N,self.N))
1213            mypde.resetRightHandSideCoefficients()
1214            self.assertFalse(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is empty after reset of right hand side coefficients")
1215      def test_setCoefficient_B_reduced_System(self):      def test_setCoefficient_B_reduced_System(self):
1216          d=self.domain.getDim()          d=self.domain.getDim()
1217          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1218          mypde.setValue(B_reduced=numpy.ones((self.N,d,self.N)))          mypde.setValue(B_reduced=numpy.ones((self.N,d,self.N)))
1219          coeff=mypde.getCoefficient("B_reduced")          coeff=mypde.getCoefficient("B_reduced")
1220          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,d,self.N),ReducedFunction(self.domain),self.N,self.N))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,d,self.N),ReducedFunction(self.domain),self.N,self.N))
1221            mypde.resetRightHandSideCoefficients()
1222            self.assertFalse(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is empty after reset of right hand side coefficients")
1223      def test_setCoefficient_C_reduced_System(self):      def test_setCoefficient_C_reduced_System(self):
1224          d=self.domain.getDim()          d=self.domain.getDim()
1225          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1226          mypde.setValue(C_reduced=numpy.ones((self.N,self.N,d)))          mypde.setValue(C_reduced=numpy.ones((self.N,self.N,d)))
1227          coeff=mypde.getCoefficient("C_reduced")          coeff=mypde.getCoefficient("C_reduced")
1228          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N,d),ReducedFunction(self.domain),self.N,self.N))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N,d),ReducedFunction(self.domain),self.N,self.N))
1229            mypde.resetRightHandSideCoefficients()
1230            self.assertFalse(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is empty after reset of right hand side coefficients")
1231      def test_setCoefficient_D_System_reduced(self):      def test_setCoefficient_D_System_reduced(self):
1232          d=self.domain.getDim()          d=self.domain.getDim()
1233          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1234          mypde.setValue(D_reduced=numpy.ones((self.N,self.N)))          mypde.setValue(D_reduced=numpy.ones((self.N,self.N)))
1235          coeff=mypde.getCoefficient("D_reduced")          coeff=mypde.getCoefficient("D_reduced")
1236          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),ReducedFunction(self.domain),self.N,self.N))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),ReducedFunction(self.domain),self.N,self.N))
1237            mypde.resetRightHandSideCoefficients()
1238            self.assertFalse(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is empty after reset of right hand side coefficients")
1239      def test_setCoefficient_X_System_reduced(self):      def test_setCoefficient_X_System_reduced(self):
1240          d=self.domain.getDim()          d=self.domain.getDim()
1241          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1242          mypde.setValue(X_reduced=numpy.ones((self.N,d)))          mypde.setValue(X_reduced=numpy.ones((self.N,d)))
1243          coeff=mypde.getCoefficient("X_reduced")          coeff=mypde.getCoefficient("X_reduced")
1244          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,d),ReducedFunction(self.domain),self.N))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,d),ReducedFunction(self.domain),self.N))
1245            mypde.resetRightHandSideCoefficients()
1246            self.assertTrue(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty after reset of right hand side coefficients")
1247      def test_setCoefficient_Y_System_reduced(self):      def test_setCoefficient_Y_System_reduced(self):
1248          d=self.domain.getDim()          d=self.domain.getDim()
1249          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1250          mypde.setValue(Y_reduced=numpy.ones((self.N,)))          mypde.setValue(Y_reduced=numpy.ones((self.N,)))
1251          coeff=mypde.getCoefficient("Y_reduced")          coeff=mypde.getCoefficient("Y_reduced")
1252          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),ReducedFunction(self.domain),self.N))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),ReducedFunction(self.domain),self.N))
1253            mypde.resetRightHandSideCoefficients()
1254            self.assertTrue(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty after reset of right hand side coefficients")
1255      def test_setCoefficient_y_System_reduced(self):      def test_setCoefficient_y_System_reduced(self):
1256          d=self.domain.getDim()          d=self.domain.getDim()
1257          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1258          mypde.setValue(y_reduced=numpy.ones((self.N,)))          mypde.setValue(y_reduced=numpy.ones((self.N,)))
1259          coeff=mypde.getCoefficient("y_reduced")          coeff=mypde.getCoefficient("y_reduced")
1260          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),ReducedFunctionOnBoundary(self.domain),self.N))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),ReducedFunctionOnBoundary(self.domain),self.N))
1261            mypde.resetRightHandSideCoefficients()
1262            self.assertTrue(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty after reset of right hand side coefficients")
1263      def test_setCoefficient_d_reduced_System(self):      def test_setCoefficient_d_reduced_System(self):
1264          d=self.domain.getDim()          d=self.domain.getDim()
1265          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1266          mypde.setValue(d_reduced=numpy.ones((self.N,self.N)))          mypde.setValue(d_reduced=numpy.ones((self.N,self.N)))
1267          coeff=mypde.getCoefficient("d_reduced")          coeff=mypde.getCoefficient("d_reduced")
1268          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),ReducedFunctionOnBoundary(self.domain),self.N,self.N))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),ReducedFunctionOnBoundary(self.domain),self.N,self.N))
1269            mypde.resetRightHandSideCoefficients()
1270            self.assertFalse(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is empty after reset of right hand side coefficients")
1271      def test_setCoefficient_d_contact_reduced_System(self):      def test_setCoefficient_d_contact_reduced_System(self):
1272          d=self.domain.getDim()          if self.domain.supportsContactElements():
1273          mypde=LinearPDE(self.domain,debug=self.DEBUG)              d=self.domain.getDim()
1274          mypde.setValue(d_contact_reduced=numpy.ones((self.N,self.N)))              mypde=LinearPDE(self.domain,debug=self.DEBUG)
1275          coeff=mypde.getCoefficient("d_contact_reduced")              mypde.setValue(d_contact_reduced=numpy.ones((self.N,self.N)))
1276          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),ReducedFunctionOnContactZero(self.domain),self.N,self.N))              coeff=mypde.getCoefficient("d_contact_reduced")
1277                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),ReducedFunctionOnContactZero(self.domain),self.N,self.N))
1278                mypde.resetRightHandSideCoefficients()
1279                self.assertFalse(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is empty after reset of right hand side coefficients")
1280      def test_setCoefficient_y_contact_reduced_System(self):      def test_setCoefficient_y_contact_reduced_System(self):
1281          d=self.domain.getDim()          if self.domain.supportsContactElements():
1282          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)              d=self.domain.getDim()
1283          mypde.setValue(y_contact_reduced=numpy.ones((self.N,)))              mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1284          coeff=mypde.getCoefficient("y_contact_reduced")              mypde.setValue(y_contact_reduced=numpy.ones((self.N,)))
1285          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),ReducedFunctionOnContactZero(self.domain),self.N))              coeff=mypde.getCoefficient("y_contact_reduced")
1286                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),ReducedFunctionOnContactZero(self.domain),self.N))
1287                mypde.resetRightHandSideCoefficients()
1288                self.assertTrue(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty after reset of right hand side coefficients")
1289      def test_setCoefficient_r_System(self):      def test_setCoefficient_r_System(self):
1290          d=self.domain.getDim()          d=self.domain.getDim()
1291          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
1292          mypde.setValue(r=numpy.ones((self.N,)))          mypde.setValue(r=numpy.ones((self.N,)))
1293          coeff=mypde.getCoefficient("r")          coeff=mypde.getCoefficient("r")
1294          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((self.N,),Solution(self.domain),self.N))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((self.N,),Solution(self.domain),self.N))
1295            mypde.resetRightHandSideCoefficients()
1296            self.assertTrue(mypde.getCoefficient("r").isEmpty(),"r is not empty after reset of right hand side coefficients")
1297      def test_setCoefficient_q_System(self):      def test_setCoefficient_q_System(self):
1298          d=self.domain.getDim()          d=self.domain.getDim()
1299          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
1300          mypde.setValue(q=numpy.ones((self.N,)))          mypde.setValue(q=numpy.ones((self.N,)))
1301          coeff=mypde.getCoefficient("q")          coeff=mypde.getCoefficient("q")
1302          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((self.N,),Solution(self.domain),self.N))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((self.N,),Solution(self.domain),self.N))
1303            mypde.resetRightHandSideCoefficients()
1304            self.assertFalse(mypde.getCoefficient("q").isEmpty(),"q is empty after reset of right hand side coefficients")
1305      def test_setCoefficient_r_System_reducedOn(self):      def test_setCoefficient_r_System_reducedOn(self):
1306          d=self.domain.getDim()          d=self.domain.getDim()
1307          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
1308          mypde.setReducedOrderOn()          mypde.setReducedOrderOn()
1309          mypde.setValue(r=numpy.ones((self.N,)))          mypde.setValue(r=numpy.ones((self.N,)))
1310          coeff=mypde.getCoefficient("r")          coeff=mypde.getCoefficient("r")
1311          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((self.N,),ReducedSolution(self.domain),self.N))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((self.N,),ReducedSolution(self.domain),self.N))
1312            mypde.resetRightHandSideCoefficients()
1313            self.assertTrue(mypde.getCoefficient("r").isEmpty(),"r is not empty after reset of right hand side coefficients")
1314      def test_setCoefficient_q_System_reducedOn(self):      def test_setCoefficient_q_System_reducedOn(self):
1315          d=self.domain.getDim()          d=self.domain.getDim()
1316          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
1317          mypde.setReducedOrderOn()          mypde.setReducedOrderOn()
1318          mypde.setValue(q=numpy.ones((self.N,)))          mypde.setValue(q=numpy.ones((self.N,)))
1319          coeff=mypde.getCoefficient("q")          coeff=mypde.getCoefficient("q")
1320          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((self.N,),ReducedSolution(self.domain),self.N))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((self.N,),ReducedSolution(self.domain),self.N))
1321            mypde.resetRightHandSideCoefficients()
1322            self.assertFalse(mypde.getCoefficient("q").isEmpty(),"q is empty after reset of right hand side coefficients")
1323    
1324      def test_setCoefficient_A_reduced_System_using_A(self):      def test_setCoefficient_A_reduced_System_using_A(self):
1325          d=self.domain.getDim()          d=self.domain.getDim()
1326          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1327          mypde.setValue(A=Data(numpy.ones((self.N,d,self.N,d)),ReducedFunction(self.domain)))          mypde.setValue(A=Data(numpy.ones((self.N,d,self.N,d)),ReducedFunction(self.domain)))
1328          coeff=mypde.getCoefficient("A_reduced")          if self.specialInterpolationSupported():
1329          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,d,self.N,d),ReducedFunction(self.domain),self.N,self.N))              coeff_name='A'
1330                FS=Function
1331            else:
1332                coeff_name='A_reduced'
1333                FS=ReducedFunction
1334            coeff=mypde.getCoefficient(coeff_name)
1335            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,d,self.N,d),FS(self.domain),self.N,self.N))
1336            mypde.resetRightHandSideCoefficients()
1337            self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
1338    
1339      def test_setCoefficient_B_reduced_System_using_B(self):      def test_setCoefficient_B_reduced_System_using_B(self):
1340          d=self.domain.getDim()          d=self.domain.getDim()
1341          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1342          mypde.setValue(B=Data(numpy.ones((self.N,d,self.N)),ReducedFunction(self.domain)))          mypde.setValue(B=Data(numpy.ones((self.N,d,self.N)),ReducedFunction(self.domain)))
1343          coeff=mypde.getCoefficient("B_reduced")          if self.specialInterpolationSupported():
1344          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,d,self.N),ReducedFunction(self.domain),self.N,self.N))              coeff_name='B'
1345                FS=Function
1346            else:
1347                coeff_name='B_reduced'
1348                FS=ReducedFunction
1349            coeff=mypde.getCoefficient(coeff_name)
1350            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,d,self.N),FS(self.domain),self.N,self.N))
1351            mypde.resetRightHandSideCoefficients()
1352            self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
1353    
1354      def test_setCoefficient_C_reduced_System_using_C(self):      def test_setCoefficient_C_reduced_System_using_C(self):
1355          d=self.domain.getDim()          d=self.domain.getDim()
1356          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1357          mypde.setValue(C=Data(numpy.ones((self.N,self.N,d)),ReducedFunction(self.domain)))          mypde.setValue(C=Data(numpy.ones((self.N,self.N,d)),ReducedFunction(self.domain)))
1358          coeff=mypde.getCoefficient("C_reduced")          if self.specialInterpolationSupported():
1359          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N,d),ReducedFunction(self.domain),self.N,self.N))              coeff_name='C'
1360      def test_setCoefficient_D_System_reduced_using_D(self):              FS=Function
1361            else:
1362                coeff_name='C_reduced'
1363                FS=ReducedFunction
1364            coeff=mypde.getCoefficient(coeff_name)
1365            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N,d),FS(self.domain),self.N,self.N))
1366            mypde.resetRightHandSideCoefficients()
1367            self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
1368    
1369        def test_setCoefficient_D_reduced_System_using_D(self):
1370          d=self.domain.getDim()          d=self.domain.getDim()
1371          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1372          mypde.setValue(D=Data(numpy.ones((self.N,self.N)),ReducedFunction(self.domain)))          mypde.setValue(D=Data(numpy.ones((self.N,self.N)),ReducedFunction(self.domain)))
1373          coeff=mypde.getCoefficient("D_reduced")          if self.specialInterpolationSupported():
1374          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),ReducedFunction(self.domain),self.N,self.N))              coeff_name='D'
1375      def test_setCoefficient_X_System_reduced_using_X(self):              FS=Function
1376            else:
1377                coeff_name='D_reduced'
1378                FS=ReducedFunction
1379            coeff=mypde.getCoefficient(coeff_name)
1380            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),FS(self.domain),self.N,self.N))
1381            mypde.resetRightHandSideCoefficients()
1382            self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
1383    
1384        def test_setCoefficient_X_reduced_System_using_X(self):
1385          d=self.domain.getDim()          d=self.domain.getDim()
1386          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1387          mypde.setValue(X=Data(numpy.ones((self.N,d)),ReducedFunction(self.domain)))          mypde.setValue(X=Data(numpy.ones((self.N,d)),ReducedFunction(self.domain)))
1388          coeff=mypde.getCoefficient("X_reduced")          if self.specialInterpolationSupported():
1389          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,d),ReducedFunction(self.domain),self.N))              coeff_name='X'
1390      def test_setCoefficient_Y_System_reduced_using_Y(self):              FS=Function
1391            else:
1392                coeff_name='X_reduced'
1393                FS=ReducedFunction
1394            coeff=mypde.getCoefficient(coeff_name)
1395            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,d),FS(self.domain),self.N))
1396            mypde.resetRightHandSideCoefficients()
1397            self.assertTrue(mypde.getCoefficient(coeff_name).isEmpty(),"%s is not empty after reset of right hand side coefficients"%coeff_name)
1398    
1399        def test_setCoefficient_Y_reduced_System_using_Y(self):
1400          d=self.domain.getDim()          d=self.domain.getDim()
1401          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1402          mypde.setValue(Y=Data(numpy.ones((self.N,)),ReducedFunction(self.domain)))          mypde.setValue(Y=Data(numpy.ones((self.N,)),ReducedFunction(self.domain)))
1403          coeff=mypde.getCoefficient("Y_reduced")          if self.specialInterpolationSupported():
1404          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),ReducedFunction(self.domain),self.N))              coeff_name='Y'
1405                FS=Function
1406            else:
1407                coeff_name='Y_reduced'
1408                FS=ReducedFunction
1409            coeff=mypde.getCoefficient(coeff_name)
1410            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),FS(self.domain),self.N))
1411            mypde.resetRightHandSideCoefficients()
1412            self.assertTrue(mypde.getCoefficient(coeff_name).isEmpty(),"%s is not empty after reset of right hand side coefficients"%coeff_name)
1413    
1414      def test_setCoefficient_y_reduced_System_using_y(self):      def test_setCoefficient_y_reduced_System_using_y(self):
1415          d=self.domain.getDim()          d=self.domain.getDim()
1416          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1417          mypde.setValue(y=Data(numpy.ones((self.N,)),ReducedFunctionOnBoundary(self.domain)))          mypde.setValue(y=Data(numpy.ones((self.N,)),ReducedFunctionOnBoundary(self.domain)))
1418          coeff=mypde.getCoefficient("y_reduced")          if self.specialInterpolationSupported():
1419          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),ReducedFunctionOnBoundary(self.domain),self.N))              coeff_name='y'
1420                FS=FunctionOnBoundary
1421            else:
1422                coeff_name='y_reduced'
1423                FS=ReducedFunctionOnBoundary
1424            coeff=mypde.getCoefficient(coeff_name)
1425            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),FS(self.domain),self.N))
1426            mypde.resetRightHandSideCoefficients()
1427            self.assertTrue(mypde.getCoefficient(coeff_name).isEmpty(),"%s is not empty after reset of right hand side coefficients"%coeff_name)
1428    
1429      def test_setCoefficient_d_reduced_System_using_d(self):      def test_setCoefficient_d_reduced_System_using_d(self):
1430          d=self.domain.getDim()          d=self.domain.getDim()
1431          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1432          mypde.setValue(d=Data(numpy.ones((self.N,self.N)),ReducedFunctionOnBoundary(self.domain)))          mypde.setValue(d=Data(numpy.ones((self.N,self.N)),ReducedFunctionOnBoundary(self.domain)))
1433          coeff=mypde.getCoefficient("d_reduced")          if self.specialInterpolationSupported():
1434          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),ReducedFunctionOnBoundary(self.domain),self.N,self.N))              coeff_name='d'
1435                FS=FunctionOnBoundary
1436            else:
1437                coeff_name='d_reduced'
1438                FS=ReducedFunctionOnBoundary
1439            coeff=mypde.getCoefficient(coeff_name)
1440            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,self.N),FS(self.domain),self.N,self.N))
1441            mypde.resetRightHandSideCoefficients()
1442            self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
1443    
1444      def test_setCoefficient_d_contact_reduced_System_using_d_contact(self):      def test_setCoefficient_d_contact_reduced_System_using_d_contact(self):
1445          d=self.domain.getDim()          if self.domain.supportsContactElements():
1446          mypde=LinearPDE(self.domain,debug=self.DEBUG)              d=self.domain.getDim()
1447          mypde.setValue(d_contact=Data(numpy.ones((self.N,self.N)),ReducedFunctionOnContactZero(self.domain)))              mypde=LinearPDE(self.domain,debug=self.DEBUG)
1448          coeff=mypde.getCoefficient("d_contact_reduced")              mypde.setValue(d_contact=Data(numpy.ones((self.N,self.N)),ReducedFunctionOnContactZero(self.domain)))
1449          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),ReducedFunctionOnContactZero(self.domain),self.N,self.N))              if self.specialInterpolationSupported():
1450                    coeff_name='d_contact'
1451                    FS=FunctionOnContactZero
1452                else:
1453                    coeff_name='d_contact_reduced'
1454                    FS=ReducedFunctionOnContactZero
1455                coeff=mypde.getCoefficient(coeff_name)
1456                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,self.N),FS(self.domain),self.N,self.N))
1457                mypde.resetRightHandSideCoefficients()
1458                self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
1459    
1460      def test_setCoefficient_y_contact_reduced_System_using_y_contact(self):      def test_setCoefficient_y_contact_reduced_System_using_y_contact(self):
1461          d=self.domain.getDim()          if self.domain.supportsContactElements():
1462          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)              d=self.domain.getDim()
1463          mypde.setValue(y_contact=Data(numpy.ones((self.N,)),ReducedFunctionOnContactZero(self.domain)))              mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1464          coeff=mypde.getCoefficient("y_contact_reduced")              mypde.setValue(y_contact=Data(numpy.ones((self.N,)),ReducedFunctionOnContactZero(self.domain)))
1465          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),ReducedFunctionOnContactZero(self.domain),self.N))              if self.specialInterpolationSupported():
1466                    coeff_name='y_contact'
1467                    FS=FunctionOnContactZero
1468                else:
1469                    coeff_name='y_contact_reduced'
1470                    FS=ReducedFunctionOnContactZero
1471                coeff=mypde.getCoefficient(coeff_name)
1472                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),FS(self.domain),self.N))
1473                mypde.resetRightHandSideCoefficients()
1474                self.assertTrue(mypde.getCoefficient(coeff_name).isEmpty(),"%s is not empty after reset of right hand side coefficients"%coeff_name)
1475    
1476      def test_resetCoefficient_HomogeneousConstraint(self):      def test_resetCoefficient_HomogeneousConstraint(self):
1477          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1478          x=self.domain.getX()          x=self.domain.getX()
# Line 1393  class Test_LinearPDE_noLumping(Test_line Line 1480  class Test_LinearPDE_noLumping(Test_line
1480          u1=mypde.getSolution()          u1=mypde.getSolution()
1481          mypde.setValue(Y=2.)          mypde.setValue(Y=2.)
1482          u2=mypde.getSolution()          u2=mypde.getSolution()
1483          self.failUnless(self.check(u2,2*u1),'solution is wrong.')          self.assertTrue(self.check(u2,2*u1),'solution is wrong.')
1484    
1485      def test_resetCoefficient_InHomogeneousConstraint(self):      def test_resetCoefficient_InHomogeneousConstraint(self):
1486          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1487          mypde.setSymmetryOn()          mypde.setSymmetryOn()
1488      mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1489          x=self.domain.getX()          x=self.domain.getX()
1490          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.,r=1,q=whereZero(x[0]))          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.,r=1,q=whereZero(x[0]))
1491          u1=mypde.getSolution()          u1=mypde.getSolution()
1492          mypde.setValue(Y=2.,D=2)          mypde.setValue(Y=2.,D=2)
1493          u2=mypde.getSolution()          u2=mypde.getSolution()
1494          self.failUnless(self.check(u2,u1),'first solution is wrong.')          self.assertTrue(self.check(u2,u1),'first solution is wrong.')
1495          u2=mypde.getSolution()          u2=mypde.getSolution()
1496          self.failUnless(self.check(u2,u1),'first solution is wrong.')          self.assertTrue(self.check(u2,u1),'first solution is wrong.')
1497          mypde.setValue(r=2,Y=4.)          mypde.setValue(r=2,Y=4.)
1498          u2=mypde.getSolution()          u2=mypde.getSolution()
1499          self.failUnless(self.check(u2,2*u1),'second solution is wrong.')          self.assertTrue(self.check(u2,2*u1),'second solution is wrong.')
1500    
1501      def test_Status(self):      def test_Status(self):
1502          DIM=self.domain.getDim()          DIM=self.domain.getDim()
# Line 1422  class Test_LinearPDE_noLumping(Test_line Line 1509  class Test_LinearPDE_noLumping(Test_line
1509          u1_ref=x1[0]*(1.-x1[0])          u1_ref=x1[0]*(1.-x1[0])
1510          u1=mypde.getSolution()          u1=mypde.getSolution()
1511          error1=Lsup(u1-u1_ref)/Lsup(u1_ref)          error1=Lsup(u1-u1_ref)/Lsup(u1_ref)
1512          self.failUnless(mypde.getDomainStatus() == mypde.getSystemStatus(), "status of first pde does not match domain status.")          self.assertTrue(mypde.getDomainStatus() == mypde.getSystemStatus(), "status of first pde does not match domain status.")
1513            try:
1514          self.domain.setX(x*5)              self.domain.setX(x*5)
1515            except:
1516          self.failUnless(mypde.getDomainStatus() != mypde.getSystemStatus(), "status of first pde matches updated domain status.")              # setX not supported
1517                return
1518            self.assertTrue(mypde.getDomainStatus() != mypde.getSystemStatus(), "status of first pde matches updated domain status.")
1519          x2=self.domain.getX()          x2=self.domain.getX()
1520          u2_ref=x2[0]*(5.-x2[0])          u2_ref=x2[0]*(5.-x2[0])
1521          u2=mypde.getSolution()          u2=mypde.getSolution()
1522          error2=Lsup(u2-u2_ref)/Lsup(u2_ref)          error2=Lsup(u2-u2_ref)/Lsup(u2_ref)
1523          self.failUnless(error2 <= max(error1,self.RES_TOL)*10., "solution of second PDE wrong.")          self.assertTrue(error2 <= max(error1,self.RES_TOL)*10., "solution of second PDE wrong.")
1524          self.failUnless(mypde.getDomainStatus() == mypde.getSystemStatus(), "status of second pde does not match domain status.")          self.assertTrue(mypde.getDomainStatus() == mypde.getSystemStatus(), "status of second pde does not match domain status.")
1525    
1526      def test_symmetryCheckTrue_System(self):      def test_symmetryCheckTrue_System(self):
1527          d=self.domain.getDim()          d=self.domain.getDim()
# Line 1443  class Test_LinearPDE_noLumping(Test_line Line 1532  class Test_LinearPDE_noLumping(Test_line
1532          D=3*numpy.ones((self.N,self.N))          D=3*numpy.ones((self.N,self.N))
1533          d=4*numpy.ones((self.N,self.N))          d=4*numpy.ones((self.N,self.N))
1534          d_contact=5*numpy.ones((self.N,self.N))          d_contact=5*numpy.ones((self.N,self.N))
1535          mypde.setValue(A=A,B=B,C=C,D=D,d=d,d_contact=d_contact,A_reduced=-A,B_reduced=-B,C_reduced=-C,D_reduced=-D,d_reduced=-d,d_contact_reduced=-d_contact)          pars={"A":A, "B":B, "C":C, "D":D, "d":d, "A_reduced":-A, "B_reduced":-B, "C_reduced":-C, "D_reduced":-D, "d_reduced":-d}
1536          self.failUnless(mypde.checkSymmetry(verbose=False),"symmetry detected")          if self.domain.supportsContactElements():
1537                    pars["d_contact"]=d_contact
1538                    pars["d_contact_reduced"]=-d_contact
1539            mypde.setValue(**pars)
1540            self.assertTrue(mypde.checkSymmetry(verbose=False),"symmetry detected")
1541    
1542      def test_symmetryCheckFalse_A_System(self):      def test_symmetryCheckFalse_A_System(self):
1543          d=self.domain.getDim()          d=self.domain.getDim()
# Line 1452  class Test_LinearPDE_noLumping(Test_line Line 1545  class Test_LinearPDE_noLumping(Test_line
1545          A=numpy.ones((self.N,d,self.N,d))          A=numpy.ones((self.N,d,self.N,d))
1546          A[1,1,1,0]=0.          A[1,1,1,0]=0.
1547          mypde.setValue(A=A)          mypde.setValue(A=A)
1548          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1549      def test_symmetryCheckFalse_BC_System(self):      def test_symmetryCheckFalse_BC_System(self):
1550          d=self.domain.getDim()          d=self.domain.getDim()
1551          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
# Line 1460  class Test_LinearPDE_noLumping(Test_line Line 1553  class Test_LinearPDE_noLumping(Test_line
1553          B=2*numpy.ones((self.N,d,self.N))          B=2*numpy.ones((self.N,d,self.N))
1554          B[0,0,1]=1.          B[0,0,1]=1.
1555          mypde.setValue(B=B,C=C)          mypde.setValue(B=B,C=C)
1556          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1557    
1558      def test_symmetryCheckFalse_D_System(self):      def test_symmetryCheckFalse_D_System(self):
1559          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1560          D=3*numpy.ones((self.N,self.N))          D=3*numpy.ones((self.N,self.N))
1561          D[0,1]=0.          D[0,1]=0.
1562          mypde.setValue(D=D)          mypde.setValue(D=D)
1563          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1564    
1565      def test_symmetryCheckFalse_d_System(self):      def test_symmetryCheckFalse_d_System(self):
1566          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1567          d=4*numpy.ones((self.N,self.N))          d=4*numpy.ones((self.N,self.N))
1568          d[0,1]=0.          d[0,1]=0.
1569          mypde.setValue(d=d)          mypde.setValue(d=d)
1570          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1571    
1572      def test_symmetryCheckFalse_d_contact_System(self):      def test_symmetryCheckFalse_d_contact_System(self):
1573          mypde=LinearPDE(self.domain,debug=self.DEBUG)          if self.domain.supportsContactElements():
1574          d_contact=5*numpy.ones((self.N,self.N))              mypde=LinearPDE(self.domain,debug=self.DEBUG)
1575          d_contact[0,1]=0.              d_contact=5*numpy.ones((self.N,self.N))
1576          mypde.setValue(d_contact=d_contact)              d_contact[0,1]=0.
1577          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")              mypde.setValue(d_contact=d_contact)
1578                self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1579    
1580      def test_symmetryCheckFalse_A_reduced_System(self):      def test_symmetryCheckFalse_A_reduced_System(self):
1581          d=self.domain.getDim()          d=self.domain.getDim()
# Line 1489  class Test_LinearPDE_noLumping(Test_line Line 1583  class Test_LinearPDE_noLumping(Test_line
1583          A=numpy.ones((self.N,d,self.N,d))          A=numpy.ones((self.N,d,self.N,d))
1584          A[1,1,1,0]=0.          A[1,1,1,0]=0.
1585          mypde.setValue(A_reduced=A)          mypde.setValue(A_reduced=A)
1586          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1587      def test_symmetryCheckFalse_BC_reduced_System(self):      def test_symmetryCheckFalse_BC_reduced_System(self):
1588          d=self.domain.getDim()          d=self.domain.getDim()
1589          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
# Line 1497  class Test_LinearPDE_noLumping(Test_line Line 1591  class Test_LinearPDE_noLumping(Test_line
1591          B=2*numpy.ones((self.N,d,self.N))          B=2*numpy.ones((self.N,d,self.N))
1592          B[0,0,1]=1.          B[0,0,1]=1.
1593          mypde.setValue(B_reduced=B,C_reduced=C)          mypde.setValue(B_reduced=B,C_reduced=C)
1594          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1595    
1596      def test_symmetryCheckFalse_D_reduced_System(self):      def test_symmetryCheckFalse_D_reduced_System(self):
1597          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1598          D=3*numpy.ones((self.N,self.N))          D=3*numpy.ones((self.N,self.N))
1599          D[0,1]=0.          D[0,1]=0.
1600          mypde.setValue(D_reduced=D)          mypde.setValue(D_reduced=D)
1601          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1602    
1603      def test_symmetryCheckFalse_d_reduced_System(self):      def test_symmetryCheckFalse_d_reduced_System(self):
1604          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1605          d=4*numpy.ones((self.N,self.N))          d=4*numpy.ones((self.N,self.N))
1606          d[0,1]=0.          d[0,1]=0.
1607          mypde.setValue(d_reduced=d)          mypde.setValue(d_reduced=d)
1608          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1609    
1610      def test_symmetryCheckFalse_d_contact_reduced_System(self):      def test_symmetryCheckFalse_d_contact_reduced_System(self):
1611          mypde=LinearPDE(self.domain,debug=self.DEBUG)          if self.domain.supportsContactElements():
1612          d_contact=5*numpy.ones((self.N,self.N))              mypde=LinearPDE(self.domain,debug=self.DEBUG)
1613          d_contact[0,1]=0.              d_contact=5*numpy.ones((self.N,self.N))
1614          mypde.setValue(d_contact_reduced=d_contact)              d_contact[0,1]=0.
1615          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")              mypde.setValue(d_contact_reduced=d_contact)
1616                self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1617    
1618      def test_symmetryCheckTrue_Scalar(self):      def test_symmetryCheckTrue_Scalar(self):
1619          d=self.domain.getDim()          d=self.domain.getDim()
# Line 1529  class Test_LinearPDE_noLumping(Test_line Line 1624  class Test_LinearPDE_noLumping(Test_line
1624          D=3          D=3
1625          d=4          d=4
1626          d_contact=5          d_contact=5
1627          mypde.setValue(A=A,B=B,C=C,D=D,d=d,d_contact=d_contact,A_reduced=-A,B_reduced=-B,C_reduced=-C,D_reduced=-D,d_reduced=-d,d_contact_reduced=-d_contact)          pars={"A":A, "B":B, "C":C, "D":D, "d":d, "A_reduced":-A, "B_reduced":-B, "C_reduced":-C, "D_reduced":-D, "d_reduced":-d}
1628          self.failUnless(mypde.checkSymmetry(verbose=False),"symmetry detected")          if self.domain.supportsContactElements():
1629                    pars["d_contact"]=d_contact
1630                    pars["d_contact_reduced"]=-d_contact
1631            mypde.setValue(**pars)
1632            self.assertTrue(mypde.checkSymmetry(verbose=False),"symmetry detected")
1633    
1634      def test_symmetryCheckFalse_A_Scalar(self):      def test_symmetryCheckFalse_A_Scalar(self):
1635          d=self.domain.getDim()          d=self.domain.getDim()
# Line 1538  class Test_LinearPDE_noLumping(Test_line Line 1637  class Test_LinearPDE_noLumping(Test_line
1637          A=numpy.ones((d,d))          A=numpy.ones((d,d))
1638          A[1,0]=0.          A[1,0]=0.
1639          mypde.setValue(A=A)          mypde.setValue(A=A)
1640          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1641      def test_symmetryCheckFalse_BC_Scalar(self):      def test_symmetryCheckFalse_BC_Scalar(self):
1642          d=self.domain.getDim()          d=self.domain.getDim()
1643          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
# Line 1546  class Test_LinearPDE_noLumping(Test_line Line 1645  class Test_LinearPDE_noLumping(Test_line
1645          B=2*numpy.ones((d,))          B=2*numpy.ones((d,))
1646          B[0]=1.          B[0]=1.
1647          mypde.setValue(B=B,C=C)          mypde.setValue(B=B,C=C)
1648          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1649      def test_symmetryCheckFalse_A_reduced_Scalar(self):      def test_symmetryCheckFalse_A_reduced_Scalar(self):
1650          d=self.domain.getDim()          d=self.domain.getDim()
1651          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1652          A=numpy.ones((d,d))          A=numpy.ones((d,d))
1653          A[1,0]=0.          A[1,0]=0.
1654          mypde.setValue(A_reduced=A)          mypde.setValue(A_reduced=A)
1655          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1656      def test_symmetryCheckFalse_BC_reduced_Scalar(self):      def test_symmetryCheckFalse_BC_reduced_Scalar(self):
1657          d=self.domain.getDim()          d=self.domain.getDim()
1658          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
# Line 1561  class Test_LinearPDE_noLumping(Test_line Line 1660  class Test_LinearPDE_noLumping(Test_line
1660          B=2*numpy.ones((d,))          B=2*numpy.ones((d,))
1661          B[0]=1.          B[0]=1.
1662          mypde.setValue(B_reduced=B,C_reduced=C)          mypde.setValue(B_reduced=B,C_reduced=C)
1663          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1664      #      #
1665      #   solver checks (single PDE)      #   solver checks (single PDE)
1666      #      #
1667      def test_symmetryOnIterative(self):      def test_symmetryOnIterative(self):
1668          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1669          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1670      mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1671          u=mypde.getSolution()          u=mypde.getSolution()
1672          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1673      def test_symmetryOnDirect(self):      def test_symmetryOnDirect(self):
1674          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1675          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1676          mypde.getSolverOptions().setSolverMethod(SolverOptions.DIRECT)          mypde.getSolverOptions().setSolverMethod(SolverOptions.DIRECT)
1677      mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1678          u=mypde.getSolution()          u=mypde.getSolution()
1679          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1680      def test_PCG_JACOBI(self):      def test_PCG_JACOBI(self):
1681          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1682          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1683          mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)          mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
1684      mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)          mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
1685      mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1686            u=mypde.getSolution()
1687            self.assertTrue(self.check(u,1.),'solution is wrong.')
1688        def test_PCG_GAUSS_SEIDEL(self):
1689            mypde=LinearPDE(self.domain,debug=self.DEBUG)
1690            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1691            mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
1692            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
1693            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1694          u=mypde.getSolution()          u=mypde.getSolution()
1695          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1696        def test_PCG_AMG(self):
1697            if self.order!=2:
1698                if getEscriptParamInt('DISABLE_AMG', 0):
1699                    print("AMG test disabled on MPI build")
1700                    return
1701                mypde=LinearPDE(self.domain,debug=self.DEBUG)
1702                mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1703                mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
1704                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
1705                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1706                u=mypde.getSolution()
1707                self.assertTrue(self.check(u,1.),'solution is wrong.')
1708      def test_PCG_ILU0(self):      def test_PCG_ILU0(self):
1709          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1710          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1711      mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)          mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
1712      mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
1713      mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1714          u=mypde.getSolution()          u=mypde.getSolution()
1715          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1716      def test_PCG_RILU(self):      def test_PCG_RILU(self):
1717          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1718          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1719      mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)          mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
1720      mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)          mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)
1721      mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1722          u=mypde.getSolution()          u=mypde.getSolution()
1723          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1724      def test_PCG_REC_ILU(self):      def test_PCG_REC_ILU(self):
1725          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1726          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1727      mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)          mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
1728      mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)          mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)
1729      mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1730          u=mypde.getSolution()          u=mypde.getSolution()
1731          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1732      def test_DIRECT(self):      def test_DIRECT(self):
1733          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1734          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1735          mypde.getSolverOptions().setSolverMethod(SolverOptions.DIRECT)          mypde.getSolverOptions().setSolverMethod(SolverOptions.DIRECT)
1736      mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1737          u=mypde.getSolution()          u=mypde.getSolution()
1738          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1739      def test_BICGSTAB_JACOBI(self):      def test_BICGSTAB_JACOBI(self):
1740          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1741      mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)          mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)
1742      mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)          mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
1743            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1744            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1745            u=mypde.getSolution()
1746            self.assertTrue(self.check(u,1.),'solution is wrong.')
1747        def test_BICGSTAB_GAUSS_SEIDEL(self):
1748            mypde=LinearPDE(self.domain,debug=self.DEBUG)
1749            mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)
1750            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
1751          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1752          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1753          u=mypde.getSolution()          u=mypde.getSolution()
1754          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1755        def test_BICGSTAB_AMG(self):
1756            if self.order!=2:
1757                if getEscriptParamInt('DISABLE_AMG', 0):
1758                    print("AMG test disabled on MPI build")
1759                    return    
1760                mypde=LinearPDE(self.domain,debug=self.DEBUG)
1761                mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)
1762                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
1763                mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1764                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1765                u=mypde.getSolution()
1766                self.assertTrue(self.check(u,1.),'solution is wrong.')
1767      def test_BICGSTAB_ILU0(self):      def test_BICGSTAB_ILU0(self):
1768          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1769          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1770      mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)          mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)
1771      mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
1772          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1773          u=mypde.getSolution()          u=mypde.getSolution()
1774          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1775      def test_BICGSTAB_RILU(self):      def test_BICGSTAB_RILU(self):
1776          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1777          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1778      mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)          mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)
1779      mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)          mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)
1780          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1781          u=mypde.getSolution()          u=mypde.getSolution()
1782          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1783      def test_BICGSTAB_REC_ILU(self):      def test_BICGSTAB_REC_ILU(self):
1784          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1785          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1786      mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)          mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)
1787      mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)          mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)
1788          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1789          u=mypde.getSolution()          u=mypde.getSolution()
1790          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1791      def test_MINRES_JACOBI(self):      def test_MINRES_JACOBI(self):
1792          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1793      mypde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)          mypde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
1794      mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)          mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
1795          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1796          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1797          u=mypde.getSolution()          u=mypde.getSolution()
1798          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1799        def test_MINRES_GAUSS_SEIDEL(self):
1800            mypde=LinearPDE(self.domain,debug=self.DEBUG)
1801            mypde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
1802            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
1803            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1804            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1805            u=mypde.getSolution()
1806            self.assertTrue(self.check(u,1.),'solution is wrong.')
1807        def test_MINRES_AMG(self):
1808            if self.order!=2:
1809                if getEscriptParamInt('DISABLE_AMG',0):
1810                    print("AMG test disabled on MPI build")
1811                    return                
1812                mypde=LinearPDE(self.domain,debug=self.DEBUG)
1813                mypde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
1814                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
1815                mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1816                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1817                u=mypde.getSolution()
1818                self.assertTrue(self.check(u,1.),'solution is wrong.')
1819      def test_MINRES_ILU0(self):      def test_MINRES_ILU0(self):
1820          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1821          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1822      mypde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)          mypde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
1823      mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
1824          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1825          u=mypde.getSolution()          u=mypde.getSolution()
1826          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1827      def test_MINRES_RILU(self):      def test_MINRES_RILU(self):
1828          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1829          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1830      mypde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)          mypde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
1831      mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)          mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)
1832          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1833          u=mypde.getSolution()          u=mypde.getSolution()
1834          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1835      def test_MINRES_REC_ILU(self):      def test_MINRES_REC_ILU(self):
1836          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1837          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1838      mypde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)          mypde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
1839      mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)          mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)
1840          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1841          u=mypde.getSolution()          u=mypde.getSolution()
1842          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1843      def test_TFQMR_JACOBI(self):      def test_TFQMR_JACOBI(self):
1844          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1845      mypde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)          mypde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
1846      mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)          mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
1847          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1848          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1849          u=mypde.getSolution()          u=mypde.getSolution()
1850          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1851        def test_TFQMR_GAUSS_SEIDEL(self):
1852            mypde=LinearPDE(self.domain,debug=self.DEBUG)
1853            mypde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
1854            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
1855            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1856            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1857            u=mypde.getSolution()
1858            self.assertTrue(self.check(u,1.),'solution is wrong.')
1859        def test_TFQMR_AMG(self):
1860            if self.order!=2:
1861                if getEscriptParamInt('DISABLE_AMG', 0):
1862                    print("AMG test disabled on MPI build")
1863                    return
1864                mypde=LinearPDE(self.domain,debug=self.DEBUG)
1865                mypde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
1866                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
1867                mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1868                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1869                u=mypde.getSolution()
1870                self.assertTrue(self.check(u,1.),'solution is wrong.')
1871      def test_TFQMR_ILU0(self):      def test_TFQMR_ILU0(self):
1872          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1873          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1874      mypde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)          mypde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
1875      mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
1876          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1877          u=mypde.getSolution()          u=mypde.getSolution()
1878          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1879      def test_TFQMR_RILU(self):      def test_TFQMR_RILU(self):
1880          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1881          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1882      mypde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)          mypde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
1883      mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)          mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)
1884          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1885          u=mypde.getSolution()          u=mypde.getSolution()
1886          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1887      def test_TFQMR_REC_ILU(self):      def test_TFQMR_REC_ILU(self):
1888          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1889          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1890      mypde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)          mypde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
1891      mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)          mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)
1892          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1893          u=mypde.getSolution()          u=mypde.getSolution()
1894          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1895      def test_PRES20_JACOBI(self):      def test_PRES20_JACOBI(self):
1896          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1897          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1898      mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)          mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)
1899      mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)          mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
1900          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1901          u=mypde.getSolution()          u=mypde.getSolution()
1902          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1903        def test_PRES20_GAUSS_SEIDEL(self):
1904            mypde=LinearPDE(self.domain,debug=self.DEBUG)
1905            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1906            mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)
1907            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
1908            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1909            u=mypde.getSolution()
1910            self.assertTrue(self.check(u,1.),'solution is wrong.')
1911        def test_PRES20_AMG(self):
1912            if self.order!=2:
1913                if getEscriptParamInt('DISABLE_AMG', 0):
1914                    print("AMG test disabled on MPI build")
1915                    return
1916                mypde=LinearPDE(self.domain,debug=self.DEBUG)
1917                mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1918                mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)
1919                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
1920                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1921                u=mypde.getSolution()
1922                self.assertTrue(self.check(u,1.),'solution is wrong.')
1923      def test_PRES20_ILU0(self):      def test_PRES20_ILU0(self):
1924          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1925          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1926      mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)          mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)
1927      mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
1928          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1929          u=mypde.getSolution()          u=mypde.getSolution()
1930          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1931      def test_PRES20_RILU(self):      def test_PRES20_RILU(self):
1932          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1933          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1934      mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)          mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)
1935      mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)          mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)
1936          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1937          u=mypde.getSolution()          u=mypde.getSolution()
1938          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1939      def test_PRES20_REC_ILU(self):      def test_PRES20_REC_ILU(self):
1940          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1941          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1942      mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)          mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)
1943      mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)          mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)
1944          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1945          u=mypde.getSolution()          u=mypde.getSolution()
1946          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1947      def test_GMRESnoRestart_JACOBI(self):      def test_GMRESnoRestart_JACOBI(self):
1948          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1949          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1950      mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
1951      mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)          mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
1952      mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1953      mypde.getSolverOptions().setTruncation(50)          mypde.getSolverOptions().setTruncation(50)
1954          u=mypde.getSolution()          u=mypde.getSolution()
1955          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1956        def test_GMRESnoRestart_GAUSS_SEIDEL(self):
1957            mypde=LinearPDE(self.domain,debug=self.DEBUG)
1958            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1959            mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
1960            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
1961            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1962            mypde.getSolverOptions().setTruncation(50)
1963            u=mypde.getSolution()
1964            self.assertTrue(self.check(u,1.),'solution is wrong.')
1965        def test_GMRESnoRestart_AMG(self):
1966            if self.order!=2:
1967                if getEscriptParamInt('DISABLE_AMG', 0):
1968                    print("AMG test disabled on MPI build")
1969                    return    
1970                mypde=LinearPDE(self.domain,debug=self.DEBUG)
1971                mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1972                mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
1973                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
1974                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1975                mypde.getSolverOptions().setTruncation(50)
1976                u=mypde.getSolution()
1977                self.assertTrue(self.check(u,1.),'solution is wrong.')
1978      def test_GMRESnoRestart_ILU0(self):      def test_GMRESnoRestart_ILU0(self):
1979          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1980          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1981      mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
1982      mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
1983      mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1984      mypde.getSolverOptions().setTruncation(50)                                  mypde.getSolverOptions().setTruncation(50)                        
1985          u=mypde.getSolution()          u=mypde.getSolution()
1986          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1987      def test_GMRESnoRestart_RILU(self):      def test_GMRESnoRestart_RILU(self):
1988          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1989          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1990      mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
1991      mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)          mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)
1992      mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1993      mypde.getSolverOptions().setTruncation(50)          mypde.getSolverOptions().setTruncation(50)
1994          u=mypde.getSolution()          u=mypde.getSolution()
1995          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
1996      def test_GMRESnoRestart_REC_ILU(self):      def test_GMRESnoRestart_REC_ILU(self):
1997          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1998          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1999      mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2000      mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)          mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)
2001      mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2002      mypde.getSolverOptions().setTruncation(50)          mypde.getSolverOptions().setTruncation(50)
2003          u=mypde.getSolution()          u=mypde.getSolution()
2004          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
2005      def test_GMRES_JACOBI(self):      def test_GMRES_JACOBI(self):
2006          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2007          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2008      mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2009      mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)          mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
2010          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2011          u=mypde.getSolution()          u=mypde.getSolution()
2012          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
2013        def test_GMRES_GAUSS_SEIDEL(self):
2014            mypde=LinearPDE(self.domain,debug=self.DEBUG)
2015            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2016            mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2017            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
2018            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2019            u=mypde.getSolution()
2020            self.assertTrue(self.check(u,1.),'solution is wrong.')
2021        def test_GMRES_AMG(self):
2022            if self.order!=2:
2023                if getEscriptParamInt('DISABLE_AMG', 0):
2024                    print("AMG test disabled on MPI build")
2025                    return    
2026                mypde=LinearPDE(self.domain,debug=self.DEBUG)
2027                mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2028                mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2029                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
2030                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2031                u=mypde.getSolution()
2032                self.assertTrue(self.check(u,1.),'solution is wrong.')        
2033      def test_GMRES_ILU0(self):      def test_GMRES_ILU0(self):
2034          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2035          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2036      mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2037      mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
2038          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2039          u=mypde.getSolution()          u=mypde.getSolution()
2040          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
2041      def test_GMRES_RILU(self):      def test_GMRES_RILU(self):
2042          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2043          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2044      mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2045      mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)          mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)
2046          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2047          u=mypde.getSolution()          u=mypde.getSolution()
2048          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
2049      def test_GMRES_REC_ILU(self):      def test_GMRES_REC_ILU(self):
2050          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2051          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2052      mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2053      mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)          mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)
2054          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2055          u=mypde.getSolution()          u=mypde.getSolution()
2056          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
2057      def test_GMRES_truncation_restart_JACOBI(self):      def test_GMRES_truncation_restart_JACOBI(self):
2058          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2059          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2060      mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2061      mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)          mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
2062      mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2063      mypde.getSolverOptions().setTruncation(10)          mypde.getSolverOptions().setTruncation(10)
2064      mypde.getSolverOptions().setRestart(20)          mypde.getSolverOptions().setRestart(20)
2065          u=mypde.getSolution()          u=mypde.getSolution()
2066          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
2067        def test_GMRES_truncation_restart_GAUSS_SEIDEL(self):
2068            mypde=LinearPDE(self.domain,debug=self.DEBUG)
2069            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2070            mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2071            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
2072            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2073            mypde.getSolverOptions().setTruncation(10)
2074            mypde.getSolverOptions().setRestart(20)
2075            u=mypde.getSolution()
2076            self.assertTrue(self.check(u,1.),'solution is wrong.')
2077        def test_GMRES_truncation_restart_AMG(self):
2078            if self.order!=2:
2079                if getEscriptParamInt('DISABLE_AMG', 0):
2080                    print("AMG test disabled on MPI build")
2081                    return    
2082                mypde=LinearPDE(self.domain,debug=self.DEBUG)
2083                mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2084                mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2085                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
2086                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2087                mypde.getSolverOptions().setTruncation(10)
2088                mypde.getSolverOptions().setRestart(20)
2089                u=mypde.getSolution()
2090                self.assertTrue(self.check(u,1.),'solution is wrong.')
2091      def test_GMRES_truncation_restart_ILU0(self):      def test_GMRES_truncation_restart_ILU0(self):
2092          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2093          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2094      mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2095      mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
2096      mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2097      mypde.getSolverOptions().setTruncation(10)          mypde.getSolverOptions().setTruncation(10)
2098      mypde.getSolverOptions().setRestart(20)          mypde.getSolverOptions().setRestart(20)
2099          u=mypde.getSolution()          u=mypde.getSolution()
2100          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
2101      def test_GMRES_truncation_restart_RILU(self):      def test_GMRES_truncation_restart_RILU(self):
2102          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2103          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2104      mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2105      mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)          mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)
2106      mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2107      mypde.getSolverOptions().setTruncation(10)          mypde.getSolverOptions().setTruncation(10)
2108      mypde.getSolverOptions().setRestart(20)          mypde.getSolverOptions().setRestart(20)
2109          u=mypde.getSolution()          u=mypde.getSolution()
2110          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
2111      def test_GMRES_truncation_restart_REC_ILU(self):      def test_GMRES_truncation_restart_REC_ILU(self):
2112          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2113          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2114      mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2115      mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)          mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)
2116      mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2117      mypde.getSolverOptions().setTruncation(10)          mypde.getSolverOptions().setTruncation(10)
2118      mypde.getSolverOptions().setRestart(20)          mypde.getSolverOptions().setRestart(20)
2119          u=mypde.getSolution()          u=mypde.getSolution()
2120          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
2121      #      #
2122      #   solver checks (PDE system)      #   solver checks (PDE system)
2123      #      #
# Line 1868  class Test_LinearPDE_noLumping(Test_line Line 2133  class Test_LinearPDE_noLumping(Test_line
2133          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2134          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2135          u=mypde.getSolution()          u=mypde.getSolution()
2136          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
2137      def test_symmetryOnDirect_System(self):      def test_symmetryOnDirect_System(self):
2138          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2139          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 1882  class Test_LinearPDE_noLumping(Test_line Line 2147  class Test_LinearPDE_noLumping(Test_line
2147          mypde.getSolverOptions().setSolverMethod(SolverOptions.DIRECT)          mypde.getSolverOptions().setSolverMethod(SolverOptions.DIRECT)
2148          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2149          u=mypde.getSolution()          u=mypde.getSolution()
2150          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
2151      def test_PCG_JACOBI_System(self):      def test_PCG_JACOBI_System(self):
2152          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2153          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 1894  class Test_LinearPDE_noLumping(Test_line Line 2159  class Test_LinearPDE_noLumping(Test_line
2159          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2160          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2161          mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)          mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
2162      mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)          mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
2163            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2164            u=mypde.getSolution()
2165            self.assertTrue(self.check(u,1.),'solution is wrong.')
2166        def test_PCG_GAUSS_SEIDEL_System(self):
2167            A=Tensor4(0.,Function(self.domain))
2168            D=Tensor(1.,Function(self.domain))
2169            Y=Vector(self.domain.getDim(),Function(self.domain))
2170            for i in range(self.domain.getDim()):
2171                A[i,:,i,:]=kronecker(self.domain)
2172                D[i,i]+=i
2173                Y[i]+=i
2174            mypde=LinearPDE(self.domain,debug=self.DEBUG)
2175            mypde.setValue(A=A,D=D,Y=Y)
2176            mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
2177            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
2178          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2179          u=mypde.getSolution()          u=mypde.getSolution()
2180          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
2181        def test_PCG_AMG_System(self):
2182            if self.order!=2:
2183                if getEscriptParamInt('DISABLE_AMG', 0):
2184                    print("AMG test disabled on MPI build")
2185                    return
2186                A=Tensor4(0.,Function(self.domain))
2187                D=Tensor(1.,Function(self.domain))
2188                Y=Vector(self.domain.getDim(),Function(self.domain))
2189                for i in range(self.domain.getDim()):
2190                    A[i,:,i,:]=kronecker(self.domain)
2191                    D[i,i]+=i
2192                    Y[i]+=i
2193                mypde=LinearPDE(self.domain,debug=self.DEBUG)
2194                mypde.setValue(A=A,D=D,Y=Y)
2195                mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
2196                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
2197                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2198                u=mypde.getSolution()
2199                self.assertTrue(self.check(u,1.),'solution is wrong.')
2200      def test_PCG_ILU0_System(self):      def test_PCG_ILU0_System(self):
2201          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2202          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 1909  class Test_LinearPDE_noLumping(Test_line Line 2208  class Test_LinearPDE_noLumping(Test_line
2208          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2209          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2210          mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)          mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
2211      mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
2212          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2213          u=mypde.getSolution()          u=mypde.getSolution()
2214          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
2215      def test_DIRECT_System(self):      def test_DIRECT_System(self):
2216          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2217          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 1926  class Test_LinearPDE_noLumping(Test_line Line 2225  class Test_LinearPDE_noLumping(Test_line
2225          mypde.getSolverOptions().setSolverMethod(SolverOptions.DIRECT)          mypde.getSolverOptions().setSolverMethod(SolverOptions.DIRECT)
2226          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2227          u=mypde.getSolution()          u=mypde.getSolution()
2228          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
2229      def test_BICGSTAB_JACOBI_System(self):      def test_BICGSTAB_JACOBI_System(self):
2230          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2231          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 1937  class Test_LinearPDE_noLumping(Test_line Line 2236  class Test_LinearPDE_noLumping(Test_line
2236              Y[i]+=i              Y[i]+=i
2237          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2238          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2239      mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)          mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)
2240      mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)          mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
2241          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2242          u=mypde.getSolution()          u=mypde.getSolution()
2243          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
2244        def test_BICGSTAB_GAUSS_SEIDEL_System(self):
2245            A=Tensor4(0.,Function(self.domain))
2246            D=Tensor(1.,Function(self.domain))
2247            Y=Vector(self.domain.getDim(),Function(self.domain))
2248            for i in range(self.domain.getDim()):
2249                A[i,:,i,:]=kronecker(self.domain)
2250                D[i,i]+=i
2251                Y[i]+=i
2252            mypde=LinearPDE(self.domain,debug=self.DEBUG)
2253            mypde.setValue(A=A,D=D,Y=Y)
2254            mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)
2255            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
2256            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2257            u=mypde.getSolution()
2258            self.assertTrue(self.check(u,1.),'solution is wrong.')
2259        def test_BICGSTAB_AMG_System(self):
2260            if self.order!=2:
2261                if getEscriptParamInt('DISABLE_AMG', 0):
2262                    print("AMG test disabled on MPI build")
2263                    return    
2264                A=Tensor4(0.,Function(self.domain))
2265                D=Tensor(1.,Function(self.domain))
2266                Y=Vector(self.domain.getDim(),Function(self.domain))
2267                for i in range(self.domain.getDim()):
2268                    A[i,:,i,:]=kronecker(self.domain)
2269                    D[i,i]+=i
2270                    Y[i]+=i
2271                mypde=LinearPDE(self.domain,debug=self.DEBUG)
2272                mypde.setValue(A=A,D=D,Y=Y)
2273                mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)
2274                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
2275                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2276                u=mypde.getSolution()
2277                self.assertTrue(self.check(u,1.),'solution is wrong.')
2278      def test_BICGSTAB_ILU0_System(self):      def test_BICGSTAB_ILU0_System(self):
2279          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2280          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 1952  class Test_LinearPDE_noLumping(Test_line Line 2285  class Test_LinearPDE_noLumping(Test_line
2285              Y[i]+=i              Y[i]+=i
2286          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2287          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2288      mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)          mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)
2289      mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
2290          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2291          u=mypde.getSolution()          u=mypde.getSolution()
2292          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
2293      def test_PRES20_JACOBI_System(self):      def test_PRES20_JACOBI_System(self):
2294          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2295          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 1967  class Test_LinearPDE_noLumping(Test_line Line 2300  class Test_LinearPDE_noLumping(Test_line
2300              Y[i]+=i              Y[i]+=i
2301          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2302          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2303      mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)          mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)
2304      mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)          mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
2305          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2306          u=mypde.getSolution()          u=mypde.getSolution()
2307          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
2308        def test_PRES20_GAUSS_SEIDEL_System(self):
2309            A=Tensor4(0.,Function(self.domain))
2310            D=Tensor(1.,Function(self.domain))
2311            Y=Vector(self.domain.getDim(),Function(self.domain))
2312            for i in range(self.domain.getDim()):
2313                A[i,:,i,:]=kronecker(self.domain)
2314                D[i,i]+=i
2315                Y[i]+=i
2316            mypde=LinearPDE(self.domain,debug=self.DEBUG)
2317            mypde.setValue(A=A,D=D,Y=Y)
2318            mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)
2319            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
2320            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2321            u=mypde.getSolution()
2322            self.assertTrue(self.check(u,1.),'solution is wrong.')
2323        def test_PRES20_AMG_System(self):
2324            if self.order!=2:
2325                if getEscriptParamInt('DISABLE_AMG', 0):
2326                    print("AMG test disabled on MPI build")
2327                    return
2328                A=Tensor4(0.,Function(self.domain))
2329                D=Tensor(1.,Function(self.domain))
2330                Y=Vector(self.domain.getDim(),Function(self.domain))
2331                for i in range(self.domain.getDim()):
2332                    A[i,:,i,:]=kronecker(self.domain)
2333                    D[i,i]+=i
2334                    Y[i]+=i
2335                mypde=LinearPDE(self.domain,debug=self.DEBUG)
2336                mypde.setValue(A=A,D=D,Y=Y)
2337                mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)
2338                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
2339                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2340                u=mypde.getSolution()
2341                self.assertTrue(self.check(u,1.),'solution is wrong.')
2342      def test_PRES20_ILU0_System(self):      def test_PRES20_ILU0_System(self):
2343          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2344          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 1982  class Test_LinearPDE_noLumping(Test_line Line 2349  class Test_LinearPDE_noLumping(Test_line
2349              Y[i]+=i              Y[i]+=i
2350          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2351          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2352      mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)          mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)
2353      mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
2354          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2355          u=mypde.getSolution()          u=mypde.getSolution()
2356          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
2357      def test_GMRESnoRestart_JACOBI_System(self):      def test_GMRESnoRestart_JACOBI_System(self):
2358          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2359          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 1997  class Test_LinearPDE_noLumping(Test_line Line 2364  class Test_LinearPDE_noLumping(Test_line
2364              Y[i]+=i              Y[i]+=i
2365          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2366          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2367      mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2368      mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)          mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
2369          # u=mypde.getSolution(verbose=self.VERBOSE,truncation=5)          # u=mypde.getSolution(verbose=self.VERBOSE,truncation=5)
2370          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2371          u=mypde.getSolution()          u=mypde.getSolution()
2372          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
2373        def test_GMRESnoRestart_GAUSS_SEIDEL_System(self):
2374            A=Tensor4(0.,Function(self.domain))
2375            D=Tensor(1.,Function(self.domain))
2376            Y=Vector(self.domain.getDim(),Function(self.domain))
2377            for i in range(self.domain.getDim()):
2378                A[i,:,i,:]=kronecker(self.domain)
2379                D[i,i]+=i
2380                Y[i]+=i
2381            mypde=LinearPDE(self.domain,debug=self.DEBUG)
2382            mypde.setValue(A=A,D=D,Y=Y)
2383            mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2384            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
2385            # u=mypde.getSolution(verbose=self.VERBOSE,truncation=5)
2386            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2387            u=mypde.getSolution()
2388            self.assertTrue(self.check(u,1.),'solution is wrong.')
2389        def test_GMRESnoRestart_AMG_System(self):
2390            if self.order!=2:
2391                if getEscriptParamInt('DISABLE_AMG',0):
2392                    print("AMG test disabled on MPI build")
2393                    return        
2394                A=Tensor4(0.,Function(self.domain))
2395                D=Tensor(1.,Function(self.domain))
2396                Y=Vector(self.domain.getDim(),Function(self.domain))
2397                for i in range(self.domain.getDim()):
2398                    A[i,:,i,:]=kronecker(self.domain)
2399                    D[i,i]+=i
2400                    Y[i]+=i
2401                mypde=LinearPDE(self.domain,debug=self.DEBUG)
2402                mypde.setValue(A=A,D=D,Y=Y)
2403                mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2404                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
2405                # u=mypde.getSolution(verbose=self.VERBOSE,truncation=5)
2406                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2407                u=mypde.getSolution()
2408                self.assertTrue(self.check(u,1.),'solution is wrong.')
2409      def test_GMRESnoRestart_ILU0_System(self):      def test_GMRESnoRestart_ILU0_System(self):
2410          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2411          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 2013  class Test_LinearPDE_noLumping(Test_line Line 2416  class Test_LinearPDE_noLumping(Test_line
2416              Y[i]+=i              Y[i]+=i
2417          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2418          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2419      mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2420      mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
2421          # u=mypde.getSolution(verbose=self.VERBOSE,truncation=5)          # u=mypde.getSolution(verbose=self.VERBOSE,truncation=5)
2422          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2423          u=mypde.getSolution()          u=mypde.getSolution()
2424          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
2425      def test_GMRES_JACOBI_System(self):      def test_GMRES_JACOBI_System(self):
2426          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2427          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 2029  class Test_LinearPDE_noLumping(Test_line Line 2432  class Test_LinearPDE_noLumping(Test_line
2432              Y[i]+=i              Y[i]+=i
2433          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2434          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2435      mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2436      mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)          mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
2437            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2438            u=mypde.getSolution()
2439            self.assertTrue(self.check(u,1.),'solution is wrong.')
2440        def test_GMRES_GAUSS_SEIDEL_System(self):
2441            A=Tensor4(0.,Function(self.domain))
2442            D=Tensor(1.,Function(self.domain))
2443            Y=Vector(self.domain.getDim(),Function(self.domain))
2444            for i in range(self.domain.getDim()):
2445                A[i,:,i,:]=kronecker(self.domain)
2446                D[i,i]+=i
2447                Y[i]+=i
2448            mypde=LinearPDE(self.domain,debug=self.DEBUG)
2449            mypde.setValue(A=A,D=D,Y=Y)
2450            mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2451            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
2452          mypde.getSolverOptions().setVerbosity(self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2453          u=mypde.getSolution()          u=mypde.getSolution()
2454          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
2455        def test_GMRES_AMG_System(self):
2456            if self.order!=2:
2457                if getEscriptParamInt('DISABLE_AMG', 0):
2458                    print("AMG test disabled on MPI build")
2459                    return    
2460                A=Tensor4(0.,Function(self.domain))
2461                D=Tensor(1.,Function(self.domain))
2462                Y=Vector(self.domain.getDim(),Function(self.domain))
2463                for i in range(self.domain.getDim()):
2464                    A[i,:,i,:]=kronecker(self.domain)
2465                    D[i,i]+=i
2466                    Y[i]+=i
2467                mypde=LinearPDE(self.domain,debug=self.DEBUG)
2468                mypde.setValue(A=A,D=D,Y=Y)
2469                mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2470     &nbs