/[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 155 - (show annotations)
Wed Nov 9 02:02:19 2005 UTC (13 years, 11 months ago) by jgs
File MIME type: text/x-python
File size: 3886 byte(s)
move all directories from trunk/esys2 into trunk and remove esys2

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(190,190,1)
28
29 print "\nGenerate mesh: finley.Rectangle(151,151,1)=>"
30 # mydomain=finley.Rectangle(151,151,1)
31 # mydomain=finley.Rectangle(128,128,1)
32
33 print "\nSetup domain and functions"
34 print "--------------------------"
35
36 print "e=Function(mydomain):"
37 e=Function(mydomain)
38
39 print "n=ContinuousFunction(mydomain):"
40 n=ContinuousFunction(mydomain)
41
42 # get handles to nodes and elements 1
43
44 print "\nGet handles to nodes and elements(1)=>"
45 print "--------------------------------------"
46
47 print "u_ex=Scalar(1,n,True):"
48 u_ex=Scalar(1,n,True)
49
50 print "x=e.getX():"
51 x=e.getX()
52
53 print "norm_u_ex=u_ex.Lsup():"
54 norm_u_ex=u_ex.Lsup()
55
56 print "\nGenerate a test solution (1)"
57 print "----------------------------"
58
59 print "mypde=LinearPDE( A=[[1.,0.8],[0.4,1.]], D=alpha, Y=alpha, domain=mydomain)"
60 mypde=LinearPDE(mydomain)
61 mypde.setDebugOn()
62 mypde.setValue(A=[[1.,0.1],[0.04,1.]],D=alpha,Y=alpha)
63
64 print "mypde.checkSymmetry()"
65 print mypde.checkSymmetry()
66
67 print "\nIterative Solver (1)=>"
68 mypde.setSolverMethod(mypde.BICGSTAB)
69 u_i=mypde.getSolution(verbose=True,iter_max=3000,preconditioner=mypde.ILU0)
70
71 print "\nDirect Solver (1)=>"
72 mypde.setSolverMethod(mypde.DIRECT)
73 u_d=mypde.getSolution(verbose=True)
74
75 print "\n***************************************************************"
76 error=u_ex-u_d
77 error_norm=error.Lsup()/norm_u_ex
78 print "norm of the error for direct solver is : ",error_norm
79 if error_norm > error_tol:
80 print "### error norm exceeded maximum tolerance ###"
81 sys.exit(1)
82 error=u_ex-u_i
83 error_norm=error.Lsup()/norm_u_ex
84 print "norm of the error for iterative solver is: ",error_norm
85 if error_norm > error_tol:
86 print "### error norm exceeded maximum tolerance ###"
87 sys.exit(1)
88 print "***************************************************************"
89
90 # get handles to nodes and elements 2
91
92 print "\nGet handles to nodes and elements(2)=>"
93 print "--------------------------------------"
94
95 print "x=n.getX():"
96 x=n.getX()
97
98 print "msk=x[0].whereZero()+(x[0]-1.).whereZero()"
99 msk=x[0].whereZero()+(x[0]-1.).whereZero()
100
101 print "mypde=LinearPDE(A=[[1.,0.],[0.,1.]],q=msk,r=u_ex)"
102 mypde=LinearPDE(mydomain)
103 mypde.setDebugOn()
104 mypde.setValue(A=[[1.,0.],[0.,1.]],q=msk,r=u_ex)
105
106 print "mypde.checkSymmetry()"
107 print mypde.checkSymmetry()
108
109 # generate a test solution 2
110
111 print "\nGenerate a test solution (2)"
112 print "----------------------------"
113
114 print "\nDirect Solver (2)=>"
115
116 #mypde.setSymmetryOn()
117 mypde.setTolerance(1.e-13)
118
119 # mypde.setSymmetryOn() : is not woking yet!
120 mypde.setSolverMethod(mypde.DIRECT)
121 u_d=mypde.getSolution(verbose=True)
122
123 print "\nIterative Solver (2)=>"
124
125 mypde.setSymmetryOn()
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=error.Lsup()/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=error.Lsup()/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