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

Contents of /trunk/finley/test/python/SimpleSolve.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 416 - (show annotations)
Wed Jan 4 05:38:35 2006 UTC (13 years, 10 months ago) by gross
File MIME type: text/x-python
File size: 3886 byte(s)
changes in Paso have been worked in
1 # $Id$
2
3 import sys
4 import os
5 import unittest
6 import time
7
8 from esys.escript import *
9 from esys.escript.linearPDEs import *
10 from esys import finley
11
12 starttime = time.clock()
13
14 print "\nSimpleSolve.py"
15 print "--------------"
16
17 alpha=0.7
18 error_tol=1.e-5
19
20 # generate mesh
21
22 # print "\nGenerate mesh: finley.Rectangle(9,12,1)=>"
23 # mydomain=finley.Rectangle(140,140)
24
25 # print "\nGenerate mesh: finley.Rectangle(4,4,1)=>"
26 mydomain=finley.Rectangle(150,10,1)
27 mydomain=finley.Rectangle(250,250,1)
28 # mydomain=finley.Rectangle(190,190,1)
29
30 print "\nGenerate mesh: finley.Rectangle(151,151,1)=>"
31 # mydomain=finley.Rectangle(151,151,1)
32 # mydomain=finley.Rectangle(128,128,1)
33
34 print "\nSetup domain and functions"
35 print "--------------------------"
36
37 print "e=Function(mydomain):"
38 e=Function(mydomain)
39
40 print "n=ContinuousFunction(mydomain):"
41 n=ContinuousFunction(mydomain)
42
43 # get handles to nodes and elements 1
44
45 print "\nGet handles to nodes and elements(1)=>"
46 print "--------------------------------------"
47
48 print "u_ex=Scalar(1,n,True):"
49 u_ex=Scalar(1,n,True)
50
51 print "x=e.getX():"
52 x=e.getX()
53
54 print "norm_u_ex=Lsup(u_ex):"
55 norm_u_ex=Lsup(u_ex)
56
57 print "\nGenerate a test solution (1)"
58 print "----------------------------"
59
60 print "mypde=LinearPDE( A=[[1.,0.8],[0.4,1.]], D=alpha, Y=alpha, domain=mydomain)"
61 mypde=LinearPDE(mydomain)
62 mypde.setDebugOn()
63 mypde.setValue(A=[[1.,0.1],[0.04,1.]],D=alpha,Y=alpha)
64
65 print "mypde.checkSymmetry()"
66 print mypde.checkSymmetry()
67
68 print "\nIterative Solver (1)=>"
69 mypde.setSolverMethod(mypde.BICGSTAB)
70 u_i=mypde.getSolution(verbose=True,iter_max=3000,preconditioner=mypde.ILU0)
71
72 print "\nDirect Solver (1)=>"
73 mypde.setSolverMethod(mypde.DIRECT)
74 u_d=mypde.getSolution(verbose=True)
75
76 print "\n***************************************************************"
77 error=u_ex-u_d
78 error_norm=Lsup(error)/norm_u_ex
79 print "norm of the error for direct solver is : ",error_norm
80 if error_norm > error_tol:
81 print "### error norm exceeded maximum tolerance ###"
82 sys.exit(1)
83 error=u_ex-u_i
84 error_norm=Lsup(error)/norm_u_ex
85 print "norm of the error for iterative solver is: ",error_norm
86 if error_norm > error_tol:
87 print "### error norm exceeded maximum tolerance ###"
88 sys.exit(1)
89 print "***************************************************************"
90
91 # get handles to nodes and elements 2
92
93 print "\nGet handles to nodes and elements(2)=>"
94 print "--------------------------------------"
95
96 print "x=n.getX():"
97 x=n.getX()
98
99 print " msk=whereZero(x[0])+whereZero(x[0]-1.)"
100 msk=whereZero(x[0])+whereZero(x[0]-1.)
101
102 print "mypde=LinearPDE(A=[[1.,0.],[0.,1.]],q=msk,r=u_ex)"
103 mypde=LinearPDE(mydomain)
104 mypde.setDebugOn()
105 mypde.setValue(A=[[1.,0.],[0.,1.]],q=msk,r=u_ex)
106
107 print "mypde.checkSymmetry()"
108 print mypde.checkSymmetry()
109
110 # generate a test solution 2
111
112 print "\nGenerate a test solution (2)"
113 print "----------------------------"
114
115 print "\nDirect Solver (2)=>"
116
117 mypde.setSymmetryOn()
118 mypde.setTolerance(1.e-13)
119
120 # mypde.setSymmetryOn() : is not woking yet!
121 mypde.setSolverMethod(mypde.DIRECT)
122 u_d=mypde.getSolution(verbose=True)
123
124 print "\nIterative Solver (2)=>"
125
126 mypde.setSolverMethod(mypde.ITERATIVE)
127 u_i=mypde.getSolution(verbose=True,iter_max=3000)
128
129 print "\n******************************************************************"
130 error=u_ex-u_d
131 error_norm=Lsup(error)/norm_u_ex
132 print "norm of the error for direct solver is : ",error_norm
133 if error_norm > error_tol:
134 print "### error norm exceeded maximum tolerance ###"
135 sys.exit(1)
136 error=u_ex-u_i
137 error_norm=Lsup(error)/norm_u_ex
138 print "norm of the error for iterative solver is: ",error_norm
139 if error_norm > error_tol:
140 print "### error norm exceeded maximum tolerance ###"
141 sys.exit(1)
142 print "******************************************************************"
143
144 print "\n-----"
145 print "Done."
146 print "-----"
147
148 stoptime = time.clock()
149 elapsed = stoptime - starttime
150 print "\nElapsed time: ", elapsed, "\n"
151
152 sys.exit(0)

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26