/[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 104 - (show annotations)
Fri Dec 17 07:43:12 2004 UTC (14 years, 4 months ago) by jgs
File MIME type: text/x-python
File size: 2906 byte(s)
*** empty log message ***

1 # $Id$
2
3 import sys
4 import os
5 import unittest
6
7 from esys.escript import *
8 from esys.linearPDEs import *
9 from esys import finley
10
11 print "\nSimpleSolve.py"
12 print "--------------"
13
14 alpha=0.01
15
16 # generate mesh
17
18 # print "\nGenerate mesh: finley.Rectangle(9,12,1)=>"
19 # mydomain=finley.Rectangle(1,1)
20
21 # print "\nGenerate mesh: finley.Rectangle(4,4,1)=>"
22 # mydomain=finley.Rectangle(4,4,1)
23
24 print "\nGenerate mesh: finley.Rectangle(151,151,1)=>"
25 mydomain=finley.Rectangle(151,151,1)
26
27 print "\nSetup domain and functions"
28 print "--------------------------"
29
30 print "e=Function(mydomain):"
31 e=Function(mydomain)
32
33 print "n=ContinuousFunction(mydomain):"
34 n=ContinuousFunction(mydomain)
35
36 # get handles to nodes and elements 1
37
38 print "\nGet handles to nodes and elements(1)=>"
39 print "--------------------------------------"
40
41 print "u_ex=Scalar(1,n,True):"
42 u_ex=Scalar(1,n,True)
43
44 print "x=e.getX():"
45 x=e.getX()
46
47 print "norm_u_ex=u_ex.Lsup():"
48 norm_u_ex=u_ex.Lsup()
49
50 print "mypde=LinearPDE( A=[[1.,0.8],[0.4,1.]], D=alpha, Y=alpha, domain=mydomain)"
51 mypde=LinearPDE(mydomain)
52 mypde.setDebugOn()
53 mypde.setValue(A=[[1.,0.8],[0.4,1.]],D=alpha,Y=alpha)
54 mypde.getOperator().saveMM("t.msh")
55
56 # generate a test solution 1
57
58 print "\nGenerate a test solution (1)"
59 print "----------------------------"
60
61 print "\nIterative Solver (1)=>"
62
63 u_i=mypde.getSolution()
64
65 print "\nDirect Solver (1)=>"
66
67 mypde.setSolverMethod(DIRECT)
68 u_d=mypde.getSolution()
69
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.setValue(A=[[1.,0.],[0.,1.]],q=msk,r=u_ex)
92 mypde.setDebugOn()
93
94 # generate a test solution 2
95
96 print "\nGenerate a test solution (2)"
97 print "----------------------------"
98
99 print "\nDirect Solver (2)=>"
100
101 # mypde.setSymmetryOn() : is not woking yet!
102 mypde.setSolverMethod(mypde.DIRECT)
103 u_d=mypde.getSolution()
104
105 print "\nIterative Solver (2)=>"
106
107 #mypde.setSymmetryOn()
108 mypde.setSolverMethod(DEFAULT_METHOD)
109 u_i=mypde.getSolution()
110
111 print "\n******************************************************************"
112 error=u_ex-u_d
113 print "norm of the error for direct solver is : ",error.Lsup()/norm_u_ex
114 error=u_ex-u_i
115 print "norm of the error for iterative solver is: ",error.Lsup()/norm_u_ex
116 print "******************************************************************"
117
118 print "\n-----"
119 print "Done."
120 print "-----"

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26