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

Annotation of /trunk/finley/test/python/RecTest.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2549 - (hide annotations)
Mon Jul 20 06:43:47 2009 UTC (10 years, 3 months ago) by jfenwick
File MIME type: text/x-python
File size: 3629 byte(s)
Remainder of copyright date fixes
1 ksteube 1809
2     ########################################################
3 ksteube 1312 #
4 jfenwick 2548 # Copyright (c) 2003-2009 by University of Queensland
5 ksteube 1809 # Earth Systems Science Computational Center (ESSCC)
6     # http://www.uq.edu.au/esscc
7 ksteube 1312 #
8 ksteube 1809 # Primary Business: Queensland, Australia
9     # Licensed under the Open Software License version 3.0
10     # http://www.opensource.org/licenses/osl-3.0.php
11 ksteube 1312 #
12 ksteube 1809 ########################################################
13 jgs 153
14 jfenwick 2549 __copyright__="""Copyright (c) 2003-2009 by University of Queensland
15 ksteube 1809 Earth Systems Science Computational Center (ESSCC)
16     http://www.uq.edu.au/esscc
17     Primary Business: Queensland, Australia"""
18 elspeth 617 __license__="""Licensed under the Open Software License version 3.0
19 ksteube 1809 http://www.opensource.org/licenses/osl-3.0.php"""
20 jfenwick 2344 __url__="https://launchpad.net/escript-finley"
21 ksteube 1809
22 jgs 82 import sys
23     import os
24     import unittest
25    
26 jgs 104 from esys.escript import *
27 jgs 149 from esys.escript.linearPDEs import *
28     from esys import finley
29 jgs 104
30 jfenwick 2455 import numpy
31 jgs 153
32 jfenwick 2455 Pi=numpy.pi
33 jgs 153 numElements=15
34     sml=0.1
35 jgs 82
36     #
37     # this test the assemblage of problems with periodic boundary conditions:
38     #
39 jgs 153
40 jgs 82 #
41     # test solution is u_ex=sin(2*Pi*n*x0)*...*sin(2*Pi*n*x_dim)
42     #
43     def TheTest(msh,constraints,reduce):
44     n=ContinuousFunction(msh)
45     x=n.getX()
46     # set a test problem:
47     u_ex=Scalar(1,what=n)
48     for i in range(msh.getDim()):
49     u_ex*=sin(2*Pi*x[i])
50 jgs 104 mypde=LinearPDE(msh)
51 jgs 102 mypde.setSymmetryOn()
52     mypde.setDebugOn()
53     mypde.setReducedOrderTo(reduce)
54 jfenwick 2455 mypde.setValue(A=numpy.identity(msh.getDim()),D=sml,Y=(sml+4*Pi**2*msh.getDim())*u_ex,q=constraints,r=u_ex)
55 jgs 82 return Lsup(mypde.getSolution()-u_ex)/Lsup(u_ex)
56    
57     max_error=0
58 jgs 102 max_text="none"
59 jgs 82 for onElements in [False,True]:
60     if onElements==True:
61     onElmtext=", with elements on faces"
62     else:
63     onElmtext=""
64     for order in [1,2]:
65     for reduce in [False,True]:
66     if reduce==True:
67     redtext=",reduced"
68     else:
69     redtext=""
70     for dim in [2,3]:
71     if dim==2:
72 jgs 102 for i0 in [True,True]:
73     for i1 in [True,True]:
74 jgs 153 msh=finley.Rectangle(numElements,numElements,order,periodic0=i0,periodic1=i1,useElementsOnFace=onElements)
75 jgs 82 n=ContinuousFunction(msh)
76     x=n.getX()
77     c=Scalar(0,what=n)
78     if i0==False:
79 gross 301 c+=whereZero(x[0])+whereZero(x[0]-1.)
80 jgs 82 if i1==False:
81 gross 301 c+=whereZero(x[1])+whereZero(x[1]-1.)
82 jgs 82 error=TheTest(msh,c,reduce)
83     text="Rectangle order = %d%s%s, periodic0= %d, periodic1= %d: error= %f"%(order,redtext,onElmtext,i0,i1,error)
84     print "@@ ",text
85     if error>max_error:
86     max_error=error
87     max_text=text
88     elif dim==3:
89     for i0 in [True,False]:
90     for i1 in [True,False]:
91     for i2 in [True,False]:
92 jgs 153 msh=finley.Brick(numElements,numElements,numElements,order,periodic0=i0,periodic1=i1,periodic2=i2,useElementsOnFace=onElements)
93 jgs 82 n=ContinuousFunction(msh)
94     x=n.getX()
95     c=Scalar(0,what=n)
96     if i0==False:
97 gross 301 c+=whereZero(x[0])+whereZero(x[0]-1.)
98 jgs 82 if i1==False:
99 gross 301 c+=whereZero(x[1])+whereZero(x[1]-1.)
100 jgs 82 if i2==False:
101 gross 301 c+=whereZero(x[2])+whereZero(x[2]-1.)
102 jgs 82 error=TheTest(msh,c,reduce)
103     text="Brick order = %d%s%s, periodic0= %d, periodic1= %d, periodic2= %d: error= %f"%(order,redtext,onElmtext,i0,i1,i2,error)
104     print "@@ ",text
105     if error>max_error:
106     max_error=error
107     max_text=text
108    
109     print "@@@@ maximum error for :",max_text

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26