/[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 2470 by gross, Thu Jun 11 08:32:32 2009 UTC trunk/escriptcore/test/python/test_linearPDEs.py revision 4783 by sshaw, Thu Mar 20 06:23:23 2014 UTC
# Line 1  Line 1 
1    # -*- coding: utf-8 -*-
2    
3  ########################################################  ##############################################################################
4  #  #
5  # Copyright (c) 2003-2008 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 LinearPDE,IllegalCoefficientValue,Poisson, IllegalCoefficientFunctionSpace, TransportPDE, IllegalCoefficient, Helmholtz, LameEquation, SolverOptions
35  import numpy  import numpy
36  import unittest  import unittest
# Line 55  class Test_linearPDEs(unittest.TestCase) Line 40  class Test_linearPDEs(unittest.TestCase)
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(mypde, ["d_contact", "y_contact"])
72                self.checkIfReducedNotEmpty(mypde, ["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.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 nor 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.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.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
391    
392      def test_SolverOptions(self):      def test_SolverOptions(self):
393          so=SolverOptions()          so=SolverOptions()
394            
395            self.assertTrue(so.getSmoother() == 28, "initial Smoother is wrong.")
396            self.assertRaises(ValueError,so.setSmoother,-1)
397            so.setSmoother(so.GAUSS_SEIDEL)
398            self.assertTrue(so.getSmoother() == 28, "Gauss-Seidel smoother is not set.")
399            so.setSmoother(so.JACOBI)
400            self.assertTrue(so.getSmoother() == 10, "Jacobi smoother is not set.")
401    
402            self.assertTrue(so.getLevelMax() == 100, "initial  LevelMax is wrong.")
403            self.assertRaises(ValueError,so.setLevelMax,-1)
404            so.setLevelMax(20)
405            self.assertTrue(so.getLevelMax() == 20, "LevelMax is wrong.")
406    
407          self.failUnless(so.getLevelMax() == 10, "initial  LevelMax is wrong.")          self.assertTrue(so.getCoarseningThreshold() == 0.25, "initial  CoarseningThreshold is wrong.")
408          self.failUnlessRaises(ValueError,so.setLevelMax,-1)          self.assertRaises(ValueError,so.setCoarseningThreshold,-1)
         so.setLevelMax(3)  
         self.failUnless(so.getLevelMax() == 3, "LevelMax is wrong.")  
   
         self.failUnless(so.getCoarseningThreshold() == 0.05, "initial  CoarseningThreshold is wrong.")  
         self.failUnlessRaises(ValueError,so.setCoarseningThreshold,-1)  
409          so.setCoarseningThreshold(0.1)          so.setCoarseningThreshold(0.1)
410          self.failUnless(so.getCoarseningThreshold() == 0.1, "CoarseningThreshold is wrong.")          self.assertTrue(so.getCoarseningThreshold() == 0.1, "CoarseningThreshold is wrong.")
411            
412            self.assertTrue(so.getMinCoarseMatrixSize() == 500, "initial  Minimum Coarse Matrix Size is wrong.")
413            self.assertRaises(ValueError,so.setMinCoarseMatrixSize,-1)
414            so.setMinCoarseMatrixSize(1000)
415            self.assertTrue(so.getMinCoarseMatrixSize() == 1000, "Minimum Coarse Matrix Size is wrong.")
416    
417          self.failUnless(so.getNumSweeps() == 2, "initial  Sweeps is wrong.")          self.assertTrue(so.getNumSweeps() == 1, "initial  Sweeps is wrong.")
418          self.failUnlessRaises(ValueError,so.setNumSweeps,-1)          self.assertRaises(ValueError,so.setNumSweeps,-1)
419          so.setNumSweeps(3)          so.setNumSweeps(3)
420          self.failUnless(so.getNumSweeps() == 3, "Sweeps is wrong.")          self.assertTrue(so.getNumSweeps() == 3, "Sweeps is wrong.")
421    
422          self.failUnless(so.getNumPreSweeps() == 2, "initial  PreSweeps is wrong.")          self.assertTrue(so.getNumPreSweeps() == 1, "initial  PreSweeps is wrong.")
423          self.failUnlessRaises(ValueError,so.setNumPreSweeps,-1)          self.assertRaises(ValueError,so.setNumPreSweeps,-1)
424          so.setNumPreSweeps(4)          so.setNumPreSweeps(4)
425          self.failUnless(so.getNumPreSweeps() == 4, "PreSweeps is wrong.")          self.assertTrue(so.getNumPreSweeps() == 4, "PreSweeps is wrong.")
426    
427          self.failUnless(so.getNumPostSweeps() == 2, "initial  PreSweeps is wrong.")          self.assertTrue(so.getNumPostSweeps() == 1, "initial  PostSweeps is wrong.")
428          self.failUnlessRaises(ValueError,so.setNumPostSweeps,-1)          self.assertRaises(ValueError,so.setNumPostSweeps,-1)
429          so.setNumPostSweeps(5)          so.setNumPostSweeps(5)
430          self.failUnless(so.getNumPostSweeps() == 5, "PreSweeps is wrong.")          self.assertTrue(so.getNumPostSweeps() == 5, "PostSweeps is wrong.")
431    
432          self.failUnless(so.getTolerance() == 1.e-8, "initial Tolerance is wrong.")          self.assertTrue(so.getTolerance() == 1.e-8, "initial Tolerance is wrong.")
433          self.failUnlessRaises(ValueError,so.setTolerance,-1)          self.assertRaises(ValueError,so.setTolerance,-1)
434          so.setTolerance(0.2)          so.setTolerance(0.2)
435          self.failUnless(so.getTolerance() == 0.2, "Tolerance is wrong.")          self.assertTrue(so.getTolerance() == 0.2, "Tolerance is wrong.")
436    
437          self.failUnless(so.getAbsoluteTolerance() == 0., "initial  AbsoluteTolerance is wrong.")          self.assertTrue(so.getAbsoluteTolerance() == 0., "initial  AbsoluteTolerance is wrong.")
438          self.failUnlessRaises(ValueError,so.setAbsoluteTolerance,-1)          self.assertRaises(ValueError,so.setAbsoluteTolerance,-1)
439          so.setAbsoluteTolerance(0.3)          so.setAbsoluteTolerance(0.3)
440          self.failUnless(so.getAbsoluteTolerance() == 0.3, "AbsoluteTolerance is wrong.")          self.assertTrue(so.getAbsoluteTolerance() == 0.3, "AbsoluteTolerance is wrong.")
441    
442          self.failUnless(so.getInnerTolerance() == 0.9, "initial InnerTolerance is wrong.")          self.assertTrue(so.getInnerTolerance() == 0.9, "initial InnerTolerance is wrong.")
443          self.failUnlessRaises(ValueError,so.setInnerTolerance,-1)          self.assertRaises(ValueError,so.setInnerTolerance,-1)
444          so.setInnerTolerance(0.4)          so.setInnerTolerance(0.4)
445          self.failUnless(so.getInnerTolerance() == 0.4, "InnerTolerance is wrong.")          self.assertTrue(so.getInnerTolerance() == 0.4, "InnerTolerance is wrong.")
446    
447          self.failUnless(so.getDropTolerance() == 0.01, "initial DropTolerance is wrong.")          self.assertTrue(so.getDropTolerance() == 0.01, "initial DropTolerance is wrong.")
448          self.failUnlessRaises(ValueError,so.setDropTolerance,-1)          self.assertRaises(ValueError,so.setDropTolerance,-1)
449          so.setDropTolerance(0.5)          so.setDropTolerance(0.5)
450          self.failUnless(so.getDropTolerance() == 0.5, "DropDropTolerance is wrong.")          self.assertTrue(so.getDropTolerance() == 0.5, "DropDropTolerance is wrong.")
451    
452          self.failUnless(so.getDropStorage() == 2., "initial DropStorage is wrong.")          self.assertTrue(so.getDropStorage() == 2., "initial DropStorage is wrong.")
453          self.failUnlessRaises(ValueError,so.setDropStorage,-1)          self.assertRaises(ValueError,so.setDropStorage,-1)
454          so.setDropStorage(10)          so.setDropStorage(10)
455          self.failUnless(so.getDropStorage() == 10, "DropStorage is wrong.")          self.assertTrue(so.getDropStorage() == 10, "DropStorage is wrong.")
456                    
457          self.failUnless(so.getRelaxationFactor() == 0.3, "initial RelaxationFactor is wrong.")          self.assertTrue(so.getRelaxationFactor() == 0.3, "initial RelaxationFactor is wrong.")
458          self.failUnlessRaises(ValueError,so.setRelaxationFactor,-1)          self.assertRaises(ValueError,so.setRelaxationFactor,-1)
459          so.setRelaxationFactor(0.1)          so.setRelaxationFactor(0.1)
460          self.failUnless(so.getRelaxationFactor() == 0.1, "Relaxation is wrong.")          self.assertTrue(so.getRelaxationFactor() == 0.1, "Relaxation is wrong.")
461    
462    
463          self.failUnless(so.getIterMax() == 10000, "initial IterMax is wrong.")          self.assertTrue(so.getIterMax() == 100000, "initial IterMax is wrong.")
464          self.failUnlessRaises(ValueError,so.setIterMax,0)          self.assertRaises(ValueError,so.setIterMax,0)
465          so.setIterMax(11)          so.setIterMax(11)
466          self.failUnless(so.getIterMax() == 11, "IterMax is wrong.")          self.assertTrue(so.getIterMax() == 11, "IterMax is wrong.")
467    
468          self.failUnless(so.getInnerIterMax() == 10, "initial InnerIterMax is wrong.")          self.assertTrue(so.getInnerIterMax() == 10, "initial InnerIterMax is wrong.")
469          self.failUnlessRaises(ValueError,so.setInnerIterMax,0)          self.assertRaises(ValueError,so.setInnerIterMax,0)
470          so.setInnerIterMax(12)          so.setInnerIterMax(12)
471          self.failUnless(so.getInnerIterMax() == 12, "InnerIterMax is wrong.")          self.assertTrue(so.getInnerIterMax() == 12, "InnerIterMax is wrong.")
472    
473          self.failUnless(so.getTruncation() == 20, "initial Truncation is wrong.")          self.assertTrue(so.getTruncation() == 20, "initial Truncation is wrong.")
474          self.failUnlessRaises(ValueError,so.setTruncation,0)          self.assertRaises(ValueError,so.setTruncation,0)
475          so.setTruncation(13)          so.setTruncation(13)
476          self.failUnless(so.getTruncation() == 13, "Truncation is wrong.")          self.assertTrue(so.getTruncation() == 13, "Truncation is wrong.")
477    
478          self.failUnless(so.getRestart() == None, "initial Truncation is wrong.")          self.assertTrue(so.getRestart() == None, "initial Truncation is wrong.")
479          self.failUnlessRaises(ValueError,so.setTruncation,0)          self.assertRaises(ValueError,so.setTruncation,0)
480          so.setRestart(14)          so.setRestart(14)
481          self.failUnless(so.getRestart() == 14, "Truncation is wrong.")          self.assertTrue(so.getRestart() == 14, "Truncation is wrong.")
482          so.setRestart(None)          so.setRestart(None)
483          self.failUnless(so.getRestart() == None, "Truncation is wrong.")          self.assertTrue(so.getRestart() == None, "Truncation is wrong.")
484                
485          self.failUnless(not so.isSymmetric(), "initial symmetry flag is wrong.")          self.assertTrue(not so.isVerbose(), "initial verbosity flag is wrong.")
486            so.setVerbosityOn()
487            self.assertTrue(so.isVerbose(), "verbosity (1) flag is wrong.")
488            so.setVerbosityOff()
489            self.assertTrue(not so.isVerbose(), "verbosity (2) flag is wrong.")
490            so.setVerbosity(verbose=True)
491            self.assertTrue(so.isVerbose(), "verbosity (3) flag is wrong.")
492            so.setVerbosity(verbose=False)
493            self.assertTrue(not so.isVerbose(), "verbosity (4) flag is wrong.")
494    
495            self.assertTrue(not so.isSymmetric(), "initial symmetry flag is wrong.")
496          so.setSymmetryOn()          so.setSymmetryOn()
497          self.failUnless(so.isSymmetric(), "symmetry (1) flag is wrong.")          self.assertTrue(so.isSymmetric(), "symmetry (1) flag is wrong.")
498          so.setSymmetryOff()          so.setSymmetryOff()
499          self.failUnless(not so.isSymmetric(), "symmetry (2) flag is wrong.")          self.assertTrue(not so.isSymmetric(), "symmetry (2) flag is wrong.")
500          so.setSymmetry(flag=True)          so.setSymmetry(flag=True)
501          self.failUnless(so.isSymmetric(), "symmetry (3) flag is wrong.")          self.assertTrue(so.isSymmetric(), "symmetry (3) flag is wrong.")
502          so.setSymmetry(flag=False)          so.setSymmetry(flag=False)
503          self.failUnless(not so.isSymmetric(), "symmetry (4) flag is wrong.")          self.assertTrue(not so.isSymmetric(), "symmetry (4) flag is wrong.")
504    
505          self.failUnless(so.adaptInnerTolerance(), "initial InnerToleranceAdaption flag is wrong.")          self.assertTrue(so.adaptInnerTolerance(), "initial InnerToleranceAdaption flag is wrong.")
506          so.setInnerToleranceAdaptionOn()          so.setInnerToleranceAdaptionOn()
507          self.failUnless(so.adaptInnerTolerance(), "InnerToleranceAdaption (1) flag is wrong.")          self.assertTrue(so.adaptInnerTolerance(), "InnerToleranceAdaption (1) flag is wrong.")
508          so.setInnerToleranceAdaptionOff()          so.setInnerToleranceAdaptionOff()
509          self.failUnless(not so.adaptInnerTolerance(), "InnerToleranceAdaption (2) flag is wrong.")          self.assertTrue(not so.adaptInnerTolerance(), "InnerToleranceAdaption (2) flag is wrong.")
510          so.setInnerToleranceAdaption(adapt=True)          so.setInnerToleranceAdaption(adapt=True)
511          self.failUnless(so.adaptInnerTolerance(), "InnerToleranceAdaption (3) flag is wrong.")          self.assertTrue(so.adaptInnerTolerance(), "InnerToleranceAdaption (3) flag is wrong.")
512          so.setInnerToleranceAdaption(adapt=False)          so.setInnerToleranceAdaption(adapt=False)
513          self.failUnless(not so.adaptInnerTolerance(), "InnerToleranceAdaption (4) flag is wrong.")          self.assertTrue(not so.adaptInnerTolerance(), "InnerToleranceAdaption (4) flag is wrong.")
514            
515          self.failUnless(not so.acceptConvergenceFailure(), "initial acceptConvergenceFailure flag is wrong.")          self.assertTrue(not so.acceptConvergenceFailure(), "initial acceptConvergenceFailure flag is wrong.")
516          so.setAcceptanceConvergenceFailureOn()          so.setAcceptanceConvergenceFailureOn()
517          self.failUnless(so.acceptConvergenceFailure(), "acceptConvergenceFailure (1) flag is wrong.")          self.assertTrue(so.acceptConvergenceFailure(), "acceptConvergenceFailure (1) flag is wrong.")
518          so.setAcceptanceConvergenceFailureOff()          so.setAcceptanceConvergenceFailureOff()
519          self.failUnless(not so.acceptConvergenceFailure(), "acceptConvergenceFailure (2) flag is wrong.")          self.assertTrue(not so.acceptConvergenceFailure(), "acceptConvergenceFailure (2) flag is wrong.")
520          so.setAcceptanceConvergenceFailure(accept=True)          so.setAcceptanceConvergenceFailure(accept=True)
521          self.failUnless(so.acceptConvergenceFailure(), "acceptConvergenceFailure (3) flag is wrong.")          self.assertTrue(so.acceptConvergenceFailure(), "acceptConvergenceFailure (3) flag is wrong.")
522          so.setAcceptanceConvergenceFailure(accept=False)          so.setAcceptanceConvergenceFailure(accept=False)
523          self.failUnless(not so.acceptConvergenceFailure(), "acceptConvergenceFailure (4) flag is wrong.")            self.assertTrue(not so.acceptConvergenceFailure(), "acceptConvergenceFailure (4) flag is wrong.")  
524                    
525          self.failUnless(so.getReordering() == 30, "initial Reordering is wrong.")          self.assertTrue(so.getReordering() == 30, "initial Reordering is wrong.")
526          self.failUnlessRaises(ValueError,so.setReordering,-1)          self.assertRaises(ValueError,so.setReordering,-1)
527          so.setReordering(so.NO_REORDERING)          so.setReordering(so.NO_REORDERING)
528          self.failUnless(so.getReordering() == 17, "NO_REORDERING is not set.")          self.assertTrue(so.getReordering() == 17, "NO_REORDERING is not set.")
529          so.setReordering(so.MINIMUM_FILL_IN)          so.setReordering(so.MINIMUM_FILL_IN)
530          self.failUnless(so.getReordering() == 18, "MINIMUM_FILL_IN is not set.")          self.assertTrue(so.getReordering() == 18, "MINIMUM_FILL_IN is not set.")
531          so.setReordering(so.NESTED_DISSECTION)          so.setReordering(so.NESTED_DISSECTION)
532          self.failUnless(so.getReordering() == 19, "NESTED_DISSECTION is not set.")          self.assertTrue(so.getReordering() == 19, "NESTED_DISSECTION is not set.")
533          so.setReordering(so.DEFAULT_REORDERING)          so.setReordering(so.DEFAULT_REORDERING)
534          self.failUnless(so.getReordering() == 30, "DEFAULT_REORDERING is not set.")          self.assertTrue(so.getReordering() == 30, "DEFAULT_REORDERING is not set.")
535                    
536          self.failUnless(so.getPackage() == 0, "initial solver package is wrong.")          self.assertTrue(so.getPackage() == 0, "initial solver package is wrong.")
537          self.failUnlessRaises(ValueError,so.setPackage,-1)          self.assertRaises(ValueError,so.setPackage,-1)
538          so.setPackage(so.PASO)          so.setPackage(so.PASO)
539          self.failUnless(so.getPackage() == 21, "PASO is not set.")          self.assertTrue(so.getPackage() == 21, "PASO is not set.")
540          so.setPackage(so.SUPER_LU)          so.setPackage(so.SUPER_LU)
541          self.failUnless(so.getPackage() == 31, "SUPER_LU is not set.")          self.assertTrue(so.getPackage() == 31, "SUPER_LU is not set.")
542          so.setPackage(so.PASTIX)          so.setPackage(so.PASTIX)
543          self.failUnless(so.getPackage() == 32, "PASTIX is not set.")          self.assertTrue(so.getPackage() == 32, "PASTIX is not set.")
544          so.setPackage(so.MKL)          so.setPackage(so.MKL)
545          self.failUnless(so.getPackage() == 15, "MKL is not set.")          self.assertTrue(so.getPackage() == 15, "MKL is not set.")
546          so.setPackage(so.UMFPACK)          so.setPackage(so.UMFPACK)
547          self.failUnless(so.getPackage() == 16, "UMFPACK is not set.")          self.assertTrue(so.getPackage() == 16, "UMFPACK is not set.")
548          so.setPackage(so.TRILINOS)          so.setPackage(so.TRILINOS)
549          self.failUnless(so.getPackage() == 24, "TRILINOS is not set.")          self.assertTrue(so.getPackage() == 24, "TRILINOS is not set.")
550    
551          self.failUnless(so.getSolverMethod() == 0, "initial SolverMethod is wrong.")          self.assertTrue(so.getSolverMethod() == 0, "initial SolverMethod is wrong.")
552          self.failUnlessRaises(ValueError,so.setSolverMethod,-1)          self.assertRaises(ValueError,so.setSolverMethod,-1)
553          so.setSolverMethod(so.DIRECT)          so.setSolverMethod(so.DIRECT)
554          self.failUnless(so.getSolverMethod() == 1, "DIRECT is not set.")          self.assertTrue(so.getSolverMethod() == 1, "DIRECT is not set.")
555          so.setSolverMethod(so.CHOLEVSKY)          so.setSolverMethod(so.CHOLEVSKY)
556          self.failUnless(so.getSolverMethod() == 2, "CHOLEVSKY is not set.")          self.assertTrue(so.getSolverMethod() == 2, "CHOLEVSKY is not set.")
557          so.setSolverMethod(so.PCG)          so.setSolverMethod(so.PCG)
558          self.failUnless(so.getSolverMethod() == 3, "PCG is not set.")          self.assertTrue(so.getSolverMethod() == 3, "PCG is not set.")
559          so.setSolverMethod(so.CR)          so.setSolverMethod(so.CR)
560          self.failUnless(so.getSolverMethod() == 4, "CR is not set.")          self.assertTrue(so.getSolverMethod() == 4, "CR is not set.")
561          so.setSolverMethod(so.CGS)          so.setSolverMethod(so.CGS)
562          self.failUnless(so.getSolverMethod() == 5, "CGS is not set.")          self.assertTrue(so.getSolverMethod() == 5, "CGS is not set.")
563          so.setSolverMethod(so.BICGSTAB)          so.setSolverMethod(so.BICGSTAB)
564          self.failUnless(so.getSolverMethod() == 6, "BICGSTAB is not set.")          self.assertTrue(so.getSolverMethod() == 6, "BICGSTAB is not set.")
         so.setSolverMethod(so.SSOR)  
         self.failUnless(so.getSolverMethod() == 7, "SSOR is not set.")  
565          so.setSolverMethod(so.GMRES)          so.setSolverMethod(so.GMRES)
566          self.failUnless(so.getSolverMethod() == 11, "GMRES is not set.")          self.assertTrue(so.getSolverMethod() == 11, "GMRES is not set.")
567          so.setSolverMethod(so.PRES20)          so.setSolverMethod(so.PRES20)
568          self.failUnless(so.getSolverMethod() == 12, "PRES20 is not set.")          self.assertTrue(so.getSolverMethod() == 12, "PRES20 is not set.")
569          so.setSolverMethod(so.LUMPING)          so.setSolverMethod(so.LUMPING)
570          self.failUnless(so.getSolverMethod() == 13, "LUMPING is not set.")          self.assertTrue(so.getSolverMethod() == 13, "LUMPING is not set.")
571          so.setSolverMethod(so.ITERATIVE)          so.setSolverMethod(so.ITERATIVE)
572          self.failUnless(so.getSolverMethod() == 20, "ITERATIVE is not set.")          self.assertTrue(so.getSolverMethod() == 20, "ITERATIVE is not set.")
         so.setSolverMethod(so.AMG)  
         self.failUnless(so.getSolverMethod() == 22, "AMG is not set.")  
573          so.setSolverMethod(so.NONLINEAR_GMRES)          so.setSolverMethod(so.NONLINEAR_GMRES)
574          self.failUnless(so.getSolverMethod() == 25, "NONLINEAR_GMRES is not set.")          self.assertTrue(so.getSolverMethod() == 25, "NONLINEAR_GMRES is not set.")
575          so.setSolverMethod(so.TFQMR)          so.setSolverMethod(so.TFQMR)
576          self.failUnless(so.getSolverMethod() == 26, "TFQMR is not set.")          self.assertTrue(so.getSolverMethod() == 26, "TFQMR is not set.")
577          so.setSolverMethod(so.MINRES)          so.setSolverMethod(so.MINRES)
578          self.failUnless(so.getSolverMethod() == 27, "MINRES is not set.")          self.assertTrue(so.getSolverMethod() == 27, "MINRES is not set.")
         so.setSolverMethod(so.GAUSS_SEIDEL)  
         self.failUnless(so.getSolverMethod() == 28, "GAUSS_SEIDEL is not set.")  
579          so.setSolverMethod(so.DEFAULT)          so.setSolverMethod(so.DEFAULT)
580          self.failUnless(so.getSolverMethod() == 0, "DEFAULT is not set.")          self.assertTrue(so.getSolverMethod() == 0, "DEFAULT is not set.")
581    
582          self.failUnless(so.getPreconditioner() == 10, "initial Preconditioner is wrong.")          self.assertTrue(so.getPreconditioner() == 10, "initial Preconditioner is wrong.")
583          self.failUnlessRaises(ValueError,so.setPreconditioner,-1)          self.assertRaises(ValueError,so.setPreconditioner,-1)
584          so.setPreconditioner(so.ILU0)          so.setPreconditioner(so.ILU0)
585          self.failUnless(so.getPreconditioner() == 8, "ILU0 is not set.")          self.assertTrue(so.getPreconditioner() == 8, "ILU0 is not set.")
         so.setPreconditioner(so.SSOR)  
         self.failUnless(so.getPreconditioner() == 7, "SSOR is not set.")  
586          so.setPreconditioner(so.ILUT)          so.setPreconditioner(so.ILUT)
587          self.failUnless(so.getPreconditioner() == 9, "ILUT is not set.")          self.assertTrue(so.getPreconditioner() == 9, "ILUT is not set.")
588          so.setPreconditioner(so.JACOBI)          so.setPreconditioner(so.JACOBI)
589          self.failUnless(so.getPreconditioner() == 10, "JACOBI is not set.")          self.assertTrue(so.getPreconditioner() == 10, "JACOBI is not set.")
590          so.setPreconditioner(so.AMG)          if getEscriptParamInt('DISABLE_AMG', 0):
591          self.failUnless(so.getPreconditioner() == 22, "AMG is not set.")              print("AMG test disabled on MPI build")
592            else:
593                so.setPreconditioner(so.AMG)
594                self.assertTrue(so.getPreconditioner() == 22, "AMG is not set.")
595          so.setPreconditioner(so.REC_ILU)          so.setPreconditioner(so.REC_ILU)
596          self.failUnless(so.getPreconditioner() == 23, "REC_ILU is not set.")          self.assertTrue(so.getPreconditioner() == 23, "REC_ILU is not set.")
597          so.setPreconditioner(so.GAUSS_SEIDEL)          so.setPreconditioner(so.GAUSS_SEIDEL)
598          self.failUnless(so.getPreconditioner() == 28, "GAUSS_SEIDEL is not set.")          self.assertTrue(so.getPreconditioner() == 28, "GAUSS_SEIDEL is not set.")
599          so.setPreconditioner(so.RILU)          so.setPreconditioner(so.RILU)
600          self.failUnless(so.getPreconditioner() == 29, "RILU is not set.")          self.assertTrue(so.getPreconditioner() == 29, "RILU is not set.")
601            so.setPreconditioner(so.AMLI)
602            self.assertTrue(so.getPreconditioner() == 38, "AMLI is not set.")
603          so.setPreconditioner(so.NO_PRECONDITIONER)          so.setPreconditioner(so.NO_PRECONDITIONER)
604          self.failUnless(so.getPreconditioner() == 36, "NO_PRECONDITIONER is not set.")                  self.assertTrue(so.getPreconditioner() == 36, "NO_PRECONDITIONER is not set.")        
605    
606          self.failUnless(so.getCoarsening() == 0, "initial Coarseningr is wrong.")          self.assertTrue(so.getCoarsening() == 0, "initial Coarseningr is wrong.")
607          self.failUnlessRaises(ValueError,so.setCoarsening,-1)          self.assertRaises(ValueError,so.setCoarsening,-1)
608          so.setCoarsening(so.YAIR_SHAPIRA_COARSENING)          so.setCoarsening(so.YAIR_SHAPIRA_COARSENING)
609          self.failUnless(so.getCoarsening() == 33, "YAIR_SHAPIRA_COARSENING is not set.")          self.assertTrue(so.getCoarsening() == 33, "YAIR_SHAPIRA_COARSENING is not set.")
610          so.setCoarsening(so.RUGE_STUEBEN_COARSENING)          so.setCoarsening(so.RUGE_STUEBEN_COARSENING)
611          self.failUnless(so.getCoarsening() == 34, "RUGE_STUEBEN_COARSENING is not set.")          self.assertTrue(so.getCoarsening() == 34, "RUGE_STUEBEN_COARSENING is not set.")
612          so.setCoarsening(so.AGGREGATION_COARSENING)          so.setCoarsening(so.AGGREGATION_COARSENING)
613          self.failUnless(so.getCoarsening() == 35, "AGREGATION_COARSENING is not set.")          self.assertTrue(so.getCoarsening() == 35, "AGREGATION_COARSENING is not set.")
614            so.setCoarsening(so.STANDARD_COARSENING)
615            self.assertTrue(so.getCoarsening() == 39, "STANDARD_COARSENING is not set.")
616          so.setCoarsening(so.DEFAULT)          so.setCoarsening(so.DEFAULT)
617          self.failUnless(so.getCoarsening() == 0, "DEFAULT is not set.")          self.assertTrue(so.getCoarsening() == 0, "DEFAULT is not set.")
618    
619          self.failUnless(so.getDiagnostics("num_iter") == None, "initial num_iter is wrong.")          self.assertTrue(so.getDiagnostics("num_iter") == None, "initial num_iter is wrong.")
620          self.failUnless(so.getDiagnostics("num_inner_iter") == None, "initial num_inner_iter is wrong.")          self.assertTrue(so.getDiagnostics("num_inner_iter") == None, "initial num_inner_iter is wrong.")
621          self.failUnless(so.getDiagnostics("time") == None, "initial time is wrong.")          self.assertTrue(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(so.getDiagnostics("set_up_time") == None, "initial set_up_time is wrong.")
623          self.failUnless(so.getDiagnostics("residual_norm") == None, "initial residual_norm is wrong.")          self.assertTrue(so.getDiagnostics("residual_norm") == None, "initial residual_norm is wrong.")
624          self.failUnless(so.getDiagnostics("converged") == None, "initial converged is wrong.")          self.assertTrue(so.getDiagnostics("converged") == None, "initial converged is wrong.")
625          self.failUnless(so.hasConverged() == None, "initial convergence flag is wrong.")          self.assertTrue(so.hasConverged() == None, "initial convergence flag is wrong.")
626          self.failUnless(so.getDiagnostics("cum_num_inner_iter") == 0, "initial cum_num_inner_iter is wrong.")          self.assertTrue(so.getDiagnostics("cum_num_inner_iter") == 0, "initial cum_num_inner_iter is wrong.")
627          self.failUnless(so.getDiagnostics("cum_num_iter") == 0, "initial cum_num_iter is wrong.")          self.assertTrue(so.getDiagnostics("cum_num_iter") == 0, "initial cum_num_iter is wrong.")
628          self.failUnless(so.getDiagnostics("cum_time") ==0, "initial cum_time is wrong.")          self.assertTrue(so.getDiagnostics("cum_time") ==0, "initial cum_time is wrong.")
629          self.failUnless(so.getDiagnostics("cum_set_up_time") == 0, "initial cum_set_up_time is wrong.")          self.assertTrue(so.getDiagnostics("cum_set_up_time") == 0, "initial cum_set_up_time is wrong.")
630    
631          so._updateDiagnostics("num_iter",1)          so._updateDiagnostics("num_iter",1)
632          so._updateDiagnostics("num_inner_iter",2)          so._updateDiagnostics("num_inner_iter",2)
# Line 865  class Test_LinearPDE_noLumping(Test_line Line 635  class Test_LinearPDE_noLumping(Test_line
635          so._updateDiagnostics("residual_norm",5)          so._updateDiagnostics("residual_norm",5)
636          so._updateDiagnostics("converged",True)          so._updateDiagnostics("converged",True)
637    
638          self.failUnless(so.getDiagnostics("num_iter") == 1, "num_iter is wrong.")          self.assertTrue(so.getDiagnostics("num_iter") == 1, "num_iter is wrong.")
639          self.failUnless(so.getDiagnostics("num_inner_iter") == 2, "num_inner_iter is wrong.")          self.assertTrue(so.getDiagnostics("num_inner_iter") == 2, "num_inner_iter is wrong.")
640          self.failUnless(so.getDiagnostics("time") == 3, "time is wrong.")          self.assertTrue(so.getDiagnostics("time") == 3, "time is wrong.")
641          self.failUnless(so.getDiagnostics("set_up_time") == 4, "set_up_time is wrong.")          self.assertTrue(so.getDiagnostics("set_up_time") == 4, "set_up_time is wrong.")
642          self.failUnless(so.getDiagnostics("residual_norm") == 5, "residual_norm is wrong.")          self.assertTrue(so.getDiagnostics("residual_norm") == 5, "residual_norm is wrong.")
643          self.failUnless(so.getDiagnostics("converged"), "converged is wrong.")          self.assertTrue(so.getDiagnostics("converged"), "converged is wrong.")
644          self.failUnless(so.hasConverged(), "convergence flag is wrong.")          self.assertTrue(so.hasConverged(), "convergence flag is wrong.")
645          self.failUnless(so.getDiagnostics("cum_num_inner_iter") == 2, "cum_num_inner_iter is wrong.")          self.assertTrue(so.getDiagnostics("cum_num_inner_iter") == 2, "cum_num_inner_iter is wrong.")
646          self.failUnless(so.getDiagnostics("cum_num_iter") == 1, "cum_num_iter is wrong.")          self.assertTrue(so.getDiagnostics("cum_num_iter") == 1, "cum_num_iter is wrong.")
647          self.failUnless(so.getDiagnostics("cum_time") ==3, "cum_time is wrong.")          self.assertTrue(so.getDiagnostics("cum_time") ==3, "cum_time is wrong.")
648          self.failUnless(so.getDiagnostics("cum_set_up_time") == 4, "cum_set_up_time is wrong.")            self.assertTrue(so.getDiagnostics("cum_set_up_time") == 4, "cum_set_up_time is wrong.")  
649                    
650          so.resetDiagnostics()          so.resetDiagnostics()
651          self.failUnless(so.getDiagnostics("num_iter") == None, "initial num_iter is wrong.")          self.assertTrue(so.getDiagnostics("num_iter") == None, "initial num_iter is wrong.")
652          self.failUnless(so.getDiagnostics("num_inner_iter") == None, "initial num_inner_iter is wrong.")          self.assertTrue(so.getDiagnostics("num_inner_iter") == None, "initial num_inner_iter is wrong.")
653          self.failUnless(so.getDiagnostics("time") == None, "initial time is wrong.")          self.assertTrue(so.getDiagnostics("time") == None, "initial time is wrong.")
654          self.failUnless(so.getDiagnostics("set_up_time") == None, "initial set_up_time is wrong.")          self.assertTrue(so.getDiagnostics("set_up_time") == None, "initial set_up_time is wrong.")
655          self.failUnless(so.getDiagnostics("residual_norm") == None, "initial residual_norm is wrong.")          self.assertTrue(so.getDiagnostics("residual_norm") == None, "initial residual_norm is wrong.")
656          self.failUnless(so.getDiagnostics("converged") == None, "initial converged is wrong.")          self.assertTrue(so.getDiagnostics("converged") == None, "initial converged is wrong.")
657          self.failUnless(so.hasConverged() == None, "initial convergence flag is wrong")                self.assertTrue(so.hasConverged() == None, "initial convergence flag is wrong")      
658          self.failUnless(so.getDiagnostics("cum_num_inner_iter") == 2, "cum_num_inner_iter is wrong.")          self.assertTrue(so.getDiagnostics("cum_num_inner_iter") == 2, "cum_num_inner_iter is wrong.")
659          self.failUnless(so.getDiagnostics("cum_num_iter") == 1, "cum_num_iter is wrong.")          self.assertTrue(so.getDiagnostics("cum_num_iter") == 1, "cum_num_iter is wrong.")
660          self.failUnless(so.getDiagnostics("cum_time") ==3, "cum_time is wrong.")          self.assertTrue(so.getDiagnostics("cum_time") ==3, "cum_time is wrong.")
661          self.failUnless(so.getDiagnostics("cum_set_up_time") == 4, "cum_set_up_time is wrong.")          self.assertTrue(so.getDiagnostics("cum_set_up_time") == 4, "cum_set_up_time is wrong.")
662    
663          so._updateDiagnostics("num_iter",10)          so._updateDiagnostics("num_iter",10)
664          so._updateDiagnostics("num_inner_iter",20)          so._updateDiagnostics("num_inner_iter",20)
# Line 897  class Test_LinearPDE_noLumping(Test_line Line 667  class Test_LinearPDE_noLumping(Test_line
667          so._updateDiagnostics("residual_norm",50)          so._updateDiagnostics("residual_norm",50)
668          so._updateDiagnostics("converged",False)          so._updateDiagnostics("converged",False)
669    
670          self.failUnless(so.getDiagnostics("num_iter") == 10, "num_iter is wrong.")          self.assertTrue(so.getDiagnostics("num_iter") == 10, "num_iter is wrong.")
671          self.failUnless(so.getDiagnostics("num_inner_iter") == 20, "num_inner_iter is wrong.")          self.assertTrue(so.getDiagnostics("num_inner_iter") == 20, "num_inner_iter is wrong.")
672          self.failUnless(so.getDiagnostics("time") == 30, "time is wrong.")          self.assertTrue(so.getDiagnostics("time") == 30, "time is wrong.")
673          self.failUnless(so.getDiagnostics("set_up_time") == 40, "set_up_time is wrong.")          self.assertTrue(so.getDiagnostics("set_up_time") == 40, "set_up_time is wrong.")
674          self.failUnless(so.getDiagnostics("residual_norm") == 50, "residual_norm is wrong.")          self.assertTrue(so.getDiagnostics("residual_norm") == 50, "residual_norm is wrong.")
675          self.failUnless(not so.getDiagnostics("converged"), "converged is wrong.")          self.assertTrue(not so.getDiagnostics("converged"), "converged is wrong.")
676          self.failUnless(not so.hasConverged(), "convergence flag is wrong.")          self.assertTrue(not so.hasConverged(), "convergence flag is wrong.")
677          self.failUnless(so.getDiagnostics("cum_num_inner_iter") == 22, "cum_num_inner_iter is wrong.")          self.assertTrue(so.getDiagnostics("cum_num_inner_iter") == 22, "cum_num_inner_iter is wrong.")
678          self.failUnless(so.getDiagnostics("cum_num_iter") == 11, "cum_num_iter is wrong.")          self.assertTrue(so.getDiagnostics("cum_num_iter") == 11, "cum_num_iter is wrong.")
679          self.failUnless(so.getDiagnostics("cum_time") ==33, "cum_time is wrong.")          self.assertTrue(so.getDiagnostics("cum_time") ==33, "cum_time is wrong.")
680          self.failUnless(so.getDiagnostics("cum_set_up_time") == 44, "cum_set_up_time is wrong.")            self.assertTrue(so.getDiagnostics("cum_set_up_time") == 44, "cum_set_up_time is wrong.")  
681    
682          so.resetDiagnostics(all=True)          so.resetDiagnostics(all=True)
683          self.failUnless(so.getDiagnostics("num_iter") == None, "initial num_iter is wrong.")          self.assertTrue(so.getDiagnostics("num_iter") == None, "initial num_iter is wrong.")
684          self.failUnless(so.getDiagnostics("num_inner_iter") == None, "initial num_inner_iter is wrong.")          self.assertTrue(so.getDiagnostics("num_inner_iter") == None, "initial num_inner_iter is wrong.")
685          self.failUnless(so.getDiagnostics("time") == None, "initial time is wrong.")          self.assertTrue(so.getDiagnostics("time") == None, "initial time is wrong.")
686          self.failUnless(so.getDiagnostics("set_up_time") == None, "initial set_up_time is wrong.")          self.assertTrue(so.getDiagnostics("set_up_time") == None, "initial set_up_time is wrong.")
687          self.failUnless(so.getDiagnostics("residual_norm") == None, "initial residual_norm is wrong.")          self.assertTrue(so.getDiagnostics("residual_norm") == None, "initial residual_norm is wrong.")
688          self.failUnless(so.getDiagnostics("converged") == None, "initial converged is wrong.")          self.assertTrue(so.getDiagnostics("converged") == None, "initial converged is wrong.")
689          self.failUnless(so.hasConverged() == None, "initial convergence flag is wrong.")          self.assertTrue(so.hasConverged() == None, "initial convergence flag is wrong.")
690          self.failUnless(so.getDiagnostics("cum_num_inner_iter") == 0, "initial cum_num_inner_iter is wrong.")          self.assertTrue(so.getDiagnostics("cum_num_inner_iter") == 0, "initial cum_num_inner_iter is wrong.")
691          self.failUnless(so.getDiagnostics("cum_num_iter") == 0, "initial cum_num_iter is wrong.")          self.assertTrue(so.getDiagnostics("cum_num_iter") == 0, "initial cum_num_iter is wrong.")
692          self.failUnless(so.getDiagnostics("cum_time") ==0, "initial cum_time is wrong.")          self.assertTrue(so.getDiagnostics("cum_time") ==0, "initial cum_time is wrong.")
693          self.failUnless(so.getDiagnostics("cum_set_up_time") == 0, "initial cum_set_up_time is wrong.")          self.assertTrue(so.getDiagnostics("cum_set_up_time") == 0, "initial cum_set_up_time is wrong.")
694                    
695      def test_setCoefficient_WithIllegalFunctionSpace(self):      def test_setCoefficient_WithIllegalFunctionSpace(self):
696          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
697          self.failUnlessRaises(IllegalCoefficientFunctionSpace,mypde.setValue, C=Vector(0.,FunctionOnBoundary(self.domain)))          self.assertRaises(IllegalCoefficientFunctionSpace, mypde.setValue, C=Vector(0.,FunctionOnBoundary(self.domain)))
698    
699      def test_setCoefficient_WithWrongName(self):      def test_setCoefficient_WithWrongName(self):
700          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
701          self.failUnlessRaises(IllegalCoefficient, mypde.setValue, ROMA=0.)          self.assertRaises(IllegalCoefficient, mypde.setValue, ROMA=0.)
702    
703      def test_resetCoefficient_WithWrongShape(self):      def test_resetCoefficient_WithWrongShape(self):
704          mypde=LinearPDE(self.domain,numEquations=2,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numEquations=2,debug=self.DEBUG)
705          self.failUnlessRaises(IllegalCoefficientValue, mypde.setValue, C=0.)          self.assertRaises(IllegalCoefficientValue, mypde.setValue, C=0.)
706    
707      def test_reducedOn(self):      def test_reducedOn(self):
708          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
# Line 940  class Test_LinearPDE_noLumping(Test_line Line 710  class Test_LinearPDE_noLumping(Test_line
710          mypde.setReducedOrderOn()          mypde.setReducedOrderOn()
711          mypde.setValue(A=kronecker(self.domain),D=x[0],Y=x[0])          mypde.setValue(A=kronecker(self.domain),D=x[0],Y=x[0])
712          u=mypde.getSolution()          u=mypde.getSolution()
713          self.failUnless(self.check(u,1.),'solution is wrong.')          self.assertTrue(self.check(u,1.),'solution is wrong.')
714    
715      def test_attemptToChangeOrderAfterDefinedCoefficient(self):      def test_attemptToChangeOrderAfterDefinedCoefficient(self):
716          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
717          mypde.setValue(D=1.)          mypde.setValue(D=1.)
718          self.failUnlessRaises(RuntimeError,mypde.setReducedOrderOn)          self.assertRaises(RuntimeError, mypde.setReducedOrderOn)
719    
720      def test_reducedOnConfig(self):      def test_reducedOnConfig(self):
721          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
722          mypde.setReducedOrderOn()          mypde.setReducedOrderOn()
723          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.")
724      #      #
725      #  set coefficients for scalars:      #  set coefficients for scalars:
726      #      #
# Line 959  class Test_LinearPDE_noLumping(Test_line Line 729  class Test_LinearPDE_noLumping(Test_line
729          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
730          mypde.setValue(A=numpy.ones((d,d)))          mypde.setValue(A=numpy.ones((d,d)))
731          coeff=mypde.getCoefficient("A")          coeff=mypde.getCoefficient("A")
732          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))
733    
734            mypde.resetRightHandSideCoefficients()
735            self.assertFalse(mypde.getCoefficient("A").isEmpty(),"A is empty after reset of right hand side coefficients")
736    
737      def test_setCoefficient_B_Scalar(self):      def test_setCoefficient_B_Scalar(self):
738          d=self.domain.getDim()          d=self.domain.getDim()
739          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
740          mypde.setValue(B=numpy.ones((d,)))          mypde.setValue(B=numpy.ones((d,)))
741          coeff=mypde.getCoefficient("B")          coeff=mypde.getCoefficient("B")
742          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))
743    
744            mypde.resetRightHandSideCoefficients()
745            self.assertFalse(mypde.getCoefficient("B").isEmpty(),"B is empty after reset of right hand side coefficients")
746    
747      def test_setCoefficient_C_Scalar(self):      def test_setCoefficient_C_Scalar(self):
748          d=self.domain.getDim()          d=self.domain.getDim()
749          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
750          mypde.setValue(C=numpy.ones((d,)))          mypde.setValue(C=numpy.ones((d,)))
751          coeff=mypde.getCoefficient("C")          coeff=mypde.getCoefficient("C")
752          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))
753    
754            mypde.resetRightHandSideCoefficients()
755            self.assertFalse(mypde.getCoefficient("C").isEmpty(),"C is empty after reset of right hand side coefficients")
756    
757      def test_setCoefficient_D_Scalar(self):      def test_setCoefficient_D_Scalar(self):
758          d=self.domain.getDim()          d=self.domain.getDim()
759          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
760          mypde.setValue(D=1.)          mypde.setValue(D=1.)
761          coeff=mypde.getCoefficient("D")          coeff=mypde.getCoefficient("D")
762          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))
763    
764            mypde.resetRightHandSideCoefficients()
765            self.assertFalse(mypde.getCoefficient("D").isEmpty(),"D is empty after reset of right hand side coefficients")
766    
767      def test_setCoefficient_X_Scalar(self):      def test_setCoefficient_X_Scalar(self):
768          d=self.domain.getDim()          d=self.domain.getDim()
769          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
770          mypde.setValue(X=numpy.ones((d,)))          mypde.setValue(X=numpy.ones((d,)))
771          coeff=mypde.getCoefficient("X")          coeff=mypde.getCoefficient("X")
772          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))
773    
774            mypde.resetRightHandSideCoefficients()
775            self.assertTrue(mypde.getCoefficient("X").isEmpty(),"X is not empty after reset of right hand side coefficients")
776    
777      def test_setCoefficient_Y_Scalar(self):      def test_setCoefficient_Y_Scalar(self):
778          d=self.domain.getDim()          d=self.domain.getDim()
779          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
780          mypde.setValue(Y=1.)          mypde.setValue(Y=1.)
781          coeff=mypde.getCoefficient("Y")          coeff=mypde.getCoefficient("Y")
782          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((),Function(self.domain),1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),Function(self.domain),1))
783    
784            mypde.resetRightHandSideCoefficients()
785            self.assertTrue(mypde.getCoefficient("Y").isEmpty(),"Y is not empty after reset of right hand side coefficients")
786    
787      def test_setCoefficient_y_Scalar(self):      def test_setCoefficient_y_Scalar(self):
788          d=self.domain.getDim()          d=self.domain.getDim()
789          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
790          mypde.setValue(y=1.)          mypde.setValue(y=1.)
791          coeff=mypde.getCoefficient("y")          coeff=mypde.getCoefficient("y")
792          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((),FunctionOnBoundary(self.domain),1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),FunctionOnBoundary(self.domain),1))
793    
794            mypde.resetRightHandSideCoefficients()
795            self.assertTrue(mypde.getCoefficient("y").isEmpty(),"y is not empty after reset of right hand side coefficients")
796    
797      def test_setCoefficient_d_Scalar(self):      def test_setCoefficient_d_Scalar(self):
798          d=self.domain.getDim()          d=self.domain.getDim()
799          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
800          mypde.setValue(d=1.)          mypde.setValue(d=1.)
801          coeff=mypde.getCoefficient("d")          coeff=mypde.getCoefficient("d")
802          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))
803    
804            mypde.resetRightHandSideCoefficients()
805            self.assertFalse(mypde.getCoefficient("d").isEmpty(),"d is empty after reset of right hand side coefficients")
806    
807      def test_setCoefficient_d_contact_Scalar(self):      def test_setCoefficient_d_contact_Scalar(self):
808          d=self.domain.getDim()          if self.domain.supportsContactElements():
809          mypde=LinearPDE(self.domain,debug=self.DEBUG)              d=self.domain.getDim()
810          mypde.setValue(d_contact=1.)              mypde=LinearPDE(self.domain,debug=self.DEBUG)
811          coeff=mypde.getCoefficient("d_contact")              mypde.setValue(d_contact=1.)
812          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((),FunctionOnContactZero(self.domain),1,1))              coeff=mypde.getCoefficient("d_contact")
813                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),FunctionOnContactZero(self.domain),1,1))
814    
815                mypde.resetRightHandSideCoefficients()
816                self.assertFalse(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is empty after reset of right hand side coefficients")
817    
818      def test_setCoefficient_y_contact_Scalar(self):      def test_setCoefficient_y_contact_Scalar(self):
819          d=self.domain.getDim()          d=self.domain.getDim()
820          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          if self.domain.supportsContactElements():
821          mypde.setValue(y_contact=1.)              mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
822          coeff=mypde.getCoefficient("y_contact")              mypde.setValue(y_contact=1.)
823          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((),FunctionOnContactZero(self.domain),1))              coeff=mypde.getCoefficient("y_contact")
824                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),FunctionOnContactZero(self.domain),1))
825    
826                mypde.resetRightHandSideCoefficients()
827                self.assertTrue(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty after reset of right hand side coefficients")
828    
829      def test_setCoefficient_A_reduced_Scalar(self):      def test_setCoefficient_A_reduced_Scalar(self):
830          d=self.domain.getDim()          d=self.domain.getDim()
831          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
832          mypde.setValue(A_reduced=numpy.ones((d,d)))          mypde.setValue(A_reduced=numpy.ones((d,d)))
833          coeff=mypde.getCoefficient("A_reduced")          coeff=mypde.getCoefficient("A_reduced")
834          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))
835    
836            mypde.resetRightHandSideCoefficients()
837            self.assertFalse(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is empty after reset of right hand side coefficients")
838    
839      def test_setCoefficient_B_reduced_Scalar(self):      def test_setCoefficient_B_reduced_Scalar(self):
840          d=self.domain.getDim()          d=self.domain.getDim()
841          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
842          mypde.setValue(B_reduced=numpy.ones((d,)))          mypde.setValue(B_reduced=numpy.ones((d,)))
843          coeff=mypde.getCoefficient("B_reduced")          coeff=mypde.getCoefficient("B_reduced")
844          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))
845    
846            mypde.resetRightHandSideCoefficients()
847            self.assertFalse(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is empty after reset of right hand side coefficients")
848    
849      def test_setCoefficient_C_reduced_Scalar(self):      def test_setCoefficient_C_reduced_Scalar(self):
850          d=self.domain.getDim()          d=self.domain.getDim()
851          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
852          mypde.setValue(C_reduced=numpy.ones((d,)))          mypde.setValue(C_reduced=numpy.ones((d,)))
853          coeff=mypde.getCoefficient("C_reduced")          coeff=mypde.getCoefficient("C_reduced")
854          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))
855    
856            mypde.resetRightHandSideCoefficients()
857            self.assertFalse(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is empty after reset of right hand side coefficients")
858    
859      def test_setCoefficient_D_reduced_Scalar(self):      def test_setCoefficient_D_reduced_Scalar(self):
860          d=self.domain.getDim()          d=self.domain.getDim()
861          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
862          mypde.setValue(D_reduced=1.)          mypde.setValue(D_reduced=1.)
863          coeff=mypde.getCoefficient("D_reduced")          coeff=mypde.getCoefficient("D_reduced")
864          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))
865    
866            mypde.resetRightHandSideCoefficients()
867            self.assertFalse(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is empty after reset of right hand side coefficients")
868    
869      def test_setCoefficient_X_reduced_Scalar(self):      def test_setCoefficient_X_reduced_Scalar(self):
870          d=self.domain.getDim()          d=self.domain.getDim()
871          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
872          mypde.setValue(X_reduced=numpy.ones((d,)))          mypde.setValue(X_reduced=numpy.ones((d,)))
873          coeff=mypde.getCoefficient("X_reduced")          coeff=mypde.getCoefficient("X_reduced")
874          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))
875    
876            mypde.resetRightHandSideCoefficients()
877            self.assertTrue(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty after reset of right hand side coefficients")
878    
879      def test_setCoefficient_Y_reduced_Scalar(self):      def test_setCoefficient_Y_reduced_Scalar(self):
880          d=self.domain.getDim()          d=self.domain.getDim()
881          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
882          mypde.setValue(Y_reduced=1.)          mypde.setValue(Y_reduced=1.)
883          coeff=mypde.getCoefficient("Y_reduced")          coeff=mypde.getCoefficient("Y_reduced")
884          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((),ReducedFunction(self.domain),1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),ReducedFunction(self.domain),1))
885    
886            mypde.resetRightHandSideCoefficients()
887            self.assertTrue(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty after reset of right hand side coefficients")
888    
889      def test_setCoefficient_y_reduced_Scalar(self):      def test_setCoefficient_y_reduced_Scalar(self):
890          d=self.domain.getDim()          d=self.domain.getDim()
891          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
892          mypde.setValue(y_reduced=1.)          mypde.setValue(y_reduced=1.)
893          coeff=mypde.getCoefficient("y_reduced")          coeff=mypde.getCoefficient("y_reduced")
894          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((),ReducedFunctionOnBoundary(self.domain),1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),ReducedFunctionOnBoundary(self.domain),1))
895    
896            mypde.resetRightHandSideCoefficients()
897            self.assertTrue(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty after reset of right hand side coefficients")
898    
899      def test_setCoefficient_d_reduced_Scalar(self):      def test_setCoefficient_d_reduced_Scalar(self):
900          d=self.domain.getDim()          d=self.domain.getDim()
901          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
902          mypde.setValue(d_reduced=1.)          mypde.setValue(d_reduced=1.)
903          coeff=mypde.getCoefficient("d_reduced")          coeff=mypde.getCoefficient("d_reduced")
904          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))
905    
906            mypde.resetRightHandSideCoefficients()
907            self.assertFalse(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is not empty after reset of right hand side coefficients")
908    
909      def test_setCoefficient_d_contact_reduced_Scalar(self):      def test_setCoefficient_d_contact_reduced_Scalar(self):
910          d=self.domain.getDim()          if self.domain.supportsContactElements():
911          mypde=LinearPDE(self.domain,debug=self.DEBUG)              d=self.domain.getDim()
912          mypde.setValue(d_contact_reduced=1.)              mypde=LinearPDE(self.domain,debug=self.DEBUG)
913          coeff=mypde.getCoefficient("d_contact_reduced")              mypde.setValue(d_contact_reduced=1.)
914          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((),ReducedFunctionOnContactZero(self.domain),1,1))              coeff=mypde.getCoefficient("d_contact_reduced")
915                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),ReducedFunctionOnContactZero(self.domain),1,1))
916    
917                mypde.resetRightHandSideCoefficients()
918                self.assertFalse(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is not empty after reset of right hand side coefficients")
919      def test_setCoefficient_y_contact_reduced_Scalar(self):      def test_setCoefficient_y_contact_reduced_Scalar(self):
920          d=self.domain.getDim()          if self.domain.supportsContactElements():
921          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)              d=self.domain.getDim()
922          mypde.setValue(y_contact_reduced=1.)              mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
923          coeff=mypde.getCoefficient("y_contact_reduced")              mypde.setValue(y_contact_reduced=1.)
924          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((),ReducedFunctionOnContactZero(self.domain),1))              coeff=mypde.getCoefficient("y_contact_reduced")
925                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),ReducedFunctionOnContactZero(self.domain),1))
926    
927                mypde.resetRightHandSideCoefficients()
928                self.assertTrue(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty after reset of right hand side coefficients")
929    
930      def test_setCoefficient_r_Scalar(self):      def test_setCoefficient_r_Scalar(self):
931          d=self.domain.getDim()          d=self.domain.getDim()
932          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
933          mypde.setValue(r=1.)          mypde.setValue(r=1.)
934          coeff=mypde.getCoefficient("r")          coeff=mypde.getCoefficient("r")
935          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions()),((),Solution(self.domain),1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((),Solution(self.domain),1))
936    
937            mypde.resetRightHandSideCoefficients()
938            self.assertTrue(mypde.getCoefficient("r").isEmpty(),"r is not empty after reset of right hand side coefficients")
939    
940      def test_setCoefficient_q_Scalar(self):      def test_setCoefficient_q_Scalar(self):
941          d=self.domain.getDim()          d=self.domain.getDim()
942          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
943          mypde.setValue(q=1.)          mypde.setValue(q=1.)
944          coeff=mypde.getCoefficient("q")          coeff=mypde.getCoefficient("q")
945          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions()),((),Solution(self.domain),1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((),Solution(self.domain),1))
946    
947            mypde.resetRightHandSideCoefficients()
948            self.assertFalse(mypde.getCoefficient("q").isEmpty(),"q is empty after reset of right hand side coefficients")
949    
950      def test_setCoefficient_r_Scalar_reducedOn(self):      def test_setCoefficient_r_Scalar_reducedOn(self):
951          d=self.domain.getDim()          d=self.domain.getDim()
952          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
953          mypde.setReducedOrderOn()          mypde.setReducedOrderOn()
954          mypde.setValue(r=1.)          mypde.setValue(r=1.)
955          coeff=mypde.getCoefficient("r")          coeff=mypde.getCoefficient("r")
956          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions()),((),ReducedSolution(self.domain),1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((),ReducedSolution(self.domain),1))
957    
958            mypde.resetRightHandSideCoefficients()
959            self.assertTrue(mypde.getCoefficient("r").isEmpty(),"r is not empty after reset of right hand side coefficients")
960    
961      def test_setCoefficient_q_Scalar_reducedOn(self):      def test_setCoefficient_q_Scalar_reducedOn(self):
962          d=self.domain.getDim()          d=self.domain.getDim()
963          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
964          mypde.setReducedOrderOn()          mypde.setReducedOrderOn()
965          mypde.setValue(q=1.)          mypde.setValue(q=1.)
966          coeff=mypde.getCoefficient("q")          coeff=mypde.getCoefficient("q")
967          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions()),((),ReducedSolution(self.domain),1))          self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions()),((),ReducedSolution(self.domain),1))
968    
969            mypde.resetRightHandSideCoefficients()
970            self.assertFalse(mypde.getCoefficient("q").isEmpty(),"q is empty after reset of right hand side coefficients")
971    
972      def test_setCoefficient_A_reduced_Scalar_usingA(self):      def test_setCoefficient_A_reduced_Scalar_usingA(self):
973          d=self.domain.getDim()          d=self.domain.getDim()
974          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
975          mypde.setValue(A=Data(numpy.ones((d,d)),ReducedFunction(self.domain)))          mypde.setValue(A=Data(numpy.ones((d,d)),ReducedFunction(self.domain)))
976          coeff=mypde.getCoefficient("A_reduced")          if self.specialInterpolationSupported():
977          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((d,d),ReducedFunction(self.domain),1,1))              coeff_name='A'
978                FS=Function
979            else:
980                coeff_name='A_reduced'
981                FS=ReducedFunction
982            coeff=mypde.getCoefficient(coeff_name)
983            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,d),FS(self.domain),1,1))
984            mypde.resetRightHandSideCoefficients()
985            self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
986    
987      def test_setCoefficient_B_reduced_Scalar_usingB(self):      def test_setCoefficient_B_reduced_Scalar_usingB(self):
988          d=self.domain.getDim()          d=self.domain.getDim()
989          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
990          mypde.setValue(B=Data(numpy.ones((d,)),ReducedFunction(self.domain)))          mypde.setValue(B=Data(numpy.ones((d,)),ReducedFunction(self.domain)))
991          coeff=mypde.getCoefficient("B_reduced")          if self.specialInterpolationSupported():
992          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1,1))              coeff_name='B'
993                FS=Function
994            else:
995                coeff_name='B_reduced'
996                FS=ReducedFunction
997            coeff=mypde.getCoefficient(coeff_name)
998            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,),FS(self.domain),1,1))
999            mypde.resetRightHandSideCoefficients()
1000            self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
1001    
1002      def test_setCoefficient_C_reduced_Scalar_usingC(self):      def test_setCoefficient_C_reduced_Scalar_usingC(self):
1003          d=self.domain.getDim()          d=self.domain.getDim()
1004          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1005          mypde.setValue(C=Data(numpy.ones((d,)),ReducedFunction(self.domain)))          mypde.setValue(C=Data(numpy.ones((d,)),ReducedFunction(self.domain)))
1006          coeff=mypde.getCoefficient("C_reduced")          if self.specialInterpolationSupported():
1007          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((d,),ReducedFunction(self.domain),1,1))              coeff_name='C'
1008                FS=Function
1009            else:
1010                coeff_name='C_reduced'
1011                FS=ReducedFunction
1012            coeff=mypde.getCoefficient(coeff_name)
1013            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((d,),FS(self.domain),1,1))
1014            mypde.resetRightHandSideCoefficients()
1015            self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
1016    
1017      def test_setCoefficient_D_reduced_Scalar_usingD(self):      def test_setCoefficient_D_reduced_Scalar_usingD(self):
1018          d=self.domain.getDim()          d=self.domain.getDim()
1019          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1020          mypde.setValue(D=Scalar(1.,ReducedFunction(self.domain)))          mypde.setValue(D=Scalar(1.,ReducedFunction(self.domain)))
1021          coeff=mypde.getCoefficient("D_reduced")          if self.specialInterpolationSupported():
1022          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((),ReducedFunction(self.domain),1,1))              coeff_name='D'
1023                FS=Function
1024            else:
1025                coeff_name='D_reduced'
1026                FS=ReducedFunction
1027            coeff=mypde.getCoefficient(coeff_name)
1028            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),FS(self.domain),1,1))
1029            mypde.resetRightHandSideCoefficients()
1030            self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
1031    
1032      def test_setCoefficient_X_reduced_Scalar_usingX(self):      def test_setCoefficient_X_reduced_Scalar_usingX(self):
1033          d=self.domain.getDim()          d=self.domain.getDim()
1034          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1035          mypde.setValue(X_reduced=Data(numpy.ones((d,)),ReducedFunction(self.domain)))          mypde.setValue(X_reduced=Data(numpy.ones((d,)),ReducedFunction(self.domain)))
1036          coeff=mypde.getCoefficient("X_reduced")          coeff=mypde.getCoefficient('X_reduced')
1037          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))
1038            mypde.resetRightHandSideCoefficients()
1039            self.assertTrue(mypde.getCoefficient('X_reduced').isEmpty(),"X_reduced is not empty after reset of right hand side coefficients")
1040    
1041      def test_setCoefficient_Y_reduced_Scalar_usingY(self):      def test_setCoefficient_Y_reduced_Scalar_usingY(self):
1042          d=self.domain.getDim()          d=self.domain.getDim()
1043          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1044          mypde.setValue(Y=Scalar(1.,ReducedFunction(self.domain)))          mypde.setValue(Y=Scalar(1.,ReducedFunction(self.domain)))
1045          coeff=mypde.getCoefficient("Y_reduced")          if self.specialInterpolationSupported():
1046          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((),ReducedFunction(self.domain),1))              coeff_name='Y'
1047                FS=Function
1048            else:
1049                coeff_name='Y_reduced'
1050                FS=ReducedFunction
1051            coeff=mypde.getCoefficient(coeff_name)
1052            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),FS(self.domain),1))
1053            mypde.resetRightHandSideCoefficients()
1054            self.assertTrue(mypde.getCoefficient(coeff_name).isEmpty(),"%s is not empty after reset of right hand side coefficients"%coeff_name)
1055    
1056      def test_setCoefficient_y_reduced_Scalar_using_y(self):      def test_setCoefficient_y_reduced_Scalar_using_y(self):
1057          d=self.domain.getDim()          d=self.domain.getDim()
1058          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1059          mypde.setValue(y=Scalar(1.,ReducedFunctionOnBoundary(self.domain)))          mypde.setValue(y=Scalar(1.,ReducedFunctionOnBoundary(self.domain)))
1060          coeff=mypde.getCoefficient("y_reduced")          if self.specialInterpolationSupported():
1061          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((),ReducedFunctionOnBoundary(self.domain),1))              coeff_name='y'
1062                FS=FunctionOnBoundary
1063            else:
1064                coeff_name='y_reduced'
1065                FS=ReducedFunctionOnBoundary
1066            coeff=mypde.getCoefficient(coeff_name)
1067            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),FS(self.domain),1))
1068            mypde.resetRightHandSideCoefficients()
1069            self.assertTrue(mypde.getCoefficient(coeff_name).isEmpty(),"%s is not empty after reset of right hand side coefficients"%coeff_name)
1070    
1071      def test_setCoefficient_d_reduced_Scalar_using_d(self):      def test_setCoefficient_d_reduced_Scalar_using_d(self):
1072          d=self.domain.getDim()          d=self.domain.getDim()
1073          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1074          mypde.setValue(d=Scalar(1.,ReducedFunctionOnBoundary(self.domain)))          mypde.setValue(d=Scalar(1.,ReducedFunctionOnBoundary(self.domain)))
1075          coeff=mypde.getCoefficient("d_reduced")          if self.specialInterpolationSupported():
1076          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((),ReducedFunctionOnBoundary(self.domain),1,1))              coeff_name='d'
1077                FS=FunctionOnBoundary
1078            else:
1079                coeff_name='d_reduced'
1080                FS=ReducedFunctionOnBoundary
1081            coeff=mypde.getCoefficient(coeff_name)
1082            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),FS(self.domain),1,1))
1083            mypde.resetRightHandSideCoefficients()
1084            self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
1085    
1086      def test_setCoefficient_d_contact_reduced_Scalar_using_d_contact(self):      def test_setCoefficient_d_contact_reduced_Scalar_using_d_contact(self):
1087          d=self.domain.getDim()          if self.domain.supportsContactElements():
1088          mypde=LinearPDE(self.domain,debug=self.DEBUG)              d=self.domain.getDim()
1089          mypde.setValue(d_contact=Scalar(1.,ReducedFunctionOnContactZero(self.domain)))              mypde=LinearPDE(self.domain,debug=self.DEBUG)
1090          coeff=mypde.getCoefficient("d_contact_reduced")              mypde.setValue(d_contact=Scalar(1.,ReducedFunctionOnContactZero(self.domain)))
1091          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((),ReducedFunctionOnContactZero(self.domain),1,1))              if self.specialInterpolationSupported():
1092                    coeff_name='d_contact'
1093                    FS=FunctionOnContactZero
1094                else:
1095                    coeff_name='d_contact_reduced'
1096                    FS=ReducedFunctionOnContactZero
1097                coeff=mypde.getCoefficient(coeff_name)
1098                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((),FS(self.domain),1,1))
1099                mypde.resetRightHandSideCoefficients()
1100                self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
1101    
1102      def test_setCoefficient_y_contact_reduced_Scalar_using_y_contact(self):      def test_setCoefficient_y_contact_reduced_Scalar_using_y_contact(self):
1103          d=self.domain.getDim()          if self.domain.supportsContactElements():
1104          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)              d=self.domain.getDim()
1105          mypde.setValue(y_contact=Scalar(1.,ReducedFunctionOnContactZero(self.domain)))              mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1106          coeff=mypde.getCoefficient("y_contact_reduced")              mypde.setValue(y_contact=Scalar(1.,ReducedFunctionOnContactZero(self.domain)))
1107          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((),ReducedFunctionOnContactZero(self.domain),1))              if self.specialInterpolationSupported():
1108                    coeff_name='y_contact'
1109                    FS=FunctionOnContactZero
1110                else:
1111                    coeff_name='y_contact_reduced'
1112                    FS=ReducedFunctionOnContactZero
1113                coeff=mypde.getCoefficient(coeff_name)
1114                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((),FS(self.domain),1))
1115                mypde.resetRightHandSideCoefficients()
1116                self.assertTrue(mypde.getCoefficient(coeff_name).isEmpty(),"%s is not empty after reset of right hand side coefficients"%coeff_name)
1117    
1118      #      #
1119      #  set coefficients for systems:      #  set coefficients for systems:
1120      #      #
# Line 1169  class Test_LinearPDE_noLumping(Test_line Line 1123  class Test_LinearPDE_noLumping(Test_line
1123          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1124          mypde.setValue(A=numpy.ones((self.N,d,self.N,d)))          mypde.setValue(A=numpy.ones((self.N,d,self.N,d)))
1125          coeff=mypde.getCoefficient("A")          coeff=mypde.getCoefficient("A")
1126          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))
1127            mypde.resetRightHandSideCoefficients()
1128            self.assertFalse(mypde.getCoefficient("A").isEmpty(),"A is empty after reset of right hand side coefficients")
1129            
1130      def test_setCoefficient_B_System(self):      def test_setCoefficient_B_System(self):
1131          d=self.domain.getDim()          d=self.domain.getDim()
1132          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1133          mypde.setValue(B=numpy.ones((self.N,d,self.N)))          mypde.setValue(B=numpy.ones((self.N,d,self.N)))
1134          coeff=mypde.getCoefficient("B")          coeff=mypde.getCoefficient("B")
1135          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))
1136            mypde.resetRightHandSideCoefficients()
1137            self.assertFalse(mypde.getCoefficient("B").isEmpty(),"B is empty after reset of right hand side coefficients")
1138      def test_setCoefficient_C_System(self):      def test_setCoefficient_C_System(self):
1139          d=self.domain.getDim()          d=self.domain.getDim()
1140          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1141          mypde.setValue(C=numpy.ones((self.N,self.N,d)))          mypde.setValue(C=numpy.ones((self.N,self.N,d)))
1142          coeff=mypde.getCoefficient("C")          coeff=mypde.getCoefficient("C")
1143          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))
1144            mypde.resetRightHandSideCoefficients()
1145            self.assertFalse(mypde.getCoefficient("C").isEmpty(),"C is empty after reset of right hand side coefficients")
1146      def test_setCoefficient_D_System(self):      def test_setCoefficient_D_System(self):
1147          d=self.domain.getDim()          d=self.domain.getDim()
1148          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1149          mypde.setValue(D=numpy.ones((self.N,self.N)))          mypde.setValue(D=numpy.ones((self.N,self.N)))
1150          coeff=mypde.getCoefficient("D")          coeff=mypde.getCoefficient("D")
1151          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))
1152            mypde.resetRightHandSideCoefficients()
1153            self.assertFalse(mypde.getCoefficient("D").isEmpty(),"D is empty after reset of right hand side coefficients")
1154      def test_setCoefficient_X_System(self):      def test_setCoefficient_X_System(self):
1155          d=self.domain.getDim()          d=self.domain.getDim()
1156          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1157          mypde.setValue(X=numpy.ones((self.N,d)))          mypde.setValue(X=numpy.ones((self.N,d)))
1158          coeff=mypde.getCoefficient("X")          coeff=mypde.getCoefficient("X")
1159          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))
1160            mypde.resetRightHandSideCoefficients()
1161            self.assertTrue(mypde.getCoefficient("X").isEmpty(),"X is not empty after reset of right hand side coefficients")
1162      def test_setCoefficient_Y_System(self):      def test_setCoefficient_Y_System(self):
1163          d=self.domain.getDim()          d=self.domain.getDim()
1164          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1165          mypde.setValue(Y=numpy.ones((self.N,)))          mypde.setValue(Y=numpy.ones((self.N,)))
1166          coeff=mypde.getCoefficient("Y")          coeff=mypde.getCoefficient("Y")
1167          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))
1168            mypde.resetRightHandSideCoefficients()
1169            self.assertTrue(mypde.getCoefficient("Y").isEmpty(),"Y is not empty after reset of right hand side coefficients")
1170      def test_setCoefficient_y_System(self):      def test_setCoefficient_y_System(self):
1171          d=self.domain.getDim()          d=self.domain.getDim()
1172          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1173          mypde.setValue(y=numpy.ones((self.N,)))          mypde.setValue(y=numpy.ones((self.N,)))
1174          coeff=mypde.getCoefficient("y")          coeff=mypde.getCoefficient("y")
1175          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))
1176            mypde.resetRightHandSideCoefficients()
1177            self.assertTrue(mypde.getCoefficient("y").isEmpty(),"y is not empty after reset of right hand side coefficients")
1178      def test_setCoefficient_d_System(self):      def test_setCoefficient_d_System(self):
1179          d=self.domain.getDim()          d=self.domain.getDim()
1180          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1181          mypde.setValue(d=numpy.ones((self.N,self.N)))          mypde.setValue(d=numpy.ones((self.N,self.N)))
1182          coeff=mypde.getCoefficient("d")          coeff=mypde.getCoefficient("d")
1183          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))
1184            mypde.resetRightHandSideCoefficients()
1185            self.assertFalse(mypde.getCoefficient("d").isEmpty(),"d is empty after reset of right hand side coefficients")
1186      def test_setCoefficient_d_contact_System(self):      def test_setCoefficient_d_contact_System(self):
1187          d=self.domain.getDim()          if self.domain.supportsContactElements():
1188          mypde=LinearPDE(self.domain,debug=self.DEBUG)              d=self.domain.getDim()
1189          mypde.setValue(d_contact=numpy.ones((self.N,self.N)))              mypde=LinearPDE(self.domain,debug=self.DEBUG)
1190          coeff=mypde.getCoefficient("d_contact")              mypde.setValue(d_contact=numpy.ones((self.N,self.N)))
1191          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")
1192                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),FunctionOnContactZero(self.domain),self.N,self.N))
1193                mypde.resetRightHandSideCoefficients()
1194                self.assertFalse(mypde.getCoefficient("d_contact").isEmpty(),"d_contact is empty after reset of right hand side coefficients")
1195      def test_setCoefficient_y_contact_System(self):      def test_setCoefficient_y_contact_System(self):
1196          d=self.domain.getDim()          if self.domain.supportsContactElements():
1197          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)              d=self.domain.getDim()
1198          mypde.setValue(y_contact=numpy.ones((self.N,)))              mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1199          coeff=mypde.getCoefficient("y_contact")              mypde.setValue(y_contact=numpy.ones((self.N,)))
1200          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((self.N,),FunctionOnContactZero(self.domain),self.N))              coeff=mypde.getCoefficient("y_contact")
1201                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),FunctionOnContactZero(self.domain),self.N))
1202                mypde.resetRightHandSideCoefficients()
1203                self.assertTrue(mypde.getCoefficient("y_contact").isEmpty(),"y_contact is not empty after reset of right hand side coefficients")
1204      def test_setCoefficient_A_reduced_System(self):      def test_setCoefficient_A_reduced_System(self):
1205          d=self.domain.getDim()          d=self.domain.getDim()
1206          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1207          mypde.setValue(A_reduced=numpy.ones((self.N,d,self.N,d)))          mypde.setValue(A_reduced=numpy.ones((self.N,d,self.N,d)))
1208          coeff=mypde.getCoefficient("A_reduced")          coeff=mypde.getCoefficient("A_reduced")
1209          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))
1210            mypde.resetRightHandSideCoefficients()
1211            self.assertFalse(mypde.getCoefficient("A_reduced").isEmpty(),"A_reduced is empty after reset of right hand side coefficients")
1212      def test_setCoefficient_B_reduced_System(self):      def test_setCoefficient_B_reduced_System(self):
1213          d=self.domain.getDim()          d=self.domain.getDim()
1214          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1215          mypde.setValue(B_reduced=numpy.ones((self.N,d,self.N)))          mypde.setValue(B_reduced=numpy.ones((self.N,d,self.N)))
1216          coeff=mypde.getCoefficient("B_reduced")          coeff=mypde.getCoefficient("B_reduced")
1217          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))
1218            mypde.resetRightHandSideCoefficients()
1219            self.assertFalse(mypde.getCoefficient("B_reduced").isEmpty(),"B_reduced is empty after reset of right hand side coefficients")
1220      def test_setCoefficient_C_reduced_System(self):      def test_setCoefficient_C_reduced_System(self):
1221          d=self.domain.getDim()          d=self.domain.getDim()
1222          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1223          mypde.setValue(C_reduced=numpy.ones((self.N,self.N,d)))          mypde.setValue(C_reduced=numpy.ones((self.N,self.N,d)))
1224          coeff=mypde.getCoefficient("C_reduced")          coeff=mypde.getCoefficient("C_reduced")
1225          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))
1226            mypde.resetRightHandSideCoefficients()
1227            self.assertFalse(mypde.getCoefficient("C_reduced").isEmpty(),"C_reduced is empty after reset of right hand side coefficients")
1228      def test_setCoefficient_D_System_reduced(self):      def test_setCoefficient_D_System_reduced(self):
1229          d=self.domain.getDim()          d=self.domain.getDim()
1230          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1231          mypde.setValue(D_reduced=numpy.ones((self.N,self.N)))          mypde.setValue(D_reduced=numpy.ones((self.N,self.N)))
1232          coeff=mypde.getCoefficient("D_reduced")          coeff=mypde.getCoefficient("D_reduced")
1233          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))
1234            mypde.resetRightHandSideCoefficients()
1235            self.assertFalse(mypde.getCoefficient("D_reduced").isEmpty(),"D_reduced is empty after reset of right hand side coefficients")
1236      def test_setCoefficient_X_System_reduced(self):      def test_setCoefficient_X_System_reduced(self):
1237          d=self.domain.getDim()          d=self.domain.getDim()
1238          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1239          mypde.setValue(X_reduced=numpy.ones((self.N,d)))          mypde.setValue(X_reduced=numpy.ones((self.N,d)))
1240          coeff=mypde.getCoefficient("X_reduced")          coeff=mypde.getCoefficient("X_reduced")
1241          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))
1242            mypde.resetRightHandSideCoefficients()
1243            self.assertTrue(mypde.getCoefficient("X_reduced").isEmpty(),"X_reduced is not empty after reset of right hand side coefficients")
1244      def test_setCoefficient_Y_System_reduced(self):      def test_setCoefficient_Y_System_reduced(self):
1245          d=self.domain.getDim()          d=self.domain.getDim()
1246          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1247          mypde.setValue(Y_reduced=numpy.ones((self.N,)))          mypde.setValue(Y_reduced=numpy.ones((self.N,)))
1248          coeff=mypde.getCoefficient("Y_reduced")          coeff=mypde.getCoefficient("Y_reduced")
1249          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))
1250            mypde.resetRightHandSideCoefficients()
1251            self.assertTrue(mypde.getCoefficient("Y_reduced").isEmpty(),"Y_reduced is not empty after reset of right hand side coefficients")
1252      def test_setCoefficient_y_System_reduced(self):      def test_setCoefficient_y_System_reduced(self):
1253          d=self.domain.getDim()          d=self.domain.getDim()
1254          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1255          mypde.setValue(y_reduced=numpy.ones((self.N,)))          mypde.setValue(y_reduced=numpy.ones((self.N,)))
1256          coeff=mypde.getCoefficient("y_reduced")          coeff=mypde.getCoefficient("y_reduced")
1257          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))
1258            mypde.resetRightHandSideCoefficients()
1259            self.assertTrue(mypde.getCoefficient("y_reduced").isEmpty(),"y_reduced is not empty after reset of right hand side coefficients")
1260      def test_setCoefficient_d_reduced_System(self):      def test_setCoefficient_d_reduced_System(self):
1261          d=self.domain.getDim()          d=self.domain.getDim()
1262          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1263          mypde.setValue(d_reduced=numpy.ones((self.N,self.N)))          mypde.setValue(d_reduced=numpy.ones((self.N,self.N)))
1264          coeff=mypde.getCoefficient("d_reduced")          coeff=mypde.getCoefficient("d_reduced")
1265          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))
1266            mypde.resetRightHandSideCoefficients()
1267            self.assertFalse(mypde.getCoefficient("d_reduced").isEmpty(),"d_reduced is empty after reset of right hand side coefficients")
1268      def test_setCoefficient_d_contact_reduced_System(self):      def test_setCoefficient_d_contact_reduced_System(self):
1269          d=self.domain.getDim()          if self.domain.supportsContactElements():
1270          mypde=LinearPDE(self.domain,debug=self.DEBUG)              d=self.domain.getDim()
1271          mypde.setValue(d_contact_reduced=numpy.ones((self.N,self.N)))              mypde=LinearPDE(self.domain,debug=self.DEBUG)
1272          coeff=mypde.getCoefficient("d_contact_reduced")              mypde.setValue(d_contact_reduced=numpy.ones((self.N,self.N)))
1273          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")
1274                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),ReducedFunctionOnContactZero(self.domain),self.N,self.N))
1275                mypde.resetRightHandSideCoefficients()
1276                self.assertFalse(mypde.getCoefficient("d_contact_reduced").isEmpty(),"d_contact_reduced is empty after reset of right hand side coefficients")
1277      def test_setCoefficient_y_contact_reduced_System(self):      def test_setCoefficient_y_contact_reduced_System(self):
1278          d=self.domain.getDim()          if self.domain.supportsContactElements():
1279          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)              d=self.domain.getDim()
1280          mypde.setValue(y_contact_reduced=numpy.ones((self.N,)))              mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1281          coeff=mypde.getCoefficient("y_contact_reduced")              mypde.setValue(y_contact_reduced=numpy.ones((self.N,)))
1282          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((self.N,),ReducedFunctionOnContactZero(self.domain),self.N))              coeff=mypde.getCoefficient("y_contact_reduced")
1283                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),ReducedFunctionOnContactZero(self.domain),self.N))
1284                mypde.resetRightHandSideCoefficients()
1285                self.assertTrue(mypde.getCoefficient("y_contact_reduced").isEmpty(),"y_contact_reduced is not empty after reset of right hand side coefficients")
1286      def test_setCoefficient_r_System(self):      def test_setCoefficient_r_System(self):
1287          d=self.domain.getDim()          d=self.domain.getDim()
1288          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
1289          mypde.setValue(r=numpy.ones((self.N,)))          mypde.setValue(r=numpy.ones((self.N,)))
1290          coeff=mypde.getCoefficient("r")          coeff=mypde.getCoefficient("r")
1291          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))
1292            mypde.resetRightHandSideCoefficients()
1293            self.assertTrue(mypde.getCoefficient("r").isEmpty(),"r is not empty after reset of right hand side coefficients")
1294      def test_setCoefficient_q_System(self):      def test_setCoefficient_q_System(self):
1295          d=self.domain.getDim()          d=self.domain.getDim()
1296          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
1297          mypde.setValue(q=numpy.ones((self.N,)))          mypde.setValue(q=numpy.ones((self.N,)))
1298          coeff=mypde.getCoefficient("q")          coeff=mypde.getCoefficient("q")
1299          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))
1300            mypde.resetRightHandSideCoefficients()
1301            self.assertFalse(mypde.getCoefficient("q").isEmpty(),"q is empty after reset of right hand side coefficients")
1302      def test_setCoefficient_r_System_reducedOn(self):      def test_setCoefficient_r_System_reducedOn(self):
1303          d=self.domain.getDim()          d=self.domain.getDim()
1304          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
1305          mypde.setReducedOrderOn()          mypde.setReducedOrderOn()
1306          mypde.setValue(r=numpy.ones((self.N,)))          mypde.setValue(r=numpy.ones((self.N,)))
1307          coeff=mypde.getCoefficient("r")          coeff=mypde.getCoefficient("r")
1308          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))
1309            mypde.resetRightHandSideCoefficients()
1310            self.assertTrue(mypde.getCoefficient("r").isEmpty(),"r is not empty after reset of right hand side coefficients")
1311      def test_setCoefficient_q_System_reducedOn(self):      def test_setCoefficient_q_System_reducedOn(self):
1312          d=self.domain.getDim()          d=self.domain.getDim()
1313          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numEquations=3,debug=self.DEBUG)
1314          mypde.setReducedOrderOn()          mypde.setReducedOrderOn()
1315          mypde.setValue(q=numpy.ones((self.N,)))          mypde.setValue(q=numpy.ones((self.N,)))
1316          coeff=mypde.getCoefficient("q")          coeff=mypde.getCoefficient("q")
1317          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))
1318            mypde.resetRightHandSideCoefficients()
1319            self.assertFalse(mypde.getCoefficient("q").isEmpty(),"q is empty after reset of right hand side coefficients")
1320    
1321      def test_setCoefficient_A_reduced_System_using_A(self):      def test_setCoefficient_A_reduced_System_using_A(self):
1322          d=self.domain.getDim()          d=self.domain.getDim()
1323          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1324          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)))
1325          coeff=mypde.getCoefficient("A_reduced")          if self.specialInterpolationSupported():
1326          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'
1327                FS=Function
1328            else:
1329                coeff_name='A_reduced'
1330                FS=ReducedFunction
1331            coeff=mypde.getCoefficient(coeff_name)
1332            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,d,self.N,d),FS(self.domain),self.N,self.N))
1333            mypde.resetRightHandSideCoefficients()
1334            self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
1335    
1336      def test_setCoefficient_B_reduced_System_using_B(self):      def test_setCoefficient_B_reduced_System_using_B(self):
1337          d=self.domain.getDim()          d=self.domain.getDim()
1338          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1339          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)))
1340          coeff=mypde.getCoefficient("B_reduced")          if self.specialInterpolationSupported():
1341          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'
1342                FS=Function
1343            else:
1344                coeff_name='B_reduced'
1345                FS=ReducedFunction
1346            coeff=mypde.getCoefficient(coeff_name)
1347            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,d,self.N),FS(self.domain),self.N,self.N))
1348            mypde.resetRightHandSideCoefficients()
1349            self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
1350    
1351      def test_setCoefficient_C_reduced_System_using_C(self):      def test_setCoefficient_C_reduced_System_using_C(self):
1352          d=self.domain.getDim()          d=self.domain.getDim()
1353          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1354          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)))
1355          coeff=mypde.getCoefficient("C_reduced")          if self.specialInterpolationSupported():
1356          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'
1357      def test_setCoefficient_D_System_reduced_using_D(self):              FS=Function
1358            else:
1359                coeff_name='C_reduced'
1360                FS=ReducedFunction
1361            coeff=mypde.getCoefficient(coeff_name)
1362            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N,d),FS(self.domain),self.N,self.N))
1363            mypde.resetRightHandSideCoefficients()
1364            self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
1365    
1366        def test_setCoefficient_D_reduced_System_using_D(self):
1367          d=self.domain.getDim()          d=self.domain.getDim()
1368          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1369          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)))
1370          coeff=mypde.getCoefficient("D_reduced")          if self.specialInterpolationSupported():
1371          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,self.N),ReducedFunction(self.domain),self.N,self.N))              coeff_name='D'
1372      def test_setCoefficient_X_System_reduced_using_X(self):              FS=Function
1373            else:
1374                coeff_name='D_reduced'
1375                FS=ReducedFunction
1376            coeff=mypde.getCoefficient(coeff_name)
1377            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(), mypde.getNumEquations()),((self.N,self.N),FS(self.domain),self.N,self.N))
1378            mypde.resetRightHandSideCoefficients()
1379            self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
1380    
1381        def test_setCoefficient_X_reduced_System_using_X(self):
1382          d=self.domain.getDim()          d=self.domain.getDim()
1383          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1384          mypde.setValue(X=Data(numpy.ones((self.N,d)),ReducedFunction(self.domain)))          mypde.setValue(X=Data(numpy.ones((self.N,d)),ReducedFunction(self.domain)))
1385          coeff=mypde.getCoefficient("X_reduced")          if self.specialInterpolationSupported():
1386          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((self.N,d),ReducedFunction(self.domain),self.N))              coeff_name='X'
1387      def test_setCoefficient_Y_System_reduced_using_Y(self):              FS=Function
1388            else:
1389                coeff_name='X_reduced'
1390                FS=ReducedFunction
1391            coeff=mypde.getCoefficient(coeff_name)
1392            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,d),FS(self.domain),self.N))
1393            mypde.resetRightHandSideCoefficients()
1394            self.assertTrue(mypde.getCoefficient(coeff_name).isEmpty(),"%s is not empty after reset of right hand side coefficients"%coeff_name)
1395    
1396        def test_setCoefficient_Y_reduced_System_using_Y(self):
1397          d=self.domain.getDim()          d=self.domain.getDim()
1398          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1399          mypde.setValue(Y=Data(numpy.ones((self.N,)),ReducedFunction(self.domain)))          mypde.setValue(Y=Data(numpy.ones((self.N,)),ReducedFunction(self.domain)))
1400          coeff=mypde.getCoefficient("Y_reduced")          if self.specialInterpolationSupported():
1401          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((self.N,),ReducedFunction(self.domain),self.N))              coeff_name='Y'
1402                FS=Function
1403            else:
1404                coeff_name='Y_reduced'
1405                FS=ReducedFunction
1406            coeff=mypde.getCoefficient(coeff_name)
1407            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),FS(self.domain),self.N))
1408            mypde.resetRightHandSideCoefficients()
1409            self.assertTrue(mypde.getCoefficient(coeff_name).isEmpty(),"%s is not empty after reset of right hand side coefficients"%coeff_name)
1410    
1411      def test_setCoefficient_y_reduced_System_using_y(self):      def test_setCoefficient_y_reduced_System_using_y(self):
1412          d=self.domain.getDim()          d=self.domain.getDim()
1413          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1414          mypde.setValue(y=Data(numpy.ones((self.N,)),ReducedFunctionOnBoundary(self.domain)))          mypde.setValue(y=Data(numpy.ones((self.N,)),ReducedFunctionOnBoundary(self.domain)))
1415          coeff=mypde.getCoefficient("y_reduced")          if self.specialInterpolationSupported():
1416          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((self.N,),ReducedFunctionOnBoundary(self.domain),self.N))              coeff_name='y'
1417                FS=FunctionOnBoundary
1418            else:
1419                coeff_name='y_reduced'
1420                FS=ReducedFunctionOnBoundary
1421            coeff=mypde.getCoefficient(coeff_name)
1422            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),FS(self.domain),self.N))
1423            mypde.resetRightHandSideCoefficients()
1424            self.assertTrue(mypde.getCoefficient(coeff_name).isEmpty(),"%s is not empty after reset of right hand side coefficients"%coeff_name)
1425    
1426      def test_setCoefficient_d_reduced_System_using_d(self):      def test_setCoefficient_d_reduced_System_using_d(self):
1427          d=self.domain.getDim()          d=self.domain.getDim()
1428          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1429          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)))
1430          coeff=mypde.getCoefficient("d_reduced")          if self.specialInterpolationSupported():
1431          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,self.N),ReducedFunctionOnBoundary(self.domain),self.N,self.N))              coeff_name='d'
1432                FS=FunctionOnBoundary
1433            else:
1434                coeff_name='d_reduced'
1435                FS=ReducedFunctionOnBoundary
1436            coeff=mypde.getCoefficient(coeff_name)
1437            self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,self.N),FS(self.domain),self.N,self.N))
1438            mypde.resetRightHandSideCoefficients()
1439            self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
1440    
1441      def test_setCoefficient_d_contact_reduced_System_using_d_contact(self):      def test_setCoefficient_d_contact_reduced_System_using_d_contact(self):
1442          d=self.domain.getDim()          if self.domain.supportsContactElements():
1443          mypde=LinearPDE(self.domain,debug=self.DEBUG)              d=self.domain.getDim()
1444          mypde.setValue(d_contact=Data(numpy.ones((self.N,self.N)),ReducedFunctionOnContactZero(self.domain)))              mypde=LinearPDE(self.domain,debug=self.DEBUG)
1445          coeff=mypde.getCoefficient("d_contact_reduced")              mypde.setValue(d_contact=Data(numpy.ones((self.N,self.N)),ReducedFunctionOnContactZero(self.domain)))
1446          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,self.N),ReducedFunctionOnContactZero(self.domain),self.N,self.N))              if self.specialInterpolationSupported():
1447                    coeff_name='d_contact'
1448                    FS=FunctionOnContactZero
1449                else:
1450                    coeff_name='d_contact_reduced'
1451                    FS=ReducedFunctionOnContactZero
1452                coeff=mypde.getCoefficient(coeff_name)
1453                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumSolutions(),mypde.getNumEquations()),((self.N,self.N),FS(self.domain),self.N,self.N))
1454                mypde.resetRightHandSideCoefficients()
1455                self.assertFalse(mypde.getCoefficient(coeff_name).isEmpty(),"%s is empty after reset of right hand side coefficients"%coeff_name)
1456    
1457      def test_setCoefficient_y_contact_reduced_System_using_y_contact(self):      def test_setCoefficient_y_contact_reduced_System_using_y_contact(self):
1458          d=self.domain.getDim()          if self.domain.supportsContactElements():
1459          mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)              d=self.domain.getDim()
1460          mypde.setValue(y_contact=Data(numpy.ones((self.N,)),ReducedFunctionOnContactZero(self.domain)))              mypde=LinearPDE(self.domain,numSolutions=3,debug=self.DEBUG)
1461          coeff=mypde.getCoefficient("y_contact_reduced")              mypde.setValue(y_contact=Data(numpy.ones((self.N,)),ReducedFunctionOnContactZero(self.domain)))
1462          self.failUnlessEqual((coeff.getShape(),coeff.getFunctionSpace(),mypde.getNumEquations()),((self.N,),ReducedFunctionOnContactZero(self.domain),self.N))              if self.specialInterpolationSupported():
1463                    coeff_name='y_contact'
1464                    FS=FunctionOnContactZero
1465                else:
1466                    coeff_name='y_contact_reduced'
1467                    FS=ReducedFunctionOnContactZero
1468                coeff=mypde.getCoefficient(coeff_name)
1469                self.assertEqual((coeff.getShape(),coeff.getFunctionSpace(), mypde.getNumEquations()),((self.N,),FS(self.domain),self.N))
1470                mypde.resetRightHandSideCoefficients()
1471                self.assertTrue(mypde.getCoefficient(coeff_name).isEmpty(),"%s is not empty after reset of right hand side coefficients"%coeff_name)
1472    
1473      def test_resetCoefficient_HomogeneousConstraint(self):      def test_resetCoefficient_HomogeneousConstraint(self):
1474          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1475          x=self.domain.getX()          x=self.domain.getX()
# Line 1378  class Test_LinearPDE_noLumping(Test_line Line 1477  class Test_LinearPDE_noLumping(Test_line
1477          u1=mypde.getSolution()          u1=mypde.getSolution()
1478          mypde.setValue(Y=2.)          mypde.setValue(Y=2.)
1479          u2=mypde.getSolution()          u2=mypde.getSolution()
1480          self.failUnless(self.check(u2,2*u1),'solution is wrong.')          self.assertTrue(self.check(u2,2*u1),'solution is wrong.')
1481    
1482      def test_resetCoefficient_InHomogeneousConstraint(self):      def test_resetCoefficient_InHomogeneousConstraint(self):
1483          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1484          mypde.setSymmetryOn()          mypde.setSymmetryOn()
1485            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1486          x=self.domain.getX()          x=self.domain.getX()
1487          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]))
1488          u1=mypde.getSolution(verbose=self.VERBOSE)          u1=mypde.getSolution()
1489          mypde.setValue(Y=2.,D=2)          mypde.setValue(Y=2.,D=2)
1490          u2=mypde.getSolution(verbose=self.VERBOSE)          u2=mypde.getSolution()
1491          self.failUnless(self.check(u2,u1),'first solution is wrong.')          self.assertTrue(self.check(u2,u1),'first solution is wrong.')
1492          u2=mypde.getSolution(verbose=self.VERBOSE)          u2=mypde.getSolution()
1493          self.failUnless(self.check(u2,u1),'first solution is wrong.')          self.assertTrue(self.check(u2,u1),'first solution is wrong.')
1494          mypde.setValue(r=2,Y=4.)          mypde.setValue(r=2,Y=4.)
1495          u2=mypde.getSolution(verbose=self.VERBOSE)          u2=mypde.getSolution()
1496          self.failUnless(self.check(u2,2*u1),'second solution is wrong.')          self.assertTrue(self.check(u2,2*u1),'second solution is wrong.')
1497    
1498        def test_Status(self):
1499            DIM=self.domain.getDim()
1500            x=self.domain.getX()
1501            mypde=LinearPDE(self.domain,debug=self.DEBUG)
1502            mypde.getSolverOptions().setSymmetryOn()
1503            mypde.getSolverOptions().setTolerance(self.RES_TOL)
1504            mypde.setValue(A=kronecker(self.domain), q=whereZero(x[0])+whereZero(x[0]-1.), Y=2.)
1505            x1=self.domain.getX()
1506            u1_ref=x1[0]*(1.-x1[0])
1507            u1=mypde.getSolution()
1508            error1=Lsup(u1-u1_ref)/Lsup(u1_ref)
1509            self.assertTrue(mypde.getDomainStatus() == mypde.getSystemStatus(), "status of first pde does not match domain status.")
1510            try:
1511                self.domain.setX(x*5)
1512            except:
1513                # setX not supported
1514                return
1515            self.assertTrue(mypde.getDomainStatus() != mypde.getSystemStatus(), "status of first pde matches updated domain status.")
1516            x2=self.domain.getX()
1517            u2_ref=x2[0]*(5.-x2[0])
1518            u2=mypde.getSolution()
1519            error2=Lsup(u2-u2_ref)/Lsup(u2_ref)
1520            self.assertTrue(error2 <= max(error1,self.RES_TOL)*10., "solution of second PDE wrong.")
1521            self.assertTrue(mypde.getDomainStatus() == mypde.getSystemStatus(), "status of second pde does not match domain status.")
1522    
1523      def test_symmetryCheckTrue_System(self):      def test_symmetryCheckTrue_System(self):
1524          d=self.domain.getDim()          d=self.domain.getDim()
# Line 1404  class Test_LinearPDE_noLumping(Test_line Line 1529  class Test_LinearPDE_noLumping(Test_line
1529          D=3*numpy.ones((self.N,self.N))          D=3*numpy.ones((self.N,self.N))
1530          d=4*numpy.ones((self.N,self.N))          d=4*numpy.ones((self.N,self.N))
1531          d_contact=5*numpy.ones((self.N,self.N))          d_contact=5*numpy.ones((self.N,self.N))
1532          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}
1533          self.failUnless(mypde.checkSymmetry(verbose=False),"symmetry detected")          if self.domain.supportsContactElements():
1534                    pars["d_contact"]=d_contact
1535                    pars["d_contact_reduced"]=-d_contact
1536            mypde.setValue(**pars)
1537            self.assertTrue(mypde.checkSymmetry(verbose=False),"symmetry detected")
1538    
1539      def test_symmetryCheckFalse_A_System(self):      def test_symmetryCheckFalse_A_System(self):
1540          d=self.domain.getDim()          d=self.domain.getDim()
# Line 1413  class Test_LinearPDE_noLumping(Test_line Line 1542  class Test_LinearPDE_noLumping(Test_line
1542          A=numpy.ones((self.N,d,self.N,d))          A=numpy.ones((self.N,d,self.N,d))
1543          A[1,1,1,0]=0.          A[1,1,1,0]=0.
1544          mypde.setValue(A=A)          mypde.setValue(A=A)
1545          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1546      def test_symmetryCheckFalse_BC_System(self):      def test_symmetryCheckFalse_BC_System(self):
1547          d=self.domain.getDim()          d=self.domain.getDim()
1548          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
# Line 1421  class Test_LinearPDE_noLumping(Test_line Line 1550  class Test_LinearPDE_noLumping(Test_line
1550          B=2*numpy.ones((self.N,d,self.N))          B=2*numpy.ones((self.N,d,self.N))
1551          B[0,0,1]=1.          B[0,0,1]=1.
1552          mypde.setValue(B=B,C=C)          mypde.setValue(B=B,C=C)
1553          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1554    
1555      def test_symmetryCheckFalse_D_System(self):      def test_symmetryCheckFalse_D_System(self):
1556          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1557          D=3*numpy.ones((self.N,self.N))          D=3*numpy.ones((self.N,self.N))
1558          D[0,1]=0.          D[0,1]=0.
1559          mypde.setValue(D=D)          mypde.setValue(D=D)
1560          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1561    
1562      def test_symmetryCheckFalse_d_System(self):      def test_symmetryCheckFalse_d_System(self):
1563          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1564          d=4*numpy.ones((self.N,self.N))          d=4*numpy.ones((self.N,self.N))
1565          d[0,1]=0.          d[0,1]=0.
1566          mypde.setValue(d=d)          mypde.setValue(d=d)
1567          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1568    
1569      def test_symmetryCheckFalse_d_contact_System(self):      def test_symmetryCheckFalse_d_contact_System(self):
1570          mypde=LinearPDE(self.domain,debug=self.DEBUG)          if self.domain.supportsContactElements():
1571          d_contact=5*numpy.ones((self.N,self.N))              mypde=LinearPDE(self.domain,debug=self.DEBUG)
1572          d_contact[0,1]=0.              d_contact=5*numpy.ones((self.N,self.N))
1573          mypde.setValue(d_contact=d_contact)              d_contact[0,1]=0.
1574          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")              mypde.setValue(d_contact=d_contact)
1575                self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1576    
1577      def test_symmetryCheckFalse_A_reduced_System(self):      def test_symmetryCheckFalse_A_reduced_System(self):
1578          d=self.domain.getDim()          d=self.domain.getDim()
# Line 1450  class Test_LinearPDE_noLumping(Test_line Line 1580  class Test_LinearPDE_noLumping(Test_line
1580          A=numpy.ones((self.N,d,self.N,d))          A=numpy.ones((self.N,d,self.N,d))
1581          A[1,1,1,0]=0.          A[1,1,1,0]=0.
1582          mypde.setValue(A_reduced=A)          mypde.setValue(A_reduced=A)
1583          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1584      def test_symmetryCheckFalse_BC_reduced_System(self):      def test_symmetryCheckFalse_BC_reduced_System(self):
1585          d=self.domain.getDim()          d=self.domain.getDim()
1586          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
# Line 1458  class Test_LinearPDE_noLumping(Test_line Line 1588  class Test_LinearPDE_noLumping(Test_line
1588          B=2*numpy.ones((self.N,d,self.N))          B=2*numpy.ones((self.N,d,self.N))
1589          B[0,0,1]=1.          B[0,0,1]=1.
1590          mypde.setValue(B_reduced=B,C_reduced=C)          mypde.setValue(B_reduced=B,C_reduced=C)
1591          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1592    
1593      def test_symmetryCheckFalse_D_reduced_System(self):      def test_symmetryCheckFalse_D_reduced_System(self):
1594          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1595          D=3*numpy.ones((self.N,self.N))          D=3*numpy.ones((self.N,self.N))
1596          D[0,1]=0.          D[0,1]=0.
1597          mypde.setValue(D_reduced=D)          mypde.setValue(D_reduced=D)
1598          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1599    
1600      def test_symmetryCheckFalse_d_reduced_System(self):      def test_symmetryCheckFalse_d_reduced_System(self):
1601          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1602          d=4*numpy.ones((self.N,self.N))          d=4*numpy.ones((self.N,self.N))
1603          d[0,1]=0.          d[0,1]=0.
1604          mypde.setValue(d_reduced=d)          mypde.setValue(d_reduced=d)
1605          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1606    
1607      def test_symmetryCheckFalse_d_contact_reduced_System(self):      def test_symmetryCheckFalse_d_contact_reduced_System(self):
1608          mypde=LinearPDE(self.domain,debug=self.DEBUG)          if self.domain.supportsContactElements():
1609          d_contact=5*numpy.ones((self.N,self.N))              mypde=LinearPDE(self.domain,debug=self.DEBUG)
1610          d_contact[0,1]=0.              d_contact=5*numpy.ones((self.N,self.N))
1611          mypde.setValue(d_contact_reduced=d_contact)              d_contact[0,1]=0.
1612          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")              mypde.setValue(d_contact_reduced=d_contact)
1613                self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1614    
1615      def test_symmetryCheckTrue_Scalar(self):      def test_symmetryCheckTrue_Scalar(self):
1616          d=self.domain.getDim()          d=self.domain.getDim()
# Line 1490  class Test_LinearPDE_noLumping(Test_line Line 1621  class Test_LinearPDE_noLumping(Test_line
1621          D=3          D=3
1622          d=4          d=4
1623          d_contact=5          d_contact=5
1624          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}
1625          self.failUnless(mypde.checkSymmetry(verbose=False),"symmetry detected")          if self.domain.supportsContactElements():
1626                    pars["d_contact"]=d_contact
1627                    pars["d_contact_reduced"]=-d_contact
1628            mypde.setValue(**pars)
1629            self.assertTrue(mypde.checkSymmetry(verbose=False),"symmetry detected")
1630    
1631      def test_symmetryCheckFalse_A_Scalar(self):      def test_symmetryCheckFalse_A_Scalar(self):
1632          d=self.domain.getDim()          d=self.domain.getDim()
# Line 1499  class Test_LinearPDE_noLumping(Test_line Line 1634  class Test_LinearPDE_noLumping(Test_line
1634          A=numpy.ones((d,d))          A=numpy.ones((d,d))
1635          A[1,0]=0.          A[1,0]=0.
1636          mypde.setValue(A=A)          mypde.setValue(A=A)
1637          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1638      def test_symmetryCheckFalse_BC_Scalar(self):      def test_symmetryCheckFalse_BC_Scalar(self):
1639          d=self.domain.getDim()          d=self.domain.getDim()
1640          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
# Line 1507  class Test_LinearPDE_noLumping(Test_line Line 1642  class Test_LinearPDE_noLumping(Test_line
1642          B=2*numpy.ones((d,))          B=2*numpy.ones((d,))
1643          B[0]=1.          B[0]=1.
1644          mypde.setValue(B=B,C=C)          mypde.setValue(B=B,C=C)
1645          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1646      def test_symmetryCheckFalse_A_reduced_Scalar(self):      def test_symmetryCheckFalse_A_reduced_Scalar(self):
1647          d=self.domain.getDim()          d=self.domain.getDim()
1648          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1649          A=numpy.ones((d,d))          A=numpy.ones((d,d))
1650          A[1,0]=0.          A[1,0]=0.
1651          mypde.setValue(A_reduced=A)          mypde.setValue(A_reduced=A)
1652          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1653      def test_symmetryCheckFalse_BC_reduced_Scalar(self):      def test_symmetryCheckFalse_BC_reduced_Scalar(self):
1654          d=self.domain.getDim()          d=self.domain.getDim()
1655          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
# Line 1522  class Test_LinearPDE_noLumping(Test_line Line 1657  class Test_LinearPDE_noLumping(Test_line
1657          B=2*numpy.ones((d,))          B=2*numpy.ones((d,))
1658          B[0]=1.          B[0]=1.
1659          mypde.setValue(B_reduced=B,C_reduced=C)          mypde.setValue(B_reduced=B,C_reduced=C)
1660          self.failUnless(not mypde.checkSymmetry(verbose=False),"symmetry detected")          self.assertTrue(not mypde.checkSymmetry(verbose=False),"symmetry detected")
1661      #      #
1662      #   solver checks (single PDE)      #   solver checks (single PDE)
1663      #      #
1664      def test_symmetryOnIterative(self):      def test_symmetryOnIterative(self):
1665          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1666          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1667          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1668          self.failUnless(self.check(u,1.),'solution is wrong.')          u=mypde.getSolution()
1669            self.assertTrue(self.check(u,1.),'solution is wrong.')
1670      def test_symmetryOnDirect(self):      def test_symmetryOnDirect(self):
1671          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1672          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1673          mypde.setSolverMethod(mypde.DIRECT)          mypde.getSolverOptions().setSolverMethod(SolverOptions.DIRECT)
1674          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1675          self.failUnless(self.check(u,1.),'solution is wrong.')          u=mypde.getSolution()
1676            self.assertTrue(self.check(u,1.),'solution is wrong.')
1677      def test_PCG_JACOBI(self):      def test_PCG_JACOBI(self):
1678          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1679          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1680          mypde.setSolverMethod(mypde.PCG,mypde.JACOBI)          mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
1681          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
1682          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1683            u=mypde.getSolution()
1684            self.assertTrue(self.check(u,1.),'solution is wrong.')
1685        def test_PCG_GAUSS_SEIDEL(self):
1686            mypde=LinearPDE(self.domain,debug=self.DEBUG)
1687            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1688            mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
1689            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
1690            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1691            u=mypde.getSolution()
1692            self.assertTrue(self.check(u,1.),'solution is wrong.')
1693        def test_PCG_AMG(self):
1694            if self.order!=2:
1695                if getEscriptParamInt('DISABLE_AMG', 0):
1696                    print("AMG test disabled on MPI build")
1697                    return
1698                mypde=LinearPDE(self.domain,debug=self.DEBUG)
1699                mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1700                mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
1701                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
1702                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1703                u=mypde.getSolution()
1704                self.assertTrue(self.check(u,1.),'solution is wrong.')
1705      def test_PCG_ILU0(self):      def test_PCG_ILU0(self):
1706          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1707          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1708          mypde.setSolverMethod(mypde.PCG,mypde.ILU0)          mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
1709          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
1710          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1711            u=mypde.getSolution()
1712            self.assertTrue(self.check(u,1.),'solution is wrong.')
1713      def test_PCG_RILU(self):      def test_PCG_RILU(self):
1714          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1715          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1716          mypde.setSolverMethod(mypde.PCG,mypde.RILU)          mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
1717          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)
1718          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1719            u=mypde.getSolution()
1720            self.assertTrue(self.check(u,1.),'solution is wrong.')
1721        def test_PCG_REC_ILU(self):
1722            mypde=LinearPDE(self.domain,debug=self.DEBUG)
1723            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1724            mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
1725            mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)
1726            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1727            u=mypde.getSolution()
1728            self.assertTrue(self.check(u,1.),'solution is wrong.')
1729      def test_DIRECT(self):      def test_DIRECT(self):
1730          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1731          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1732          mypde.setSolverMethod(mypde.DIRECT)          mypde.getSolverOptions().setSolverMethod(SolverOptions.DIRECT)
1733          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1734          self.failUnless(self.check(u,1.),'solution is wrong.')          u=mypde.getSolution()
1735            self.assertTrue(self.check(u,1.),'solution is wrong.')
1736      def test_BICGSTAB_JACOBI(self):      def test_BICGSTAB_JACOBI(self):
1737          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1738      mypde.setSolverMethod(mypde.BICGSTAB,mypde.JACOBI)          mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)
1739            mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
1740            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1741            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1742            u=mypde.getSolution()
1743            self.assertTrue(self.check(u,1.),'solution is wrong.')
1744        def test_BICGSTAB_GAUSS_SEIDEL(self):
1745            mypde=LinearPDE(self.domain,debug=self.DEBUG)
1746            mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)
1747            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
1748          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1749          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1750          self.failUnless(self.check(u,1.),'solution is wrong.')          u=mypde.getSolution()
1751            self.assertTrue(self.check(u,1.),'solution is wrong.')
1752        def test_BICGSTAB_AMG(self):
1753            if self.order!=2:
1754                if getEscriptParamInt('DISABLE_AMG', 0):
1755                    print("AMG test disabled on MPI build")
1756                    return    
1757                mypde=LinearPDE(self.domain,debug=self.DEBUG)
1758                mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)
1759                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
1760                mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1761                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1762                u=mypde.getSolution()
1763                self.assertTrue(self.check(u,1.),'solution is wrong.')
1764      def test_BICGSTAB_ILU0(self):      def test_BICGSTAB_ILU0(self):
1765          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1766          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1767      mypde.setSolverMethod(mypde.BICGSTAB,mypde.ILU0)          mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)
1768          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
1769          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1770            u=mypde.getSolution()
1771            self.assertTrue(self.check(u,1.),'solution is wrong.')
1772      def test_BICGSTAB_RILU(self):      def test_BICGSTAB_RILU(self):
1773          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1774          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1775      mypde.setSolverMethod(mypde.BICGSTAB,mypde.RILU)          mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)
1776          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)
1777          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1778            u=mypde.getSolution()
1779            self.assertTrue(self.check(u,1.),'solution is wrong.')
1780        def test_BICGSTAB_REC_ILU(self):
1781            mypde=LinearPDE(self.domain,debug=self.DEBUG)
1782            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1783            mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)
1784            mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)
1785            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1786            u=mypde.getSolution()
1787            self.assertTrue(self.check(u,1.),'solution is wrong.')
1788      def test_MINRES_JACOBI(self):      def test_MINRES_JACOBI(self):
1789          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1790      mypde.setSolverMethod(mypde.MINRES,mypde.JACOBI)          mypde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
1791            mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
1792          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1793          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1794          self.failUnless(self.check(u,1.),'solution is wrong.')          u=mypde.getSolution()
1795            self.assertTrue(self.check(u,1.),'solution is wrong.')
1796        def test_MINRES_GAUSS_SEIDEL(self):
1797            mypde=LinearPDE(self.domain,debug=self.DEBUG)
1798            mypde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
1799            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
1800            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1801            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1802            u=mypde.getSolution()
1803            self.assertTrue(self.check(u,1.),'solution is wrong.')
1804        def test_MINRES_AMG(self):
1805            if self.order!=2:
1806                if getEscriptParamInt('DISABLE_AMG',0):
1807                    print("AMG test disabled on MPI build")
1808                    return                
1809                mypde=LinearPDE(self.domain,debug=self.DEBUG)
1810                mypde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
1811                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
1812                mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1813                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1814                u=mypde.getSolution()
1815                self.assertTrue(self.check(u,1.),'solution is wrong.')
1816      def test_MINRES_ILU0(self):      def test_MINRES_ILU0(self):
1817          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1818          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1819      mypde.setSolverMethod(mypde.MINRES,mypde.ILU0)          mypde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
1820          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
1821          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1822            u=mypde.getSolution()
1823            self.assertTrue(self.check(u,1.),'solution is wrong.')
1824      def test_MINRES_RILU(self):      def test_MINRES_RILU(self):
1825          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1826          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1827      mypde.setSolverMethod(mypde.MINRES,mypde.RILU)          mypde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
1828          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)
1829          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1830            u=mypde.getSolution()
1831            self.assertTrue(self.check(u,1.),'solution is wrong.')
1832        def test_MINRES_REC_ILU(self):
1833            mypde=LinearPDE(self.domain,debug=self.DEBUG)
1834            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1835            mypde.getSolverOptions().setSolverMethod(SolverOptions.MINRES)
1836            mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)
1837            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1838            u=mypde.getSolution()
1839            self.assertTrue(self.check(u,1.),'solution is wrong.')
1840      def test_TFQMR_JACOBI(self):      def test_TFQMR_JACOBI(self):
1841          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1842      mypde.setSolverMethod(mypde.TFQMR,mypde.JACOBI)          mypde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
1843            mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
1844            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1845            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1846            u=mypde.getSolution()
1847            self.assertTrue(self.check(u,1.),'solution is wrong.')
1848        def test_TFQMR_GAUSS_SEIDEL(self):
1849            mypde=LinearPDE(self.domain,debug=self.DEBUG)
1850            mypde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
1851            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
1852          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1853          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1854          self.failUnless(self.check(u,1.),'solution is wrong.')          u=mypde.getSolution()
1855            self.assertTrue(self.check(u,1.),'solution is wrong.')
1856        def test_TFQMR_AMG(self):
1857            if self.order!=2:
1858                if getEscriptParamInt('DISABLE_AMG', 0):
1859                    print("AMG test disabled on MPI build")
1860                    return
1861                mypde=LinearPDE(self.domain,debug=self.DEBUG)
1862                mypde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
1863                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
1864                mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1865                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1866                u=mypde.getSolution()
1867                self.assertTrue(self.check(u,1.),'solution is wrong.')
1868      def test_TFQMR_ILU0(self):      def test_TFQMR_ILU0(self):
1869          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1870          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1871      mypde.setSolverMethod(mypde.TFQMR,mypde.ILU0)          mypde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
1872          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
1873          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1874            u=mypde.getSolution()
1875            self.assertTrue(self.check(u,1.),'solution is wrong.')
1876      def test_TFQMR_RILU(self):      def test_TFQMR_RILU(self):
1877          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1878          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1879      mypde.setSolverMethod(mypde.TFQMR,mypde.RILU)          mypde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
1880          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)
1881          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1882            u=mypde.getSolution()
1883            self.assertTrue(self.check(u,1.),'solution is wrong.')
1884        def test_TFQMR_REC_ILU(self):
1885            mypde=LinearPDE(self.domain,debug=self.DEBUG)
1886            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1887            mypde.getSolverOptions().setSolverMethod(SolverOptions.TFQMR)
1888            mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)
1889            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1890            u=mypde.getSolution()
1891            self.assertTrue(self.check(u,1.),'solution is wrong.')
1892      def test_PRES20_JACOBI(self):      def test_PRES20_JACOBI(self):
1893          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1894          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1895      mypde.setSolverMethod(mypde.PRES20,mypde.JACOBI)          mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)
1896          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
1897          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1898            u=mypde.getSolution()
1899            self.assertTrue(self.check(u,1.),'solution is wrong.')
1900        def test_PRES20_GAUSS_SEIDEL(self):
1901            mypde=LinearPDE(self.domain,debug=self.DEBUG)
1902            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1903            mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)
1904            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
1905            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1906            u=mypde.getSolution()
1907            self.assertTrue(self.check(u,1.),'solution is wrong.')
1908        def test_PRES20_AMG(self):
1909            if self.order!=2:
1910                if getEscriptParamInt('DISABLE_AMG', 0):
1911                    print("AMG test disabled on MPI build")
1912                    return
1913                mypde=LinearPDE(self.domain,debug=self.DEBUG)
1914                mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1915                mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)
1916                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
1917                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1918                u=mypde.getSolution()
1919                self.assertTrue(self.check(u,1.),'solution is wrong.')
1920      def test_PRES20_ILU0(self):      def test_PRES20_ILU0(self):
1921          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1922          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1923      mypde.setSolverMethod(mypde.PRES20,mypde.ILU0)          mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)
1924          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
1925          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1926            u=mypde.getSolution()
1927            self.assertTrue(self.check(u,1.),'solution is wrong.')
1928      def test_PRES20_RILU(self):      def test_PRES20_RILU(self):
1929          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1930          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1931      mypde.setSolverMethod(mypde.PRES20,mypde.RILU)          mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)
1932          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)
1933          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1934            u=mypde.getSolution()
1935            self.assertTrue(self.check(u,1.),'solution is wrong.')
1936        def test_PRES20_REC_ILU(self):
1937            mypde=LinearPDE(self.domain,debug=self.DEBUG)
1938            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1939            mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)
1940            mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)
1941            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1942            u=mypde.getSolution()
1943            self.assertTrue(self.check(u,1.),'solution is wrong.')
1944      def test_GMRESnoRestart_JACOBI(self):      def test_GMRESnoRestart_JACOBI(self):
1945          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1946          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1947      mypde.setSolverMethod(mypde.GMRES,mypde.JACOBI)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
1948          u=mypde.getSolution(verbose=self.VERBOSE,truncation=50)          mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
1949          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1950            mypde.getSolverOptions().setTruncation(50)
1951            u=mypde.getSolution()
1952            self.assertTrue(self.check(u,1.),'solution is wrong.')
1953        def test_GMRESnoRestart_GAUSS_SEIDEL(self):
1954            mypde=LinearPDE(self.domain,debug=self.DEBUG)
1955            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1956            mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
1957            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
1958            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1959            mypde.getSolverOptions().setTruncation(50)
1960            u=mypde.getSolution()
1961            self.assertTrue(self.check(u,1.),'solution is wrong.')
1962        def test_GMRESnoRestart_AMG(self):
1963            if self.order!=2:
1964                if getEscriptParamInt('DISABLE_AMG', 0):
1965                    print("AMG test disabled on MPI build")
1966                    return    
1967                mypde=LinearPDE(self.domain,debug=self.DEBUG)
1968                mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1969                mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
1970                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
1971                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1972                mypde.getSolverOptions().setTruncation(50)
1973                u=mypde.getSolution()
1974                self.assertTrue(self.check(u,1.),'solution is wrong.')
1975      def test_GMRESnoRestart_ILU0(self):      def test_GMRESnoRestart_ILU0(self):
1976          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1977          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1978      mypde.setSolverMethod(mypde.GMRES,mypde.ILU0)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
1979          u=mypde.getSolution(verbose=self.VERBOSE,truncation=50)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
1980          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1981          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setTruncation(50)                        
1982            u=mypde.getSolution()
1983            self.assertTrue(self.check(u,1.),'solution is wrong.')
1984      def test_GMRESnoRestart_RILU(self):      def test_GMRESnoRestart_RILU(self):
1985          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
1986          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1987      mypde.setSolverMethod(mypde.GMRES,mypde.RILU)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
1988          u=mypde.getSolution(verbose=self.VERBOSE,truncation=50)          mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)
1989          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1990          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setTruncation(50)
1991            u=mypde.getSolution()
1992            self.assertTrue(self.check(u,1.),'solution is wrong.')
1993        def test_GMRESnoRestart_REC_ILU(self):
1994            mypde=LinearPDE(self.domain,debug=self.DEBUG)
1995            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
1996            mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
1997            mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)
1998            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
1999            mypde.getSolverOptions().setTruncation(50)
2000            u=mypde.getSolution()
2001            self.assertTrue(self.check(u,1.),'solution is wrong.')
2002      def test_GMRES_JACOBI(self):      def test_GMRES_JACOBI(self):
2003          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2004          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2005      mypde.setSolverMethod(mypde.GMRES,mypde.JACOBI)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2006          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
2007          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2008            u=mypde.getSolution()
2009            self.assertTrue(self.check(u,1.),'solution is wrong.')
2010        def test_GMRES_GAUSS_SEIDEL(self):
2011            mypde=LinearPDE(self.domain,debug=self.DEBUG)
2012            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2013            mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2014            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
2015            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2016            u=mypde.getSolution()
2017            self.assertTrue(self.check(u,1.),'solution is wrong.')
2018        def test_GMRES_AMG(self):
2019            if self.order!=2:
2020                if getEscriptParamInt('DISABLE_AMG', 0):
2021                    print("AMG test disabled on MPI build")
2022                    return    
2023                mypde=LinearPDE(self.domain,debug=self.DEBUG)
2024                mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2025                mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2026                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
2027                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2028                u=mypde.getSolution()
2029                self.assertTrue(self.check(u,1.),'solution is wrong.')        
2030      def test_GMRES_ILU0(self):      def test_GMRES_ILU0(self):
2031          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2032          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2033      mypde.setSolverMethod(mypde.GMRES,mypde.ILU0)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2034          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
2035          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2036            u=mypde.getSolution()
2037            self.assertTrue(self.check(u,1.),'solution is wrong.')
2038      def test_GMRES_RILU(self):      def test_GMRES_RILU(self):
2039          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2040          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2041      mypde.setSolverMethod(mypde.GMRES,mypde.RILU)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2042          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)
2043          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2044            u=mypde.getSolution()
2045            self.assertTrue(self.check(u,1.),'solution is wrong.')
2046        def test_GMRES_REC_ILU(self):
2047            mypde=LinearPDE(self.domain,debug=self.DEBUG)
2048            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2049            mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2050            mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)
2051            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2052            u=mypde.getSolution()
2053            self.assertTrue(self.check(u,1.),'solution is wrong.')
2054      def test_GMRES_truncation_restart_JACOBI(self):      def test_GMRES_truncation_restart_JACOBI(self):
2055          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2056          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2057      mypde.setSolverMethod(mypde.GMRES,mypde.JACOBI)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2058          u=mypde.getSolution(verbose=self.VERBOSE,truncation=10,restart=20)          mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
2059          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2060            mypde.getSolverOptions().setTruncation(10)
2061            mypde.getSolverOptions().setRestart(20)
2062            u=mypde.getSolution()
2063            self.assertTrue(self.check(u,1.),'solution is wrong.')
2064        def test_GMRES_truncation_restart_GAUSS_SEIDEL(self):
2065            mypde=LinearPDE(self.domain,debug=self.DEBUG)
2066            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2067            mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2068            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
2069            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2070            mypde.getSolverOptions().setTruncation(10)
2071            mypde.getSolverOptions().setRestart(20)
2072            u=mypde.getSolution()
2073            self.assertTrue(self.check(u,1.),'solution is wrong.')
2074        def test_GMRES_truncation_restart_AMG(self):
2075            if self.order!=2:
2076                if getEscriptParamInt('DISABLE_AMG', 0):
2077                    print("AMG test disabled on MPI build")
2078                    return    
2079                mypde=LinearPDE(self.domain,debug=self.DEBUG)
2080                mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2081                mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2082                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
2083                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2084                mypde.getSolverOptions().setTruncation(10)
2085                mypde.getSolverOptions().setRestart(20)
2086                u=mypde.getSolution()
2087                self.assertTrue(self.check(u,1.),'solution is wrong.')
2088      def test_GMRES_truncation_restart_ILU0(self):      def test_GMRES_truncation_restart_ILU0(self):
2089          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2090          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2091      mypde.setSolverMethod(mypde.GMRES,mypde.ILU0)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2092          u=mypde.getSolution(verbose=self.VERBOSE,truncation=10,restart=20)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
2093          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2094            mypde.getSolverOptions().setTruncation(10)
2095            mypde.getSolverOptions().setRestart(20)
2096            u=mypde.getSolution()
2097            self.assertTrue(self.check(u,1.),'solution is wrong.')
2098      def test_GMRES_truncation_restart_RILU(self):      def test_GMRES_truncation_restart_RILU(self):
2099          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2100          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)          mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2101      mypde.setSolverMethod(mypde.GMRES,mypde.RILU)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2102          u=mypde.getSolution(verbose=self.VERBOSE,truncation=10,restart=20)          mypde.getSolverOptions().setPreconditioner(SolverOptions.RILU)
2103          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2104            mypde.getSolverOptions().setTruncation(10)
2105            mypde.getSolverOptions().setRestart(20)
2106            u=mypde.getSolution()
2107            self.assertTrue(self.check(u,1.),'solution is wrong.')
2108        def test_GMRES_truncation_restart_REC_ILU(self):
2109            mypde=LinearPDE(self.domain,debug=self.DEBUG)
2110            mypde.setValue(A=kronecker(self.domain),D=1.,Y=1.)
2111            mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2112            mypde.getSolverOptions().setPreconditioner(SolverOptions.REC_ILU)
2113            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2114            mypde.getSolverOptions().setTruncation(10)
2115            mypde.getSolverOptions().setRestart(20)
2116            u=mypde.getSolution()
2117            self.assertTrue(self.check(u,1.),'solution is wrong.')
2118      #      #
2119      #   solver checks (PDE system)      #   solver checks (PDE system)
2120      #      #
# Line 1702  class Test_LinearPDE_noLumping(Test_line Line 2128  class Test_LinearPDE_noLumping(Test_line
2128              Y[i]+=i              Y[i]+=i
2129          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2130          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2131          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2132          self.failUnless(self.check(u,1.),'solution is wrong.')          u=mypde.getSolution()
2133            self.assertTrue(self.check(u,1.),'solution is wrong.')
2134      def test_symmetryOnDirect_System(self):      def test_symmetryOnDirect_System(self):
2135          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2136          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 1714  class Test_LinearPDE_noLumping(Test_line Line 2141  class Test_LinearPDE_noLumping(Test_line
2141              Y[i]+=i              Y[i]+=i
2142          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2143          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2144          mypde.setSolverMethod(mypde.DIRECT)          mypde.getSolverOptions().setSolverMethod(SolverOptions.DIRECT)
2145          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2146          self.failUnless(self.check(u,1.),'solution is wrong.')          u=mypde.getSolution()
2147            self.assertTrue(self.check(u,1.),'solution is wrong.')
2148      def test_PCG_JACOBI_System(self):      def test_PCG_JACOBI_System(self):
2149          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2150          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 1727  class Test_LinearPDE_noLumping(Test_line Line 2155  class Test_LinearPDE_noLumping(Test_line
2155              Y[i]+=i              Y[i]+=i
2156          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2157          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2158          mypde.setSolverMethod(mypde.PCG,mypde.JACOBI)          mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
2159          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
2160          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2161            u=mypde.getSolution()
2162            self.assertTrue(self.check(u,1.),'solution is wrong.')
2163        def test_PCG_GAUSS_SEIDEL_System(self):
2164            A=Tensor4(0.,Function(self.domain))
2165            D=Tensor(1.,Function(self.domain))
2166            Y=Vector(self.domain.getDim(),Function(self.domain))
2167            for i in range(self.domain.getDim()):
2168                A[i,:,i,:]=kronecker(self.domain)
2169                D[i,i]+=i
2170                Y[i]+=i
2171            mypde=LinearPDE(self.domain,debug=self.DEBUG)
2172            mypde.setValue(A=A,D=D,Y=Y)
2173            mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
2174            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
2175            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2176            u=mypde.getSolution()
2177            self.assertTrue(self.check(u,1.),'solution is wrong.')
2178        def test_PCG_AMG_System(self):
2179            if self.order!=2:
2180                if getEscriptParamInt('DISABLE_AMG', 0):
2181                    print("AMG test disabled on MPI build")
2182                    return
2183                A=Tensor4(0.,Function(self.domain))
2184                D=Tensor(1.,Function(self.domain))
2185                Y=Vector(self.domain.getDim(),Function(self.domain))
2186                for i in range(self.domain.getDim()):
2187                    A[i,:,i,:]=kronecker(self.domain)
2188                    D[i,i]+=i
2189                    Y[i]+=i
2190                mypde=LinearPDE(self.domain,debug=self.DEBUG)
2191                mypde.setValue(A=A,D=D,Y=Y)
2192                mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
2193                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
2194                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2195                u=mypde.getSolution()
2196                self.assertTrue(self.check(u,1.),'solution is wrong.')
2197      def test_PCG_ILU0_System(self):      def test_PCG_ILU0_System(self):
2198          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2199          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 1740  class Test_LinearPDE_noLumping(Test_line Line 2204  class Test_LinearPDE_noLumping(Test_line
2204              Y[i]+=i              Y[i]+=i
2205          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2206          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2207          mypde.setSolverMethod(mypde.PCG,mypde.ILU0)          mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG)
2208          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
2209          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2210            u=mypde.getSolution()
2211            self.assertTrue(self.check(u,1.),'solution is wrong.')
2212      def test_DIRECT_System(self):      def test_DIRECT_System(self):
2213          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2214          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 1753  class Test_LinearPDE_noLumping(Test_line Line 2219  class Test_LinearPDE_noLumping(Test_line
2219              Y[i]+=i              Y[i]+=i
2220          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2221          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2222          mypde.setSolverMethod(mypde.DIRECT)          mypde.getSolverOptions().setSolverMethod(SolverOptions.DIRECT)
2223          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2224          self.failUnless(self.check(u,1.),'solution is wrong.')          u=mypde.getSolution()
2225            self.assertTrue(self.check(u,1.),'solution is wrong.')
2226      def test_BICGSTAB_JACOBI_System(self):      def test_BICGSTAB_JACOBI_System(self):
2227          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2228          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 1766  class Test_LinearPDE_noLumping(Test_line Line 2233  class Test_LinearPDE_noLumping(Test_line
2233              Y[i]+=i              Y[i]+=i
2234          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2235          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2236      mypde.setSolverMethod(mypde.BICGSTAB,mypde.JACOBI)          mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)
2237          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
2238          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2239            u=mypde.getSolution()
2240            self.assertTrue(self.check(u,1.),'solution is wrong.')
2241        def test_BICGSTAB_GAUSS_SEIDEL_System(self):
2242            A=Tensor4(0.,Function(self.domain))
2243            D=Tensor(1.,Function(self.domain))
2244            Y=Vector(self.domain.getDim(),Function(self.domain))
2245            for i in range(self.domain.getDim()):
2246                A[i,:,i,:]=kronecker(self.domain)
2247                D[i,i]+=i
2248                Y[i]+=i
2249            mypde=LinearPDE(self.domain,debug=self.DEBUG)
2250            mypde.setValue(A=A,D=D,Y=Y)
2251            mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)
2252            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
2253            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2254            u=mypde.getSolution()
2255            self.assertTrue(self.check(u,1.),'solution is wrong.')
2256        def test_BICGSTAB_AMG_System(self):
2257            if self.order!=2:
2258                if getEscriptParamInt('DISABLE_AMG', 0):
2259                    print("AMG test disabled on MPI build")
2260                    return    
2261                A=Tensor4(0.,Function(self.domain))
2262                D=Tensor(1.,Function(self.domain))
2263                Y=Vector(self.domain.getDim(),Function(self.domain))
2264                for i in range(self.domain.getDim()):
2265                    A[i,:,i,:]=kronecker(self.domain)
2266                    D[i,i]+=i
2267                    Y[i]+=i
2268                mypde=LinearPDE(self.domain,debug=self.DEBUG)
2269                mypde.setValue(A=A,D=D,Y=Y)
2270                mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)
2271                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
2272                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2273                u=mypde.getSolution()
2274                self.assertTrue(self.check(u,1.),'solution is wrong.')
2275      def test_BICGSTAB_ILU0_System(self):      def test_BICGSTAB_ILU0_System(self):
2276          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2277          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 1779  class Test_LinearPDE_noLumping(Test_line Line 2282  class Test_LinearPDE_noLumping(Test_line
2282              Y[i]+=i              Y[i]+=i
2283          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2284          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2285      mypde.setSolverMethod(mypde.BICGSTAB,mypde.ILU0)          mypde.getSolverOptions().setSolverMethod(SolverOptions.BICGSTAB)
2286          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
2287          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2288            u=mypde.getSolution()
2289            self.assertTrue(self.check(u,1.),'solution is wrong.')
2290      def test_PRES20_JACOBI_System(self):      def test_PRES20_JACOBI_System(self):
2291          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2292          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 1792  class Test_LinearPDE_noLumping(Test_line Line 2297  class Test_LinearPDE_noLumping(Test_line
2297              Y[i]+=i              Y[i]+=i
2298          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2299          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2300      mypde.setSolverMethod(mypde.PRES20,mypde.JACOBI)          mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)
2301          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
2302          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2303            u=mypde.getSolution()
2304            self.assertTrue(self.check(u,1.),'solution is wrong.')
2305        def test_PRES20_GAUSS_SEIDEL_System(self):
2306            A=Tensor4(0.,Function(self.domain))
2307            D=Tensor(1.,Function(self.domain))
2308            Y=Vector(self.domain.getDim(),Function(self.domain))
2309            for i in range(self.domain.getDim()):
2310                A[i,:,i,:]=kronecker(self.domain)
2311                D[i,i]+=i
2312                Y[i]+=i
2313            mypde=LinearPDE(self.domain,debug=self.DEBUG)
2314            mypde.setValue(A=A,D=D,Y=Y)
2315            mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)
2316            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
2317            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2318            u=mypde.getSolution()
2319            self.assertTrue(self.check(u,1.),'solution is wrong.')
2320        def test_PRES20_AMG_System(self):
2321            if self.order!=2:
2322                if getEscriptParamInt('DISABLE_AMG', 0):
2323                    print("AMG test disabled on MPI build")
2324                    return
2325                A=Tensor4(0.,Function(self.domain))
2326                D=Tensor(1.,Function(self.domain))
2327                Y=Vector(self.domain.getDim(),Function(self.domain))
2328                for i in range(self.domain.getDim()):
2329                    A[i,:,i,:]=kronecker(self.domain)
2330                    D[i,i]+=i
2331                    Y[i]+=i
2332                mypde=LinearPDE(self.domain,debug=self.DEBUG)
2333                mypde.setValue(A=A,D=D,Y=Y)
2334                mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)
2335                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
2336                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2337                u=mypde.getSolution()
2338                self.assertTrue(self.check(u,1.),'solution is wrong.')
2339      def test_PRES20_ILU0_System(self):      def test_PRES20_ILU0_System(self):
2340          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2341          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 1805  class Test_LinearPDE_noLumping(Test_line Line 2346  class Test_LinearPDE_noLumping(Test_line
2346              Y[i]+=i              Y[i]+=i
2347          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2348          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2349      mypde.setSolverMethod(mypde.PRES20,mypde.ILU0)          mypde.getSolverOptions().setSolverMethod(SolverOptions.PRES20)
2350          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
2351          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2352            u=mypde.getSolution()
2353            self.assertTrue(self.check(u,1.),'solution is wrong.')
2354      def test_GMRESnoRestart_JACOBI_System(self):      def test_GMRESnoRestart_JACOBI_System(self):
2355          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2356          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 1818  class Test_LinearPDE_noLumping(Test_line Line 2361  class Test_LinearPDE_noLumping(Test_line
2361              Y[i]+=i              Y[i]+=i
2362          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2363          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2364      mypde.setSolverMethod(mypde.GMRES,mypde.JACOBI)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2365            mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
2366            # u=mypde.getSolution(verbose=self.VERBOSE,truncation=5)
2367            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2368            u=mypde.getSolution()
2369            self.assertTrue(self.check(u,1.),'solution is wrong.')
2370        def test_GMRESnoRestart_GAUSS_SEIDEL_System(self):
2371            A=Tensor4(0.,Function(self.domain))
2372            D=Tensor(1.,Function(self.domain))
2373            Y=Vector(self.domain.getDim(),Function(self.domain))
2374            for i in range(self.domain.getDim()):
2375                A[i,:,i,:]=kronecker(self.domain)
2376                D[i,i]+=i
2377                Y[i]+=i
2378            mypde=LinearPDE(self.domain,debug=self.DEBUG)
2379            mypde.setValue(A=A,D=D,Y=Y)
2380            mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2381            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
2382          # u=mypde.getSolution(verbose=self.VERBOSE,truncation=5)          # u=mypde.getSolution(verbose=self.VERBOSE,truncation=5)
2383          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2384          self.failUnless(self.check(u,1.),'solution is wrong.')          u=mypde.getSolution()
2385            self.assertTrue(self.check(u,1.),'solution is wrong.')
2386        def test_GMRESnoRestart_AMG_System(self):
2387            if self.order!=2:
2388                if getEscriptParamInt('DISABLE_AMG',0):
2389                    print("AMG test disabled on MPI build")
2390                    return        
2391                A=Tensor4(0.,Function(self.domain))
2392                D=Tensor(1.,Function(self.domain))
2393                Y=Vector(self.domain.getDim(),Function(self.domain))
2394                for i in range(self.domain.getDim()):
2395                    A[i,:,i,:]=kronecker(self.domain)
2396                    D[i,i]+=i
2397                    Y[i]+=i
2398                mypde=LinearPDE(self.domain,debug=self.DEBUG)
2399                mypde.setValue(A=A,D=D,Y=Y)
2400                mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2401                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
2402                # u=mypde.getSolution(verbose=self.VERBOSE,truncation=5)
2403                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2404                u=mypde.getSolution()
2405                self.assertTrue(self.check(u,1.),'solution is wrong.')
2406      def test_GMRESnoRestart_ILU0_System(self):      def test_GMRESnoRestart_ILU0_System(self):
2407          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2408          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 1832  class Test_LinearPDE_noLumping(Test_line Line 2413  class Test_LinearPDE_noLumping(Test_line
2413              Y[i]+=i              Y[i]+=i
2414          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2415          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2416      mypde.setSolverMethod(mypde.GMRES,mypde.ILU0)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2417            mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
2418          # u=mypde.getSolution(verbose=self.VERBOSE,truncation=5)          # u=mypde.getSolution(verbose=self.VERBOSE,truncation=5)
2419          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2420          self.failUnless(self.check(u,1.),'solution is wrong.')          u=mypde.getSolution()
2421            self.assertTrue(self.check(u,1.),'solution is wrong.')
2422      def test_GMRES_JACOBI_System(self):      def test_GMRES_JACOBI_System(self):
2423          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2424          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 1846  class Test_LinearPDE_noLumping(Test_line Line 2429  class Test_LinearPDE_noLumping(Test_line
2429              Y[i]+=i              Y[i]+=i
2430          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2431          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2432      mypde.setSolverMethod(mypde.GMRES,mypde.JACOBI)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2433          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
2434          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2435            u=mypde.getSolution()
2436            self.assertTrue(self.check(u,1.),'solution is wrong.')
2437        def test_GMRES_GAUSS_SEIDEL_System(self):
2438            A=Tensor4(0.,Function(self.domain))
2439            D=Tensor(1.,Function(self.domain))
2440            Y=Vector(self.domain.getDim(),Function(self.domain))
2441            for i in range(self.domain.getDim()):
2442                A[i,:,i,:]=kronecker(self.domain)
2443                D[i,i]+=i
2444                Y[i]+=i
2445            mypde=LinearPDE(self.domain,debug=self.DEBUG)
2446            mypde.setValue(A=A,D=D,Y=Y)
2447            mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2448            mypde.getSolverOptions().setPreconditioner(SolverOptions.GAUSS_SEIDEL)
2449            mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2450            u=mypde.getSolution()
2451            self.assertTrue(self.check(u,1.),'solution is wrong.')
2452        def test_GMRES_AMG_System(self):
2453            if self.order!=2:
2454                if getEscriptParamInt('DISABLE_AMG', 0):
2455                    print("AMG test disabled on MPI build")
2456                    return    
2457                A=Tensor4(0.,Function(self.domain))
2458                D=Tensor(1.,Function(self.domain))
2459                Y=Vector(self.domain.getDim(),Function(self.domain))
2460                for i in range(self.domain.getDim()):
2461                    A[i,:,i,:]=kronecker(self.domain)
2462                    D[i,i]+=i
2463                    Y[i]+=i
2464                mypde=LinearPDE(self.domain,debug=self.DEBUG)
2465                mypde.setValue(A=A,D=D,Y=Y)
2466                mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2467                mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG)
2468                mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2469                u=mypde.getSolution()
2470                self.assertTrue(self.check(u,1.),'solution is wrong.')
2471      def test_GMRES_ILU0_System(self):      def test_GMRES_ILU0_System(self):
2472          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2473          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 1859  class Test_LinearPDE_noLumping(Test_line Line 2478  class Test_LinearPDE_noLumping(Test_line
2478              Y[i]+=i              Y[i]+=i
2479          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2480          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2481      mypde.setSolverMethod(mypde.GMRES,mypde.ILU0)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2482          u=mypde.getSolution(verbose=self.VERBOSE)          mypde.getSolverOptions().setPreconditioner(SolverOptions.ILU0)
2483          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2484            u=mypde.getSolution()
2485            self.assertTrue(self.check(u,1.),'solution is wrong.')
2486      def test_GMRES_truncation_restart_JACOBI_System(self):      def test_GMRES_truncation_restart_JACOBI_System(self):
2487          A=Tensor4(0.,Function(self.domain))          A=Tensor4(0.,Function(self.domain))
2488          D=Tensor(1.,Function(self.domain))          D=Tensor(1.,Function(self.domain))
# Line 1872  class Test_LinearPDE_noLumping(Test_line Line 2493  class Test_LinearPDE_noLumping(Test_line
2493              Y[i]+=i              Y[i]+=i
2494          mypde=LinearPDE(self.domain,debug=self.DEBUG)          mypde=LinearPDE(self.domain,debug=self.DEBUG)
2495          mypde.setValue(A=A,D=D,Y=Y)          mypde.setValue(A=A,D=D,Y=Y)
2496      mypde.setSolverMethod(mypde.GMRES,mypde.JACOBI)          mypde.getSolverOptions().setSolverMethod(SolverOptions.GMRES)
2497          u=mypde.getSolution(verbose=self.VERBOSE,truncation=10,restart=20)          mypde.getSolverOptions().setPreconditioner(SolverOptions.JACOBI)
2498          self.failUnless(self.check(u,1.),'solution is wrong.')          mypde.getSolverOptions().setVerbosity(self.VERBOSE)
2499            mypde.getSolverOptions().setTruncation(10)
2500            mypde.getSolverOptions().setRestart(20)
2501            u=mypde.getSolution()
2502            self.assertTrue(self.check(u,1.),'solution is wrong.')
2503        def test_GMRES_truncation_restart_GAUSS_SEIDEL_System(self):
2504            A=Tensor4(0.,Function(self.domain))
2505            D=Tensor(1.,Function(self.domain))
2506            Y=Vector(self.domain.getDim(),Function(self.domain))
2507            for i in range(self.domain.getDim()):
2508