/[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 1517 by artak, Wed Apr 16 00:15:44 2008 UTC revision 1518 by artak, Fri Apr 18 05:47:52 2008 UTC
# Line 144  class Test_Simple2DModels(unittest.TestC Line 144  class Test_Simple2DModels(unittest.TestC
144         error_v1=Lsup(u[1]-u0[1])/0.25         error_v1=Lsup(u[1]-u0[1])/0.25
145         zz=P1*x[0]*x[1]-p         zz=P1*x[0]*x[1]-p
146         error_p=Lsup(zz-integrate(zz))         error_p=Lsup(zz-integrate(zz))
147         # print error_p, error_v0,error_v1  
148          # print error_p, error_v0,error_v1
149         self.failUnless(error_p<TOL,"pressure error too large.")         self.failUnless(error_p<TOL,"pressure error too large.")
150         self.failUnless(error_v0<TOL,"0-velocity error too large.")         self.failUnless(error_v0<TOL,"0-velocity error too large.")
151         self.failUnless(error_v1<TOL,"1-velocity error too large.")         self.failUnless(error_v1<TOL,"1-velocity error too large.")
# Line 203  class Test_Simple2DModels(unittest.TestC Line 204  class Test_Simple2DModels(unittest.TestC
204         self.failUnless(error_v0<TOL,"0-velocity error too large.")         self.failUnless(error_v0<TOL,"0-velocity error too large.")
205         self.failUnless(error_v1<TOL,"1-velocity error too large.")         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):
# Line 241  class Test_Simple3DModels(unittest.TestC Line 406  class Test_Simple3DModels(unittest.TestC
406         self.failUnless(error_v0<TOL,"0-velocity error too large.")         self.failUnless(error_v0<TOL,"0-velocity error too large.")
407         self.failUnless(error_v1<TOL,"1-velocity error too large.")         self.failUnless(error_v1<TOL,"1-velocity error too large.")
408         self.failUnless(error_v2<TOL,"2-velocity error too large.")         self.failUnless(error_v2<TOL,"2-velocity error too large.")
409    
410     def test_StokesProblemCartesian_PCG_P_small(self):     def test_StokesProblemCartesian_PCG_P_small(self):
411         ETA=1.         ETA=1.
412         P1=1.         P1=1.
# Line 399  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))
Legend:
Removed from v.1517  
changed lines
  Added in v.1518
  ViewVC Help
Powered by ViewVC 1.1.26