/[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 617 - (hide annotations)
Wed Mar 22 02:58:17 2006 UTC (14 years ago) by elspeth
File MIME type: text/x-python
File size: 3174 byte(s)
More copyright.

1 jgs 82 # $Id$
2 jgs 153
3 elspeth 617 __copyright__=""" Copyright (c) 2006 by ACcESS MNRF
4     http://www.access.edu.au
5     Primary Business: Queensland, Australia"""
6     __license__="""Licensed under the Open Software License version 3.0
7     http://www.opensource.org/licenses/osl-3.0.php"""
8 jgs 82 import sys
9     import os
10     import unittest
11    
12 jgs 104 from esys.escript import *
13 jgs 149 from esys.escript.linearPDEs import *
14     from esys import finley
15 jgs 104
16 jgs 102 import numarray
17 jgs 153
18 jgs 149 Pi=numarray.pi
19 jgs 153 numElements=15
20     sml=0.1
21 jgs 82
22     #
23     # this test the assemblage of problems with periodic boundary conditions:
24     #
25 jgs 153
26 jgs 82 #
27     # test solution is u_ex=sin(2*Pi*n*x0)*...*sin(2*Pi*n*x_dim)
28     #
29     def TheTest(msh,constraints,reduce):
30     n=ContinuousFunction(msh)
31     x=n.getX()
32     # set a test problem:
33     u_ex=Scalar(1,what=n)
34     for i in range(msh.getDim()):
35     u_ex*=sin(2*Pi*x[i])
36 jgs 104 mypde=LinearPDE(msh)
37 jgs 102 mypde.setSymmetryOn()
38     mypde.setDebugOn()
39     mypde.setReducedOrderTo(reduce)
40 gross 301 mypde.setValue(A=numarray.identity(msh.getDim()),D=sml,Y=(sml+4*Pi**2*msh.getDim())*u_ex,q=constraints,r=u_ex)
41 jgs 82 return Lsup(mypde.getSolution()-u_ex)/Lsup(u_ex)
42    
43     max_error=0
44 jgs 102 max_text="none"
45 jgs 82 for onElements in [False,True]:
46     if onElements==True:
47     onElmtext=", with elements on faces"
48     else:
49     onElmtext=""
50     for order in [1,2]:
51     for reduce in [False,True]:
52     if reduce==True:
53     redtext=",reduced"
54     else:
55     redtext=""
56     for dim in [2,3]:
57     if dim==2:
58 jgs 102 for i0 in [True,True]:
59     for i1 in [True,True]:
60 jgs 153 msh=finley.Rectangle(numElements,numElements,order,periodic0=i0,periodic1=i1,useElementsOnFace=onElements)
61 jgs 82 n=ContinuousFunction(msh)
62     x=n.getX()
63     c=Scalar(0,what=n)
64     if i0==False:
65 gross 301 c+=whereZero(x[0])+whereZero(x[0]-1.)
66 jgs 82 if i1==False:
67 gross 301 c+=whereZero(x[1])+whereZero(x[1]-1.)
68 jgs 82 error=TheTest(msh,c,reduce)
69     text="Rectangle order = %d%s%s, periodic0= %d, periodic1= %d: error= %f"%(order,redtext,onElmtext,i0,i1,error)
70     print "@@ ",text
71     if error>max_error:
72     max_error=error
73     max_text=text
74     elif dim==3:
75     for i0 in [True,False]:
76     for i1 in [True,False]:
77     for i2 in [True,False]:
78 jgs 153 msh=finley.Brick(numElements,numElements,numElements,order,periodic0=i0,periodic1=i1,periodic2=i2,useElementsOnFace=onElements)
79 jgs 82 n=ContinuousFunction(msh)
80     x=n.getX()
81     c=Scalar(0,what=n)
82     if i0==False:
83 gross 301 c+=whereZero(x[0])+whereZero(x[0]-1.)
84 jgs 82 if i1==False:
85 gross 301 c+=whereZero(x[1])+whereZero(x[1]-1.)
86 jgs 82 if i2==False:
87 gross 301 c+=whereZero(x[2])+whereZero(x[2]-1.)
88 jgs 82 error=TheTest(msh,c,reduce)
89     text="Brick order = %d%s%s, periodic0= %d, periodic1= %d, periodic2= %d: error= %f"%(order,redtext,onElmtext,i0,i1,i2,error)
90     print "@@ ",text
91     if error>max_error:
92     max_error=error
93     max_text=text
94    
95     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