test/python/SimpleSolve.py

# $Id$

import sys
import os
import unittest

esys_root=os.getenv('ESYS_ROOT')
sys.path.append(esys_root+'/finley/lib')
sys.path.append(esys_root+'/escript/lib')
sys.path.append(esys_root+'/escript/py_src')

from escript import *
from util import *
from linearPDE import *
import finley

print "\nSimpleSolve.py"
print "--------------"

my_options = {
          "verbose" : True,
          "reordering" : NO_REORDERING,
          "tolerance" : 1.E-8,
          "final_residual" : 0.,
          "iterative_method" : PCG,
          "preconditioner" : JACOBI,
          "iter_max" :  5000,
          "iter" : 0,
          "drop_tolerance" : 1.10,
          "drop_storage" : 2.
}

alpha=0.01

print "\nOptions: ", my_options

# generate mesh

print "\nGenerate mesh: finley.Rectangle(9,12,1)=>"
mydomain=finley.Rectangle(9,12,1)

print "\nGenerate mesh: finley.Rectangle(4,4,1)=>"
mydomain=finley.Rectangle(4,4,1)

print "\nGenerate mesh: finley.Rectangle(151,151,1)=>"
mydomain=finley.Rectangle(151,151,1)

print "\nSetup domain and functions"
print "--------------------------"

print "e=Function(mydomain):"
e=Function(mydomain)

print "n=ContinuousFunction(mydomain):"
n=ContinuousFunction(mydomain)

# get handles to nodes and elements 1

print "\nGet handles to nodes and elements(1)=>"
print "--------------------------------------"

print "u_ex=Scalar(1,n,True):"
u_ex=Scalar(1,n,True)

print "x=e.getX():"
x=e.getX()

print "norm_u_ex=u_ex.Lsup():"
norm_u_ex=u_ex.Lsup()

print "mypde=linearPDE( A=[[1.,0.7],[0.7,1.]], D=alpha, Y=alpha, domain=mydomain)"
mypde=linearPDE(A=[[1.,0.7],[0.7,1.]],D=alpha,Y=alpha,domain=mydomain)
#mypde=linearPDE(A=[[1.,0.],[0.,1.]],D=alpha,Y=alpha,domain=mydomain)

# generate a test solution 1

print "\nGenerate a test solution (1)"
print "----------------------------"

print "\nDirect Solver (1)=>"

u_d=mypde.getSolution(iterative=False,**my_options)

print "\nIterative Solver (1)=>"

u_i=mypde.getSolution(iterative=True,**my_options)

print "\n***************************************************************"
error=u_ex-u_d
print "norm of the error for direct solver is   : ",error.Lsup()/norm_u_ex
error=u_ex-u_i
print "norm of the error for iterative solver is: ",error.Lsup()/norm_u_ex
print "***************************************************************"

# get handles to nodes and elements 2

print "\nGet handles to nodes and elements(2)=>"
print "--------------------------------------"

print "x=n.getX():"
x=n.getX()

print "msk=x[0].whereZero()+(x[0]-1.).whereZero()"
msk=x[0].whereZero()+(x[0]-1.).whereZero()

print "mypde=linearPDE(A=[[1.,0.],[0.,1.]],q=msk,r=u_ex)"
mypde=linearPDE(A=[[1.,0.],[0.,1.]],q=msk,r=u_ex)

# generate a test solution 2

print "\nGenerate a test solution (2)"
print "----------------------------"

print "\nDirect Solver (2)=>"

u_d=mypde.getSolution(iterative=False,**my_options)

print "\nIterative Solver (2)=>"

u_i=mypde.getSolution(iterative=True,**my_options)

print "\n******************************************************************"
error=u_ex-u_d
print "norm of the error for direct solver is   : ",error.Lsup()/norm_u_ex
error=u_ex-u_i
print "norm of the error for iterative solver is: ",error.Lsup()/norm_u_ex
print "******************************************************************"

print "\n-----"
print "Done."
print "-----"
1	# $Id$
2
3	import sys
4	import os
5	import unittest
6
7	esys_root=os.getenv('ESYS_ROOT')
8	sys.path.append(esys_root+'/finley/lib')
9	sys.path.append(esys_root+'/escript/lib')
10	sys.path.append(esys_root+'/escript/py_src')
11
12	from escript import *
13	from util import *
14	from linearPDE import *
15	import finley
16
17	print "\nSimpleSolve.py"
18	print "--------------"
19
20	my_options = {
21	"verbose" : True,
22	"reordering" : NO_REORDERING,
23	"tolerance" : 1.E-8,
24	"final_residual" : 0.,
25	"iterative_method" : PCG,
26	"preconditioner" : JACOBI,
27	"iter_max" : 5000,
28	"iter" : 0,
29	"drop_tolerance" : 1.10,
30	"drop_storage" : 2.
31	}
32
33	alpha=0.01
34
35	print "\nOptions: ", my_options
36
37	# generate mesh
38
39	print "\nGenerate mesh: finley.Rectangle(9,12,1)=>"
40	mydomain=finley.Rectangle(9,12,1)
41
42	print "\nGenerate mesh: finley.Rectangle(4,4,1)=>"
43	mydomain=finley.Rectangle(4,4,1)
44
45	print "\nGenerate mesh: finley.Rectangle(151,151,1)=>"
46	mydomain=finley.Rectangle(151,151,1)
47
48	print "\nSetup domain and functions"
49	print "--------------------------"
50
51	print "e=Function(mydomain):"
52	e=Function(mydomain)
53
54	print "n=ContinuousFunction(mydomain):"
55	n=ContinuousFunction(mydomain)
56
57	# get handles to nodes and elements 1
58
59	print "\nGet handles to nodes and elements(1)=>"
60	print "--------------------------------------"
61
62	print "u_ex=Scalar(1,n,True):"
63	u_ex=Scalar(1,n,True)
64
65	print "x=e.getX():"
66	x=e.getX()
67
68	print "norm_u_ex=u_ex.Lsup():"
69	norm_u_ex=u_ex.Lsup()
70
71	print "mypde=linearPDE( A=[[1.,0.7],[0.7,1.]], D=alpha, Y=alpha, domain=mydomain)"
72	mypde=linearPDE(A=[[1.,0.7],[0.7,1.]],D=alpha,Y=alpha,domain=mydomain)
73	#mypde=linearPDE(A=[[1.,0.],[0.,1.]],D=alpha,Y=alpha,domain=mydomain)
74
75	# generate a test solution 1
76
77	print "\nGenerate a test solution (1)"
78	print "----------------------------"
79
80	print "\nDirect Solver (1)=>"
81
82	u_d=mypde.getSolution(iterative=False,**my_options)
83
84	print "\nIterative Solver (1)=>"
85
86	u_i=mypde.getSolution(iterative=True,**my_options)
87
88	print "\n***************************************************************"
89	error=u_ex-u_d
90	print "norm of the error for direct solver is : ",error.Lsup()/norm_u_ex
91	error=u_ex-u_i
92	print "norm of the error for iterative solver is: ",error.Lsup()/norm_u_ex
93	print "***************************************************************"
94
95	# get handles to nodes and elements 2
96
97	print "\nGet handles to nodes and elements(2)=>"
98	print "--------------------------------------"
99
100	print "x=n.getX():"
101	x=n.getX()
102
103	print "msk=x[0].whereZero()+(x[0]-1.).whereZero()"
104	msk=x[0].whereZero()+(x[0]-1.).whereZero()
105
106	print "mypde=linearPDE(A=[[1.,0.],[0.,1.]],q=msk,r=u_ex)"
107	mypde=linearPDE(A=[[1.,0.],[0.,1.]],q=msk,r=u_ex)
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	u_d=mypde.getSolution(iterative=False,**my_options)
117
118	print "\nIterative Solver (2)=>"
119
120	u_i=mypde.getSolution(iterative=True,**my_options)
121
122	print "\n******************************************************************"
123	error=u_ex-u_d
124	print "norm of the error for direct solver is : ",error.Lsup()/norm_u_ex
125	error=u_ex-u_i
126	print "norm of the error for iterative solver is: ",error.Lsup()/norm_u_ex
127	print "******************************************************************"
128
129	print "\n-----"
130	print "Done."
131	print "-----"
Name	Value
svn:eol-style	native
svn:keywords	Author Date Id Revision