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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 141 - (show annotations)
Mon Jul 25 05:15:41 2005 UTC (14 years, 3 months ago) by jgs
File MIME type: text/x-python
File size: 3150 byte(s)
fixed to work with new python source/install directory structure

1 # $Id$
2
3 import sys
4 import os
5 import unittest
6
7 from escript.escript import *
8 from escript.linearPDEs import *
9 from finley import finley
10
11 print "\nSimpleSolve.py"
12 print "--------------"
13
14 alpha=0.025
15
16 # generate mesh
17
18 # print "\nGenerate mesh: finley.Rectangle(9,12,1)=>"
19 # mydomain=finley.Rectangle(140,140)
20
21 # print "\nGenerate mesh: finley.Rectangle(4,4,1)=>"
22 # mydomain=finley.Rectangle(5,5,1)
23
24 print "\nGenerate mesh: finley.Rectangle(151,151,1)=>"
25 mydomain=finley.Rectangle(151,151,1)
26 # mydomain=finley.Rectangle(128,128,1)
27
28 print "\nSetup domain and functions"
29 print "--------------------------"
30
31 print "e=Function(mydomain):"
32 e=Function(mydomain)
33
34 print "n=ContinuousFunction(mydomain):"
35 n=ContinuousFunction(mydomain)
36
37 # get handles to nodes and elements 1
38
39 print "\nGet handles to nodes and elements(1)=>"
40 print "--------------------------------------"
41
42 print "u_ex=Scalar(1,n,True):"
43 u_ex=Scalar(1,n,True)
44
45 print "x=e.getX():"
46 x=e.getX()
47
48 print "norm_u_ex=u_ex.Lsup():"
49 norm_u_ex=u_ex.Lsup()
50
51 print "\nGenerate a test solution (1)"
52 print "----------------------------"
53
54 print "mypde=LinearPDE( A=[[1.,0.8],[0.4,1.]], D=alpha, Y=alpha, domain=mydomain)"
55 mypde=LinearPDE(mydomain)
56 mypde.setDebugOn()
57 mypde.setValue(A=[[1.,0.8],[0.4,1.]],D=alpha,Y=alpha)
58
59 print "mypde.checkSymmetry()"
60 print mypde.checkSymmetry()
61
62 print "\nIterative Solver (1)=>"
63 # u_i=mypde.getSolution(preconditioner=ILU0,iter_max=3000)
64 u_i=mypde.getSolution(iter_max=3000)
65
66
67 print "\nDirect Solver (1)=>"
68 mypde.setSolverMethod(DIRECT)
69 u_d=mypde.getSolution()
70
71 print "\n***************************************************************"
72 error=u_ex-u_d
73 print "norm of the error for direct solver is : ",error.Lsup()/norm_u_ex
74 error=u_ex-u_i
75 print "norm of the error for iterative solver is: ",error.Lsup()/norm_u_ex
76 print "***************************************************************"
77
78 # get handles to nodes and elements 2
79
80 print "\nGet handles to nodes and elements(2)=>"
81 print "--------------------------------------"
82
83 print "x=n.getX():"
84 x=n.getX()
85
86 print "msk=x[0].whereZero()+(x[0]-1.).whereZero()"
87 msk=x[0].whereZero()+(x[0]-1.).whereZero()
88
89 print "mypde=LinearPDE(A=[[1.,0.],[0.,1.]],q=msk,r=u_ex)"
90 mypde=LinearPDE(mydomain)
91 mypde.setDebugOn()
92 mypde.setValue(A=[[1.,0.],[0.,1.]],q=msk,r=u_ex)
93
94 print "mypde.checkSymmetry()"
95 print mypde.checkSymmetry()
96
97 # generate a test solution 2
98
99 print "\nGenerate a test solution (2)"
100 print "----------------------------"
101
102 print "\nDirect Solver (2)=>"
103
104 #mypde.setSymmetryOn()
105 mypde.setTolerance(1.e-13)
106
107 # mypde.setSymmetryOn() : is not woking yet!
108 mypde.setSolverMethod(mypde.DIRECT)
109 u_d=mypde.getSolution()
110
111 print "\nIterative Solver (2)=>"
112
113 mypde.setSymmetryOn()
114 mypde.setSolverMethod(mypde.DEFAULT_METHOD)
115 u_i=mypde.getSolution(iter_max=3000)
116
117 print "\n******************************************************************"
118 error=u_ex-u_d
119 print "norm of the error for direct solver is : ",error.Lsup()/norm_u_ex
120 error=u_ex-u_i
121 print "norm of the error for iterative solver is: ",error.Lsup()/norm_u_ex
122 print "******************************************************************"
123
124 print "\n-----"
125 print "Done."
126 print "-----"

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26