← Back to branch summary

~vcs-imports/escript-finley/trunk

« back to all changes in this revision

Viewing changes to dudley/src/generateReferenceElementList.py

  • Committer: jfenwick
  • Date: 2010-10-11 01:48:14 UTC
  • Revision ID: svn-v4:77569008-7704-0410-b7a0-a92fef0b09fd:trunk:3259
Merging dudley and scons updates from branches

Show diffs side-by-side

added added

removed removed

Lines of Context:
dudley/src/generateReferenceElementList.py
 
1
#*******************************************************
 
2
#
 
3
# Copyright (c) 2003-2010 by University of Queensland
 
4
# Earth Systems Science Computational Center (ESSCC)
 
5
# http://www.uq.edu.au/esscc
 
6
#
 
7
# Primary Business: Queensland, Australia
 
8
# Licensed under the Open Software License version 3.0
 
9
# http://www.opensource.org/licenses/osl-3.0.php
 
10
#
 
11
#******************************************************/
 
12
#
 
13
#  this code generates the  Dudley_ReferenceElement_InfoList in ReferenceElements.c
 
14
#
 
15
GEOBASE = {
 
16
"Point": (1, 1, "Point", [0] ),
 
17
"Line":  (2, 2, "Point", [0,2, 2,1]),
 
18
"Tri":(4, 3,   "Line",  [0, 3, 5,  5, 4, 2,  3, 1, 4,  4, 5, 3]),
 
19
"Rec" :(4,  4, "Line", [0, 4, 8, 7,  7, 8, 6, 3,   4, 1, 5, 8,  8, 5, 2, 6] ),
 
20
"Tet": (8, 4, "Tri",  [0, 4, 6, 7,  4, 1, 5, 8,  6, 5, 2, 9,  7, 8, 9, 3,  4, 5, 6, 8,  5, 9, 8, 6,  9, 7, 6, 8,  7, 4, 6, 8] ),
 
21
"Hex": (8, 8, "Rec", [0, 8, 20, 11, 12, 21, 26, 24,  8, 1, 9, 20, 21, 13, 22, 26,  11, 20, 10, 3, 24, 26, 23, 15,  20, 9, 2, 10, 26, 22, 14, 23,  12, 21, 26, 24, 4, 16, 25, 19,  21, 13, 22, 26, 16, 5, 17, 25,  24, 26, 23, 15, 19, 25, 18, 7,  26, 22, 14, 23, 25, 17, 6, 18] )
 
22
 
 
23
}
 
24
 
 
25
RELEVANTGEO =  {
 
26
  "Point1" : [0],
 
27
  "Line2" : [0],
 
28
  "Line3" : [0],
 
29
  "Line4" : [0],
 
30
  "Tri3" : [0,1],
 
31
  "Tri6" : [0,1,3],
 
32
  "Tri9" : [0,1,3,4], 
 
33
  "Tri10" : [0,1,3,4],
 
34
  "Rec4" : [0,1], 
 
35
  "Rec8" : [0,1,4], 
 
36
  "Rec9" : [0,1,4],
 
37
  "Rec12" : [0,1,4,5], 
 
38
  "Rec16" : [0,1,4,5], 
 
39
  "Tet4" : [0,1,2], 
 
40
  "Tet10" : [0,1,2,4,5,6], 
 
41
  "Tet16" : [0,1,2,4,5,6,7,8,9],
 
42
  "Hex8" : [0,1,2,3], 
 
43
  "Hex20" : [0,1,2,3,8,9,10,11],
 
44
  "Hex27" : [0,1,2,3,8,9,10,11,20],
 
45
  "Hex32" : [0,1,2,3,8,9,10,11,12,13,14,15]
 
46
}
 
47
 
 
48
FACENODES= {
 
49
  "Line2" : [0], 
 
50
  "Line3" : [0], 
 
51
  "Line4" : [0], 
 
52
  "Tri3" : [0,1],
 
53
  "Tri6" : [0,1,3], 
 
54
  "Tri9" : [0,1,3,4], 
 
55
  "Tri10" : [0,1,3,4], 
 
56
  "Rec4" : [0,1], 
 
57
  "Rec8" : [0,1,4], 
 
58
  "Rec9" : [0,1,4], 
 
59
  "Rec12" : [0,1,4,5], 
 
60
  "Rec16" : [0,1,4,5], 
 
61
  "Tet4" : [0,1,2,3] ,
 
62
  "Tet10" : [0,1,2,4,5,6],
 
63
  "Tet16" : [0,1,2,4,5,6,7,8,9], 
 
64
  "Hex8" : [0,1,2,3], 
 
65
  "Hex20" : [0,1,2,3,8,9,10,11], 
 
66
  "Hex27" : [0,1,2,3,8,9,10,11,20],
 
67
  "Hex32" : [0,1,2,3,8,9,10,11,12,13,14,15]
 
68
} 
 
69
 
 
70
SHIFTNODES= {
 
71
  "Point1" : ( (0,),(-1,)),
 
72
  "Line2" : ( (1,0),(-1,)),
 
73
  "Line3" : ( (1,0,2),(-1,)),
 
74
  "Line4" : ( (1,0,3,2),(-1,)),
 
75
  "Tri3" : ( (1,2,0),(0,2,1)),
 
76
  "Tri6" : ( (1,2,0,4,5,3),(0,2,1,5,4,3)),
 
77
  "Tri9" : ( (1,2,0,5,6,7,8,3,4),(0,2,1,8,7,6,5,4,3)),
 
78
  "Tri10" : ( (1,2,0,5,6,7,8,3,4,9),(0,2,1,8,7,6,5,4,3,9)),
 
79
  "Rec4" : ( (1,2,3,0),(0,3,2,1)),
 
80
  "Rec8" : ( (1,2,3,0,5,6,7,4),(0,3,2,1,7,6,5,4)),
 
81
  "Rec9" : ( (1,2,3,0,5,6,7,4,8),(0,3,2,1,7,6,5,4,8)),
 
82
  "Rec12" : ( (1,2,3,0,6,7,8,9,10,11,4,5),(0,3,2,1,11,10,9,8,7,6,5,4)) ,
 
83
  "Rec16" : ( (1,2,3,0,6,7,8,9,10,11,4,5,13,14,15,12),(0,3,2,1,11,10,9,8,7,6,5,4,12,15,14,13)),
 
84
  "Tet4" : ( (-1,),(-1,)),
 
85
  "Tet10" : ( (-1,),(-1,)),
 
86
  "Tet16" : ( (-1,),(-1,)),
 
87
  "Hex8" : ( (-1,),(-1,)),
 
88
  "Hex20" : ( (-1,),(-1,)),
 
89
  "Hex27" : ( (-1,),(-1,)),
 
90
  "Hex32" : ( (-1,),(-1,)),
 
91
  "Line2Face" : ( (0,1,2),(-1,)),
 
92
  "Line3Face" : ( (0,1,2),(-1,)),
 
93
  "Line4Face" : ( (0,1,2),(-1,)),
 
94
  "Tri3Face" : ( (1,0,2), (-1,)),
 
95
  "Tri6Face" : ( (1,0,2,3,5,4),(-1,)),
 
96
  "Tri9Face" : ( (1,0,2,4,3,7,8,6,5),(-1,)),
 
97
  "Tri10Face" : ( (1,0,2,4,3,7,8,6,5,9),(-1,)),
 
98
  "Rec4Face" : ( (1,0,3,2),(-1,)),
 
99
  "Rec8Face" : ( (1,0,3,2,4,7,6,5),(-1,)),
 
100
  "Rec9Face" : ( (1,0,3,2,4,7,6,5,8),(-1,)),
 
101
  "Rec12Face" : ( (1,0,3,2,5,4,11,10,9,8,7,6),(-1,)),
 
102
  "Rec16Face" : ( (1,0,3,2,5,4,11,10,9,8,7,6,13,12,15,14),(-1,)),
 
103
  "Tet4Face" : ( (1,2,0,3),(0,2,1,3)) ,
 
104
  "Tet10Face" : ( (1,2,0,3,5,6,4,8,9,7),(0,2,1,3,6,7,9,8)),
 
105
  "Tet16Face" : ( (1,2,0,3,6,7,8,9,4,5,11,12,10,14,15,13),(0,2,1,3,9,8,7,6,5,4,9,8,7,6,10,12,11,13,15,14)),
 
106
  "Hex8Face" : ( (1,2,3,0,5,6,7,4),(0,3,2,1,4,7,6,5)),
 
107
  "Hex20Face" : ( (1,2,3,0,5,6,7,4,9,10,11,8,13,14,15,12,17,18,19,16),(0,3,2,1,4,7,6,5,11,10,9,8,12,15,14,13,19,18,17,16)),
 
108
  "Hex27Face" : ( (1,2,3,0,5,6,7,4,9,10,11,8,13,14,15,12,17,18,19,16,20,22,23,24,22,25,26),(0,3,2,1,4,7,6,5,11,10,9,8,12,15,14,13,19,18,17,16,20,24,23,22,21,25,26)),
 
109
  "Hex32Face" : ( (1,2,3,0,5,6,7,4,10,11,12,13,14,15,8,9,17,18,19,16,21,22,23,20,26,27,28,29,30,31,34,25), (0,3,2,1,4,7,6,5,15,14,13,12,11,10,9,8,16,19,18,17,20,23,22,21,31,30,29,28,27,26,25,24))
 
110
}
 
111
def listToArrayStr(l):
 
112
   out="{"
 
113
   for s in l:
 
114
        if len(out)==1:
 
115
          out+=" %s"%s
 
116
        else:
 
117
          out+=", %s"%s
 
118
   return out+" }"
 
119
def LENLEN(l):
 
120
   if len(l) == 1:
 
121
      if l[0] == -1:
 
122
         return -1
 
123
      else:
 
124
         return 1
 
125
   else:
 
126
      return len(l)
 
127
 
 
128
outall="Dudley_ReferenceElementInfo Dudley_ReferenceElement_InfoList[]={\n"
 
129
 
 
130
for name in ["Point1", "Line2", "Tri3", "Tet4", "Line2Face", "Tri3Face", "Tet4Face" ]:
 
131
        isFace=False
 
132
        isMacro=False
 
133
        isContact=False
 
134
        n=1
 
135
        z="Point"
 
136
        for z in [ "Point", "Line", "Tri", "Tet"] :
 
137
            if name.startswith(z):
 
138
               s=name[len(z):]
 
139
               n=s
 
140
               if s.find("Macro")>=0:
 
141
                  isMacro=True
 
142
                  if s.find("Macro")>0: n=s[:s.find("Macro")]
 
143
                  s=s[s.find("Macro")+5:]
 
144
               if s.find("Face")>=0:
 
145
                  isFace=True
 
146
                  if s.find("Face")>0: n=s[:s.find("Face")]
 
147
                  s=s[s.find("Face")+4:]
 
148
               if s.find("_Contact")>=0:
 
149
                  isContact=True
 
150
                  if s.find("_Contact")>0: n=s[:s.find("_Contact")]
 
151
                  s=s[s.find("_Contact")+8:]
 
152
               n=int(n)
 
153
               break
 
154
        numLinearNodes=GEOBASE[z][1]
 
155
        if isFace:
 
156
           Quadrature=GEOBASE[z][2]
 
157
        else:
 
158
           Quadrature=z
 
159
        Parametrization="%s%s"%(z,n)
 
160
        if isContact:
 
161
            offsets=[0, n, 2*n]
 
162
            numSides=2
 
163
        else:
 
164
            numSides=1
 
165
            offsets=[0, n]
 
166
        numNodes=offsets[-1]
 
167
        if isMacro:
 
168
           numSubElements=GEOBASE[z][0]
 
169
           BasisFunctions="%s%s"%(z,numLinearNodes)
 
170
           subElementNodes=GEOBASE[z][3]
 
171
        else:
 
172
            numSubElements=1
 
173
            BasisFunctions=Parametrization
 
174
            subElementNodes=[ i for i in xrange(numNodes) ]
 
175
 
 
176
        linearTypeId="%s%s"%(z,numLinearNodes)
 
177
        linearNodes=[ i for i in xrange(numLinearNodes) ]
 
178
        if isFace: linearTypeId+="Face"
 
179
       
 
180
        if isContact: 
 
181
             linearTypeId+="_Contact"
 
182
             linearNodes+=[ n+i for i in xrange(numLinearNodes) ]
 
183
         
 
184
        if isFace:
 
185
           relevantGeoNodes=RELEVANTGEO["%s%s"%(z,n)]
 
186
        else:
 
187
           relevantGeoNodes=[ i for i in xrange(n) ]
 
188
 
 
189
        if isContact:
 
190
             faceNodes =[-1]
 
191
             shiftNodes = [-1]
 
192
             reverseNodes = [-1]
 
193
        else:
 
194
           if isFace:
 
195
                faceNodes=FACENODES["%s%s"%(z,n)]
 
196
                shiftNodes = SHIFTNODES["%s%sFace"%(z,n)][0]
 
197
                reverseNodes = SHIFTNODES["%s%sFace"%(z,n)][1]
 
198
           else:
 
199
                faceNodes=[  i for i in xrange(numNodes) ]
 
200
                shiftNodes = SHIFTNODES["%s%s"%(z,n)][0]
 
201
                reverseNodes = SHIFTNODES["%s%s"%(z,n)][1]
 
202
 
 
203
 
 
204
        out = '{ %s, "%s", %d, %d, %d, %s, %s,\n    %s, %sQuad, %sShape, %sShape,\n    %s,\n  %s, %s,\n  %s, %s,\n    %s,\n    %s },\n'% \
 
205
               (name, name, numNodes, numSubElements, numSides, listToArrayStr(offsets), linearTypeId, listToArrayStr(linearNodes),
 
206
                  Quadrature, Parametrization, BasisFunctions, listToArrayStr(subElementNodes),
 
207
                  len(relevantGeoNodes), listToArrayStr(relevantGeoNodes),
 
208
                  LENLEN(faceNodes), listToArrayStr(faceNodes), listToArrayStr(shiftNodes), listToArrayStr(reverseNodes) )
 
209
        outall+=out
 
210
out = '{ %s, "%s", %d, %d, %d, %s, %s,\n    %s, %s, %s, %s,\n    %s,\n  %s, %s,\n  %s, %s,\n    %s,\n    %s }\n'% \
 
211
        ("NoRef", "noElement", 0, 0, 0, listToArrayStr([0]), "NoRef", listToArrayStr([0]),
 
212
         "NoQuad", "NoShape", "NoShape", listToArrayStr([0]),
 
213
         -1, listToArrayStr([0]), -1, listToArrayStr([0]), listToArrayStr([0]), listToArrayStr([0]) )
 
214
outall+=out
 
215
print outall+"\n};"