/[escript]/trunk/finley/test/python/run_models.py
ViewVC logotype

Diff of /trunk/finley/test/python/run_models.py

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

revision 1414 by gross, Thu Feb 14 10:01:43 2008 UTC revision 1562 by gross, Wed May 21 13:04:40 2008 UTC
# Line 31  except KeyError: Line 31  except KeyError:
31    
32  NE=6  NE=6
33  TOL=1.e-5  TOL=1.e-5
34  VERBOSE=False  VERBOSE=False or True
35    
36  from esys.escript import *  from esys.escript import *
37  from esys.escript.models import StokesProblemCartesian  from esys.escript.models import StokesProblemCartesian
# Line 41  class Test_Simple2DModels(unittest.TestC Line 41  class Test_Simple2DModels(unittest.TestC
41         self.domain=Rectangle(NE,NE,order=2,useFullElementOrder=True)         self.domain=Rectangle(NE,NE,order=2,useFullElementOrder=True)
42     def tearDown(self):     def tearDown(self):
43         del self.domain         del self.domain
44     def test_StokesProblemCartesian_P_0(self):     def test_StokesProblemCartesian_PCG_P_0(self):
45         ETA=1.         ETA=1.
46         P1=0.         P1=0.
47    
# Line 57  class Test_Simple2DModels(unittest.TestC Line 57  class Test_Simple2DModels(unittest.TestC
57         sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)         sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
58         u0=(1-x[0])*x[0]*[0.,1.]         u0=(1-x[0])*x[0]*[0.,1.]
59         p0=Scalar(P1,ReducedSolution(self.domain))         p0=Scalar(P1,ReducedSolution(self.domain))
60         u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100)         u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="PCG")
61                
62         error_v0=Lsup(u[0]-u0[0])         error_v0=Lsup(u[0]-u0[0])
63         error_v1=Lsup(u[1]-u0[1])/0.25         error_v1=Lsup(u[1]-u0[1])/0.25
# Line 68  class Test_Simple2DModels(unittest.TestC Line 68  class Test_Simple2DModels(unittest.TestC
68         self.failUnless(error_v0<TOL,"0-velocity error too large.")         self.failUnless(error_v0<TOL,"0-velocity error too large.")
69         self.failUnless(error_v1<TOL,"1-velocity error too large.")         self.failUnless(error_v1<TOL,"1-velocity error too large.")
70    
71     def test_StokesProblemCartesian_P_small(self):     def test_StokesProblemCartesian_PCG_P_small(self):
72         ETA=1.         ETA=1.
73         P1=1.         P1=1.
74    
# Line 84  class Test_Simple2DModels(unittest.TestC Line 84  class Test_Simple2DModels(unittest.TestC
84         sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)         sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
85         u0=(1-x[0])*x[0]*[0.,1.]         u0=(1-x[0])*x[0]*[0.,1.]
86         p0=Scalar(P1,ReducedSolution(self.domain))         p0=Scalar(P1,ReducedSolution(self.domain))
87         u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100)         u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="PCG")
88                
89         error_v0=Lsup(u[0]-u0[0])         error_v0=Lsup(u[0]-u0[0])
90         error_v1=Lsup(u[1]-u0[1])/0.25         error_v1=Lsup(u[1]-u0[1])/0.25
# Line 95  class Test_Simple2DModels(unittest.TestC Line 95  class Test_Simple2DModels(unittest.TestC
95         self.failUnless(error_v0<TOL,"0-velocity error too large.")         self.failUnless(error_v0<TOL,"0-velocity error too large.")
96         self.failUnless(error_v1<TOL,"1-velocity error too large.")         self.failUnless(error_v1<TOL,"1-velocity error too large.")
97    
98     def test_StokesProblemCartesian_P_large(self):     def test_StokesProblemCartesian_PCG_P_large(self):
99         ETA=1.         ETA=1.
100         P1=1000.         P1=1000.
101    
# Line 111  class Test_Simple2DModels(unittest.TestC Line 111  class Test_Simple2DModels(unittest.TestC
111         sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)         sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
112         u0=(1-x[0])*x[0]*[0.,1.]         u0=(1-x[0])*x[0]*[0.,1.]
113         p0=Scalar(P1,ReducedSolution(self.domain))         p0=Scalar(P1,ReducedSolution(self.domain))
114         u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100)         u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="PCG")
115                
116         error_v0=Lsup(u[0]-u0[0])         error_v0=Lsup(u[0]-u0[0])
117         error_v1=Lsup(u[1]-u0[1])/0.25         error_v1=Lsup(u[1]-u0[1])/0.25
# Line 122  class Test_Simple2DModels(unittest.TestC Line 122  class Test_Simple2DModels(unittest.TestC
122         self.failUnless(error_v0<TOL,"0-velocity error too large.")         self.failUnless(error_v0<TOL,"0-velocity error too large.")
123         self.failUnless(error_v1<TOL,"1-velocity error too large.")         self.failUnless(error_v1<TOL,"1-velocity error too large.")
124    
125       def test_StokesProblemCartesian_GMRES_P_0(self):
126           ETA=1.
127           P1=0.
128    
129           x=self.domain.getX()
130           F=-P1*x[1]*[1.,0]+(2*ETA-P1*x[0])*[0.,1.]
131           mask=whereZero(x[0])    * [1.,1.] \
132                  +whereZero(x[0]-1)  * [1.,1.] \
133                  +whereZero(x[1])    * [1.,0.] \
134                  +whereZero(x[1]-1)  * [1.,1.]
135          
136           sp=StokesProblemCartesian(self.domain)
137          
138           sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
139           u0=(1-x[0])*x[0]*[0.,1.]
140           p0=Scalar(P1,ReducedSolution(self.domain))
141           u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="GMRES",iter_restart=20)
142          
143           error_v0=Lsup(u[0]-u0[0])
144           error_v1=Lsup(u[1]-u0[1])/0.25
145           zz=P1*x[0]*x[1]-p
146           error_p=Lsup(zz-integrate(zz))
147    
148          # print error_p, error_v0,error_v1
149           self.failUnless(error_p<TOL,"pressure error too large.")
150           self.failUnless(error_v0<TOL,"0-velocity error too large.")
151           self.failUnless(error_v1<TOL,"1-velocity error too large.")
152    
153       def test_StokesProblemCartesian_GMRES_P_small(self):
154           ETA=1.
155           P1=1.
156    
157           x=self.domain.getX()
158           F=-P1*x[1]*[1.,0]+(2*ETA-P1*x[0])*[0.,1.]
159           mask=whereZero(x[0])    * [1.,1.] \
160                  +whereZero(x[0]-1)  * [1.,1.] \
161                  +whereZero(x[1])    * [1.,0.] \
162                  +whereZero(x[1]-1)  * [1.,1.]
163          
164           sp=StokesProblemCartesian(self.domain)
165          
166           sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
167           u0=(1-x[0])*x[0]*[0.,1.]
168           p0=Scalar(P1,ReducedSolution(self.domain))
169           u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="GMRES")
170          
171           error_v0=Lsup(u[0]-u0[0])
172           error_v1=Lsup(u[1]-u0[1])/0.25
173           zz=P1*x[0]*x[1]-p
174           error_p=Lsup(zz-integrate(zz))
175           # print error_p, error_v0,error_v1
176           self.failUnless(error_p<TOL,"pressure error too large.")
177           self.failUnless(error_v0<TOL,"0-velocity error too large.")
178           self.failUnless(error_v1<TOL,"1-velocity error too large.")
179    
180       def test_StokesProblemCartesian_GMRES_P_large(self):
181           ETA=1.
182           P1=1000.
183    
184           x=self.domain.getX()
185           F=-P1*x[1]*[1.,0]+(2*ETA-P1*x[0])*[0.,1.]
186           mask=whereZero(x[0])    * [1.,1.] \
187                  +whereZero(x[0]-1)  * [1.,1.] \
188                  +whereZero(x[1])    * [1.,0.] \
189                  +whereZero(x[1]-1)  * [1.,1.]
190          
191           sp=StokesProblemCartesian(self.domain)
192          
193           sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
194           u0=(1-x[0])*x[0]*[0.,1.]
195           p0=Scalar(P1,ReducedSolution(self.domain))
196           u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="GMRES")
197          
198           error_v0=Lsup(u[0]-u0[0])
199           error_v1=Lsup(u[1]-u0[1])/0.25
200           zz=P1*x[0]*x[1]-p
201           error_p=Lsup(zz-integrate(zz))/P1
202           # print error_p, error_v0,error_v1
203           self.failUnless(error_p<TOL,"pressure error too large.")
204           self.failUnless(error_v0<TOL,"0-velocity error too large.")
205           self.failUnless(error_v1<TOL,"1-velocity error too large.")
206    
207       def test_StokesProblemCartesian_MINRES_P_0(self):
208           ETA=1.
209           P1=0.
210    
211           x=self.domain.getX()
212           F=-P1*x[1]*[1.,0]+(2*ETA-P1*x[0])*[0.,1.]
213           mask=whereZero(x[0])    * [1.,1.] \
214                  +whereZero(x[0]-1)  * [1.,1.] \
215                  +whereZero(x[1])    * [1.,0.] \
216                  +whereZero(x[1]-1)  * [1.,1.]
217          
218           sp=StokesProblemCartesian(self.domain)
219          
220           sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
221           u0=(1-x[0])*x[0]*[0.,1.]
222           p0=Scalar(P1,ReducedSolution(self.domain))
223           u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="MINRES")
224          
225           error_v0=Lsup(u[0]-u0[0])
226           error_v1=Lsup(u[1]-u0[1])/0.25
227           zz=P1*x[0]*x[1]-p
228           error_p=Lsup(zz-integrate(zz))
229           # print error_p, error_v0,error_v1
230           self.failUnless(error_p<TOL,"pressure error too large.")
231           self.failUnless(error_v0<TOL,"0-velocity error too large.")
232           self.failUnless(error_v1<TOL,"1-velocity error too large.")
233    
234       def test_StokesProblemCartesian_MINRES_P_small(self):
235           ETA=1.
236           P1=1.
237    
238           x=self.domain.getX()
239           F=-P1*x[1]*[1.,0]+(2*ETA-P1*x[0])*[0.,1.]
240           mask=whereZero(x[0])    * [1.,1.] \
241                  +whereZero(x[0]-1)  * [1.,1.] \
242                  +whereZero(x[1])    * [1.,0.] \
243                  +whereZero(x[1]-1)  * [1.,1.]
244          
245           sp=StokesProblemCartesian(self.domain)
246          
247           sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
248           u0=(1-x[0])*x[0]*[0.,1.]
249           p0=Scalar(P1,ReducedSolution(self.domain))
250           u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="MINRES")
251          
252           error_v0=Lsup(u[0]-u0[0])
253           error_v1=Lsup(u[1]-u0[1])/0.25
254           zz=P1*x[0]*x[1]-p
255           error_p=Lsup(zz-integrate(zz))/P1
256           # print error_p, error_v0,error_v1
257           self.failUnless(error_p<TOL,"pressure error too large.")
258           self.failUnless(error_v0<TOL,"0-velocity error too large.")
259           self.failUnless(error_v1<TOL,"1-velocity error too large.")
260    
261    #   def test_StokesProblemCartesian_MINRES_P_large(self):
262    #       ETA=1.
263    #       P1=1000.
264    #
265    #       x=self.domain.getX()
266    #       F=-P1*x[1]*[1.,0]+(2*ETA-P1*x[0])*[0.,1.]
267    #       mask=whereZero(x[0])    * [1.,1.] \
268    #              +whereZero(x[0]-1)  * [1.,1.] \
269    #              +whereZero(x[1])    * [1.,0.] \
270    #              +whereZero(x[1]-1)  * [1.,1.]
271          
272    #       sp=StokesProblemCartesian(self.domain)
273          
274    #       sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
275    #       u0=(1-x[0])*x[0]*[0.,1.]
276    #       p0=Scalar(P1,ReducedSolution(self.domain))
277    #       u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="MINRES")
278          
279    #       error_v0=Lsup(u[0]-u0[0])
280    #       error_v1=Lsup(u[1]-u0[1])/0.25
281    #       zz=P1*x[0]*x[1]-p
282    #       error_p=Lsup(zz-integrate(zz))/P1
283           # print error_p, error_v0,error_v1
284    #       self.failUnless(error_p<TOL,"pressure error too large.")
285    #       self.failUnless(error_v0<TOL,"0-velocity error too large.")
286    #       self.failUnless(error_v1<TOL,"1-velocity error too large.")
287    
288    
289    #   def test_StokesProblemCartesian_TFQMR_P_0(self):
290    #       ETA=1.
291    #       P1=0.
292    
293    #       x=self.domain.getX()
294    #       F=-P1*x[1]*[1.,0]+(2*ETA-P1*x[0])*[0.,1.]
295    #       mask=whereZero(x[0])    * [1.,1.] \
296    #              +whereZero(x[0]-1)  * [1.,1.] \
297    #              +whereZero(x[1])    * [1.,0.] \
298    #              +whereZero(x[1]-1)  * [1.,1.]
299          
300    #       sp=StokesProblemCartesian(self.domain)
301          
302    #       sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
303    #       u0=(1-x[0])*x[0]*[0.,1.]
304    #       p0=Scalar(P1,ReducedSolution(self.domain))
305    #       u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="TFQMR")
306          
307    #       error_v0=Lsup(u[0]-u0[0])
308    #       error_v1=Lsup(u[1]-u0[1])/0.25
309    #       zz=P1*x[0]*x[1]-p
310    #       error_p=Lsup(zz-integrate(zz))
311           # print error_p, error_v0,error_v1
312    #       self.failUnless(error_p<TOL,"pressure error too large.")
313    #       self.failUnless(error_v0<TOL,"0-velocity error too large.")
314    #       self.failUnless(error_v1<TOL,"1-velocity error too large.")
315    
316       def test_StokesProblemCartesian_TFQMR_P_small(self):
317           ETA=1.
318           P1=1.
319    
320           x=self.domain.getX()
321           F=-P1*x[1]*[1.,0]+(2*ETA-P1*x[0])*[0.,1.]
322           mask=whereZero(x[0])    * [1.,1.] \
323                  +whereZero(x[0]-1)  * [1.,1.] \
324                  +whereZero(x[1])    * [1.,0.] \
325                  +whereZero(x[1]-1)  * [1.,1.]
326          
327           sp=StokesProblemCartesian(self.domain)
328          
329           sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
330           u0=(1-x[0])*x[0]*[0.,1.]
331           p0=Scalar(P1,ReducedSolution(self.domain))
332           u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="TFQMR")
333          
334           error_v0=Lsup(u[0]-u0[0])
335           error_v1=Lsup(u[1]-u0[1])/0.25
336           zz=P1*x[0]*x[1]-p
337           error_p=Lsup(zz-integrate(zz))/P1
338           # print error_p, error_v0,error_v1
339           self.failUnless(error_p<TOL,"pressure error too large.")
340           self.failUnless(error_v0<TOL,"0-velocity error too large.")
341           self.failUnless(error_v1<TOL,"1-velocity error too large.")
342    
343       def test_StokesProblemCartesian_TFQMR_P_large(self):
344           ETA=1.
345           P1=1000.
346    
347           x=self.domain.getX()
348           F=-P1*x[1]*[1.,0]+(2*ETA-P1*x[0])*[0.,1.]
349           mask=whereZero(x[0])    * [1.,1.] \
350                  +whereZero(x[0]-1)  * [1.,1.] \
351                  +whereZero(x[1])    * [1.,0.] \
352                  +whereZero(x[1]-1)  * [1.,1.]
353          
354           sp=StokesProblemCartesian(self.domain)
355          
356           sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
357           u0=(1-x[0])*x[0]*[0.,1.]
358           p0=Scalar(P1,ReducedSolution(self.domain))
359           u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="TFQMR")
360          
361           error_v0=Lsup(u[0]-u0[0])
362           error_v1=Lsup(u[1]-u0[1])/0.25
363           zz=P1*x[0]*x[1]-p
364           error_p=Lsup(zz-integrate(zz))/P1
365           # print error_p, error_v0,error_v1
366           self.failUnless(error_p<TOL,"pressure error too large.")
367           self.failUnless(error_v0<TOL,"0-velocity error too large.")
368           self.failUnless(error_v1<TOL,"1-velocity error too large.")
369    
370    
371    
372  class Test_Simple3DModels(unittest.TestCase):  class Test_Simple3DModels(unittest.TestCase):
373     def setUp(self):     def setUp(self):
374         self.domain=Brick(NE,NE,NE,order=2,useFullElementOrder=True)         self.domain=Brick(NE,NE,NE,order=2,useFullElementOrder=True)
375     def tearDown(self):     def tearDown(self):
376         del self.domain         del self.domain
377     def test_StokesProblemCartesian_P_0(self):     def test_StokesProblemCartesian_PCG_P_0(self):
378           ETA=1.
379           P1=0.
380    
381           x=self.domain.getX()
382           F=-P1*x[1]*x[2]*[1.,0.,0.]-P1*x[0]*x[2]*[0.,1.,0.]+(2*ETA*((1-x[0])*x[0]+(1-x[1])*x[1])-P1*x[0]*x[1])*[0.,0.,1.]
383           x=self.domain.getX()
384           mask=whereZero(x[0])    * [1.,1.,1.] \
385                  +whereZero(x[0]-1)  * [1.,1.,1.] \
386                  +whereZero(x[1])    * [1.,1.,1.] \
387                  +whereZero(x[1]-1)  * [1.,1.,1.] \
388                  +whereZero(x[2])    * [1.,1.,0.] \
389                  +whereZero(x[2]-1)  * [1.,1.,1.]
390          
391          
392           sp=StokesProblemCartesian(self.domain)
393          
394           sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
395           u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.]
396           p0=Scalar(P1,ReducedSolution(self.domain))
397           u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="PCG")
398          
399           error_v0=Lsup(u[0]-u0[0])
400           error_v1=Lsup(u[1]-u0[1])
401           error_v2=Lsup(u[2]-u0[2])/0.25**2
402           zz=P1*x[0]*x[1]*x[2]-p
403           error_p=Lsup(zz-integrate(zz))
404           # print error_p, error_v0,error_v1,error_v2
405           self.failUnless(error_p<TOL,"pressure error too large.")
406           self.failUnless(error_v0<TOL,"0-velocity error too large.")
407           self.failUnless(error_v1<TOL,"1-velocity error too large.")
408           self.failUnless(error_v2<TOL,"2-velocity error too large.")
409    
410       def test_StokesProblemCartesian_PCG_P_small(self):
411           ETA=1.
412           P1=1.
413    
414           x=self.domain.getX()
415           F=-P1*x[1]*x[2]*[1.,0.,0.]-P1*x[0]*x[2]*[0.,1.,0.]+(2*ETA*((1-x[0])*x[0]+(1-x[1])*x[1])-P1*x[0]*x[1])*[0.,0.,1.]
416           mask=whereZero(x[0])    * [1.,1.,1.] \
417                  +whereZero(x[0]-1)  * [1.,1.,1.] \
418                  +whereZero(x[1])    * [1.,1.,1.] \
419                  +whereZero(x[1]-1)  * [1.,1.,1.] \
420                  +whereZero(x[2])    * [1.,1.,0.] \
421                  +whereZero(x[2]-1)  * [1.,1.,1.]
422          
423          
424           sp=StokesProblemCartesian(self.domain)
425          
426           sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
427           u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.]
428           p0=Scalar(P1,ReducedSolution(self.domain))
429           u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="PCG")
430          
431           error_v0=Lsup(u[0]-u0[0])
432           error_v1=Lsup(u[1]-u0[1])
433           error_v2=Lsup(u[2]-u0[2])/0.25**2
434           zz=P1*x[0]*x[1]*x[2]-p
435           error_p=Lsup(zz-integrate(zz))/P1
436           # print error_p, error_v0,error_v1,error_v2
437           self.failUnless(error_p<TOL,"pressure error too large.")
438           self.failUnless(error_v0<TOL,"0-velocity error too large.")
439           self.failUnless(error_v1<TOL,"1-velocity error too large.")
440           self.failUnless(error_v2<TOL,"2-velocity error too large.")
441       def test_StokesProblemCartesian_PCG_P_large(self):
442           ETA=1.
443           P1=1000.
444    
445           x=self.domain.getX()
446           F=-P1*x[1]*x[2]*[1.,0.,0.]-P1*x[0]*x[2]*[0.,1.,0.]+(2*ETA*((1-x[0])*x[0]+(1-x[1])*x[1])-P1*x[0]*x[1])*[0.,0.,1.]
447           mask=whereZero(x[0])    * [1.,1.,1.] \
448                  +whereZero(x[0]-1)  * [1.,1.,1.] \
449                  +whereZero(x[1])    * [1.,1.,1.] \
450                  +whereZero(x[1]-1)  * [1.,1.,1.] \
451                  +whereZero(x[2])    * [1.,1.,0.] \
452                  +whereZero(x[2]-1)  * [1.,1.,1.]
453          
454          
455           sp=StokesProblemCartesian(self.domain)
456          
457           sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
458           u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.]
459           p0=Scalar(P1,ReducedSolution(self.domain))
460           u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="PCG")
461          
462           error_v0=Lsup(u[0]-u0[0])
463           error_v1=Lsup(u[1]-u0[1])
464           error_v2=Lsup(u[2]-u0[2])/0.25**2
465           zz=P1*x[0]*x[1]*x[2]-p
466           error_p=Lsup(zz-integrate(zz))/P1
467           # print error_p, error_v0,error_v1,error_v2
468           self.failUnless(error_p<TOL,"pressure error too large.")
469           self.failUnless(error_v0<TOL,"0-velocity error too large.")
470           self.failUnless(error_v1<TOL,"1-velocity error too large.")
471           self.failUnless(error_v2<TOL,"2-velocity error too large.")
472    
473       def test_StokesProblemCartesian_GMRES_P_0(self):
474         ETA=1.         ETA=1.
475         P1=0.         P1=0.
476    
# Line 148  class Test_Simple3DModels(unittest.TestC Line 490  class Test_Simple3DModels(unittest.TestC
490         sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)         sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
491         u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.]         u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.]
492         p0=Scalar(P1,ReducedSolution(self.domain))         p0=Scalar(P1,ReducedSolution(self.domain))
493         u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100)         u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="GMRES",iter_restart=20)
494                
495         error_v0=Lsup(u[0]-u0[0])         error_v0=Lsup(u[0]-u0[0])
496         error_v1=Lsup(u[1]-u0[1])         error_v1=Lsup(u[1]-u0[1])
# Line 160  class Test_Simple3DModels(unittest.TestC Line 502  class Test_Simple3DModels(unittest.TestC
502         self.failUnless(error_v0<TOL,"0-velocity error too large.")         self.failUnless(error_v0<TOL,"0-velocity error too large.")
503         self.failUnless(error_v1<TOL,"1-velocity error too large.")         self.failUnless(error_v1<TOL,"1-velocity error too large.")
504         self.failUnless(error_v2<TOL,"2-velocity error too large.")         self.failUnless(error_v2<TOL,"2-velocity error too large.")
505     def test_StokesProblemCartesian_P_small(self):     def test_StokesProblemCartesian_GMRES_P_small(self):
506         ETA=1.         ETA=1.
507         P1=1.         P1=1.
508    
# Line 179  class Test_Simple3DModels(unittest.TestC Line 521  class Test_Simple3DModels(unittest.TestC
521         sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)         sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
522         u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.]         u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.]
523         p0=Scalar(P1,ReducedSolution(self.domain))         p0=Scalar(P1,ReducedSolution(self.domain))
524         u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100)         u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="GMRES")
525                
526         error_v0=Lsup(u[0]-u0[0])         error_v0=Lsup(u[0]-u0[0])
527         error_v1=Lsup(u[1]-u0[1])         error_v1=Lsup(u[1]-u0[1])
# Line 191  class Test_Simple3DModels(unittest.TestC Line 533  class Test_Simple3DModels(unittest.TestC
533         self.failUnless(error_v0<TOL,"0-velocity error too large.")         self.failUnless(error_v0<TOL,"0-velocity error too large.")
534         self.failUnless(error_v1<TOL,"1-velocity error too large.")         self.failUnless(error_v1<TOL,"1-velocity error too large.")
535         self.failUnless(error_v2<TOL,"2-velocity error too large.")         self.failUnless(error_v2<TOL,"2-velocity error too large.")
536     def test_StokesProblemCartesian_P_large(self):     def test_StokesProblemCartesian_GMRES_P_large(self):
537         ETA=1.         ETA=1.
538         P1=1000.         P1=1000.
539    
# Line 210  class Test_Simple3DModels(unittest.TestC Line 552  class Test_Simple3DModels(unittest.TestC
552         sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)         sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
553         u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.]         u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.]
554         p0=Scalar(P1,ReducedSolution(self.domain))         p0=Scalar(P1,ReducedSolution(self.domain))
555         u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100)         u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="GMRES")
556                
557         error_v0=Lsup(u[0]-u0[0])         error_v0=Lsup(u[0]-u0[0])
558         error_v1=Lsup(u[1]-u0[1])         error_v1=Lsup(u[1]-u0[1])
# Line 223  class Test_Simple3DModels(unittest.TestC Line 565  class Test_Simple3DModels(unittest.TestC
565         self.failUnless(error_v1<TOL,"1-velocity error too large.")         self.failUnless(error_v1<TOL,"1-velocity error too large.")
566         self.failUnless(error_v2<TOL,"2-velocity error too large.")         self.failUnless(error_v2<TOL,"2-velocity error too large.")
567    
568       def test_StokesProblemCartesian_MINRES_P_0(self):
569           ETA=1.
570           P1=0.
571    
572           x=self.domain.getX()
573           F=-P1*x[1]*x[2]*[1.,0.,0.]-P1*x[0]*x[2]*[0.,1.,0.]+(2*ETA*((1-x[0])*x[0]+(1-x[1])*x[1])-P1*x[0]*x[1])*[0.,0.,1.]
574           x=self.domain.getX()
575           mask=whereZero(x[0])    * [1.,1.,1.] \
576                  +whereZero(x[0]-1)  * [1.,1.,1.] \
577                  +whereZero(x[1])    * [1.,1.,1.] \
578                  +whereZero(x[1]-1)  * [1.,1.,1.] \
579                  +whereZero(x[2])    * [1.,1.,0.] \
580                  +whereZero(x[2]-1)  * [1.,1.,1.]
581          
582          
583           sp=StokesProblemCartesian(self.domain)
584          
585           sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
586           u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.]
587           p0=Scalar(P1,ReducedSolution(self.domain))
588           u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="MINRES")
589          
590           error_v0=Lsup(u[0]-u0[0])
591           error_v1=Lsup(u[1]-u0[1])
592           error_v2=Lsup(u[2]-u0[2])/0.25**2
593           zz=P1*x[0]*x[1]*x[2]-p
594           error_p=Lsup(zz-integrate(zz))
595           # print error_p, error_v0,error_v1,error_v2
596           self.failUnless(error_p<TOL,"pressure error too large.")
597           self.failUnless(error_v0<TOL,"0-velocity error too large.")
598           self.failUnless(error_v1<TOL,"1-velocity error too large.")
599           self.failUnless(error_v2<TOL,"2-velocity error too large.")
600    
601       def test_StokesProblemCartesian_MINRES_P_small(self):
602           ETA=1.
603           P1=1.
604    
605           x=self.domain.getX()
606           F=-P1*x[1]*x[2]*[1.,0.,0.]-P1*x[0]*x[2]*[0.,1.,0.]+(2*ETA*((1-x[0])*x[0]+(1-x[1])*x[1])-P1*x[0]*x[1])*[0.,0.,1.]
607           mask=whereZero(x[0])    * [1.,1.,1.] \
608                  +whereZero(x[0]-1)  * [1.,1.,1.] \
609                  +whereZero(x[1])    * [1.,1.,1.] \
610                  +whereZero(x[1]-1)  * [1.,1.,1.] \
611                  +whereZero(x[2])    * [1.,1.,0.] \
612                  +whereZero(x[2]-1)  * [1.,1.,1.]
613          
614          
615           sp=StokesProblemCartesian(self.domain)
616          
617           sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
618           u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.]
619           p0=Scalar(P1,ReducedSolution(self.domain))
620           u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="MINRES")
621          
622           error_v0=Lsup(u[0]-u0[0])
623           error_v1=Lsup(u[1]-u0[1])
624           error_v2=Lsup(u[2]-u0[2])/0.25**2
625           zz=P1*x[0]*x[1]*x[2]-p
626           error_p=Lsup(zz-integrate(zz))/P1
627           # print error_p, error_v0,error_v1,error_v2
628           self.failUnless(error_p<TOL,"pressure error too large.")
629           self.failUnless(error_v0<TOL,"0-velocity error too large.")
630           self.failUnless(error_v1<TOL,"1-velocity error too large.")
631           self.failUnless(error_v2<TOL,"2-velocity error too large.")
632    
633    #   def test_StokesProblemCartesian_MINRES_P_large(self):
634    #       ETA=1.
635    #       P1=1000.
636    
637    #       x=self.domain.getX()
638    #       F=-P1*x[1]*x[2]*[1.,0.,0.]-P1*x[0]*x[2]*[0.,1.,0.]+(2*ETA*((1-x[0])*x[0]+(1-x[1])*x[1])-P1*x[0]*x[1])*[0.,0.,1.]
639    #       mask=whereZero(x[0])    * [1.,1.,1.] \
640    #              +whereZero(x[0]-1)  * [1.,1.,1.] \
641    #              +whereZero(x[1])    * [1.,1.,1.] \
642    #              +whereZero(x[1]-1)  * [1.,1.,1.] \
643    #              +whereZero(x[2])    * [1.,1.,0.] \
644    #              +whereZero(x[2]-1)  * [1.,1.,1.]
645          
646          
647    #       sp=StokesProblemCartesian(self.domain)
648          
649    #       sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
650    #       u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.]
651    #       p0=Scalar(P1,ReducedSolution(self.domain))
652    #       u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="MINRES")
653          
654    #       error_v0=Lsup(u[0]-u0[0])
655    #       error_v1=Lsup(u[1]-u0[1])
656    #       error_v2=Lsup(u[2]-u0[2])/0.25**2
657    #       zz=P1*x[0]*x[1]*x[2]-p
658    #       error_p=Lsup(zz-integrate(zz))/P1
659           # print error_p, error_v0,error_v1,error_v2
660    #       self.failUnless(error_p<TOL,"pressure error too large.")
661    #       self.failUnless(error_v0<TOL,"0-velocity error too large.")
662    #       self.failUnless(error_v1<TOL,"1-velocity error too large.")
663    #       self.failUnless(error_v2<TOL,"2-velocity error too large.")
664    
665    #   def test_StokesProblemCartesian_TFQMR_P_0(self):
666    #       ETA=1.
667    #       P1=0.
668    
669    #       x=self.domain.getX()
670    #       F=-P1*x[1]*x[2]*[1.,0.,0.]-P1*x[0]*x[2]*[0.,1.,0.]+(2*ETA*((1-x[0])*x[0]+(1-x[1])*x[1])-P1*x[0]*x[1])*[0.,0.,1.]
671    #       x=self.domain.getX()
672    #       mask=whereZero(x[0])    * [1.,1.,1.] \
673    #              +whereZero(x[0]-1)  * [1.,1.,1.] \
674    #              +whereZero(x[1])    * [1.,1.,1.] \
675    #              +whereZero(x[1]-1)  * [1.,1.,1.] \
676    #              +whereZero(x[2])    * [1.,1.,0.] \
677    #              +whereZero(x[2]-1)  * [1.,1.,1.]
678          
679          
680    #       sp=StokesProblemCartesian(self.domain)
681          
682    #       sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
683    #       u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.]
684    #       p0=Scalar(P1,ReducedSolution(self.domain))
685    #       u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="TFQMR")
686          
687    #       error_v0=Lsup(u[0]-u0[0])
688    #       error_v1=Lsup(u[1]-u0[1])
689    #       error_v2=Lsup(u[2]-u0[2])/0.25**2
690    #       zz=P1*x[0]*x[1]*x[2]-p
691    #       error_p=Lsup(zz-integrate(zz))
692           # print error_p, error_v0,error_v1,error_v2
693    #       self.failUnless(error_p<TOL,"pressure error too large.")
694    #       self.failUnless(error_v0<TOL,"0-velocity error too large.")
695    #       self.failUnless(error_v1<TOL,"1-velocity error too large.")
696    #       self.failUnless(error_v2<TOL,"2-velocity error too large.")
697    
698       def test_StokesProblemCartesian_TFQMR_P_small(self):
699           ETA=1.
700           P1=1.
701    
702           x=self.domain.getX()
703           F=-P1*x[1]*x[2]*[1.,0.,0.]-P1*x[0]*x[2]*[0.,1.,0.]+(2*ETA*((1-x[0])*x[0]+(1-x[1])*x[1])-P1*x[0]*x[1])*[0.,0.,1.]
704           mask=whereZero(x[0])    * [1.,1.,1.] \
705                  +whereZero(x[0]-1)  * [1.,1.,1.] \
706                  +whereZero(x[1])    * [1.,1.,1.] \
707                  +whereZero(x[1]-1)  * [1.,1.,1.] \
708                  +whereZero(x[2])    * [1.,1.,0.] \
709                  +whereZero(x[2]-1)  * [1.,1.,1.]
710          
711          
712           sp=StokesProblemCartesian(self.domain)
713          
714           sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
715           u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.]
716           p0=Scalar(P1,ReducedSolution(self.domain))
717           u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="TFQMR")
718          
719           error_v0=Lsup(u[0]-u0[0])
720           error_v1=Lsup(u[1]-u0[1])
721           error_v2=Lsup(u[2]-u0[2])/0.25**2
722           zz=P1*x[0]*x[1]*x[2]-p
723           error_p=Lsup(zz-integrate(zz))/P1
724           # print error_p, error_v0,error_v1,error_v2
725           self.failUnless(error_p<TOL,"pressure error too large.")
726           self.failUnless(error_v0<TOL,"0-velocity error too large.")
727           self.failUnless(error_v1<TOL,"1-velocity error too large.")
728           self.failUnless(error_v2<TOL,"2-velocity error too large.")
729    
730       def test_StokesProblemCartesian_TFQMR_P_large(self):
731           ETA=1.
732           P1=1000.
733    
734           x=self.domain.getX()
735           F=-P1*x[1]*x[2]*[1.,0.,0.]-P1*x[0]*x[2]*[0.,1.,0.]+(2*ETA*((1-x[0])*x[0]+(1-x[1])*x[1])-P1*x[0]*x[1])*[0.,0.,1.]
736           mask=whereZero(x[0])    * [1.,1.,1.] \
737                  +whereZero(x[0]-1)  * [1.,1.,1.] \
738                  +whereZero(x[1])    * [1.,1.,1.] \
739                  +whereZero(x[1]-1)  * [1.,1.,1.] \
740                  +whereZero(x[2])    * [1.,1.,0.] \
741                  +whereZero(x[2]-1)  * [1.,1.,1.]
742          
743          
744           sp=StokesProblemCartesian(self.domain)
745          
746           sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
747           u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.]
748           p0=Scalar(P1,ReducedSolution(self.domain))
749           u,p=sp.solve(u0,p0,show_details=VERBOSE, verbose=VERBOSE,max_iter=100,solver="TFQMR")
750          
751           error_v0=Lsup(u[0]-u0[0])
752           error_v1=Lsup(u[1]-u0[1])
753           error_v2=Lsup(u[2]-u0[2])/0.25**2
754           zz=P1*x[0]*x[1]*x[2]-p
755           error_p=Lsup(zz-integrate(zz))/P1
756           # print error_p, error_v0,error_v1,error_v2
757           self.failUnless(error_p<TOL,"pressure error too large.")
758           self.failUnless(error_v0<TOL,"0-velocity error too large.")
759           self.failUnless(error_v1<TOL,"1-velocity error too large.")
760           self.failUnless(error_v2<TOL,"2-velocity error too large.")
761    
762    
763  if __name__ == '__main__':  if __name__ == '__main__':
764     suite = unittest.TestSuite()     suite = unittest.TestSuite()
765     suite.addTest(unittest.makeSuite(Test_Simple2DModels))     suite.addTest(unittest.makeSuite(Test_Simple2DModels))
766     suite.addTest(unittest.makeSuite(Test_Simple3DModels))     # suite.addTest(unittest.makeSuite(Test_Simple3DModels))
767     s=unittest.TextTestRunner(verbosity=2).run(suite)     s=unittest.TextTestRunner(verbosity=2).run(suite)
768     if not s.wasSuccessful(): sys.exit(1)     if not s.wasSuccessful(): sys.exit(1)
769    

Legend:
Removed from v.1414  
changed lines
  Added in v.1562

  ViewVC Help
Powered by ViewVC 1.1.26