3 import sys, salome, geompy, smesh, SMESH, math, copy, commands
4 CWD = commands.getoutput('pwd')
7 from MacObject import *
8 import Config, GenFunctions
10 def CompositeBoxF (Pt1 , Pt2 , Pt3 , Pt4 , **args ) :
11 [Pt1 , Pt2 , Pt3 , Pt4] = GenFunctions.SortPoints([Pt1 , Pt2 , Pt3 , Pt4])
12 if args.__contains__('groups') :
13 GroupNames = args['groups']
14 else : GroupNames = [None, None, None, None]
15 # Create a full NonOrtho box just to inherit, globally, the mesh parameters of bounding objects
16 dummy = MacObject('NonOrtho',[Pt1,Pt2,Pt3,Pt4],['auto'],publish=0)
17 # Save the existing number of segments on each direction
18 ExistingSeg0 = Config.ListObj[-1].DirectionalMeshParams
19 Convention = [2,3,0,1]
20 ExistingSegments = [ExistingSeg0[Convention[i]] for i in range(4)]
21 # Save Boundary lengths on each direction
22 BoundaryLengths = [IntLen(dummy.DirBoundaries(i)) for i in range(4) ]
23 # Calculate global mesh element size on each direction
24 GlobalDelta = [1.*BoundaryLengths[i]/ExistingSegments[i] for i in range(4) ]
25 print "GlobalDelta :",GlobalDelta
26 # Sort the connection list for the full Box
27 [(X0,Y0),(DX,DY)] = dummy.GeoPar
28 ObjIDLists = SortObjLists(Config.Connections[-1],X0 , Y0 , DX , DY )
29 [Xmin,Xmax,Ymin,Ymax] = dummy.Boundaries() # Used for groups determination
35 if not(args.__contains__('recursive')) :
38 Config.Criterion = GetCriterion(ObjIDLists)
39 for index, ObjList in enumerate(ObjIDLists) :
40 if not (ObjList[0] == -1 or Config.Count >= Config.Criterion):
41 if not(args.__contains__('recursive')) :
42 Config.DirIndex = index
43 if index > 1 : Config.RefPts = [Pt2, Pt3]
44 elif index == 0 : Config.RefPts = [Pt1, Pt2]
45 else : Config.RefPts = [Pt4, Pt3]
47 if len(ObjList)>1 : flag = 1
50 ToLook0 = [2,3,0,1][index]
51 ToLook1 = [3,2,1,0][index]
52 CommonSide = FindCommonSide(Config.ListObj[ObjID].DirBoundaries(ToLook1),dummy.DirBoundaries(ToLook0))
53 ToLook2 = [1,0,3,2][index]
54 RealSegments = Config.ListObj[ObjID].DirectionalMeshParams[ToLook2]*IntLen(CommonSide)/IntLen(Config.ListObj[ObjID].DirBoundaries(ToLook1))
55 LocalDelta = 1.*IntLen(CommonSide)/RealSegments
56 print "Direction:", ["West","East","South","North"][ToLook2]
57 print "IntLen(CommonSide):",IntLen(CommonSide)
58 print "RealSegments:",RealSegments
59 print "LocalDelta:",LocalDelta
60 if flag and Config.Count < Config.Criterion:
62 if abs(CommonSide[0] - Ymin)<1e-7 : SouthGR = GroupNames[0]
64 if abs(CommonSide[1] - Ymax)<1e-7 : NorthGR = GroupNames[1]
67 NDelta = Config.ListObj[ObjID].DirectionalMeshParams[ToLook2]* (LocalDelta-GlobalDelta[Convention[index]])
68 [Pt1,Pt2] = Config.RefPts
69 Coef = [1.,-1.][index]
70 Vref1 = [Coef*(Pt2[0]-Pt1[0]),Coef*(Pt2[1]-Pt1[1])]
71 Vref2 = NormalizeVector([Pt2[0]-Pt3[0],Pt2[1]-Pt3[1]])
72 Ptref = Config.ListObj[ObjID].PtCoor[[2,3][index]]
73 NewPt = ExtrapPoint (Ptref,Vref1,Vref2,NDelta)
74 CompositeBoxF (Pt1, Pt2, NewPt, Ptref, recursive=1, groups = [SouthGR,NorthGR]+GroupNames[2:4])
76 if abs(CommonSide[0] - Ymin)<1e-7 : SouthGR = GroupNames[0]
78 if abs(CommonSide[1] - Ymax)<1e-7 : NorthGR = GroupNames[1]
81 NDelta = Config.ListObj[ObjID].DirectionalMeshParams[ToLook2]* (LocalDelta-GlobalDelta[Convention[index]])
82 [Pt4,Pt3] = Config.RefPts
84 Vref1 = [Coef*(Pt4[0]-Pt3[0]),Coef*(Pt4[1]-Pt3[1])]
85 Vref2 = NormalizeVector([Pt1[0]-Pt4[0],Pt1[1]-Pt4[1]])
86 Ptref = Config.ListObj[ObjID].PtCoor[0]
87 NewPt = ExtrapPoint (Ptref,Vref1,Vref2,NDelta)
88 CompositeBoxF (NewPt, Ptref, Pt3, Pt4, recursive=1, groups = [SouthGR,NorthGR]+GroupNames[2:4])
90 if abs(CommonSide[0] - Xmin)<1e-7 : WestGR = GroupNames[2]
92 if abs(CommonSide[1] - Xmax)<1e-7 : EastGR = GroupNames[3]
95 NDelta = Config.ListObj[ObjID].DirectionalMeshParams[ToLook2]* (LocalDelta-GlobalDelta[Convention[index]])
96 [Pt2,Pt3] = Config.RefPts
97 Coef = [1.,-1.][index-2]
98 Vref1 = [Coef*(Pt3[0]-Pt2[0]),Coef*(Pt3[1]-Pt2[1])]
99 Vref2 = NormalizeVector([Pt3[0]-Pt4[0],Pt3[1]-Pt4[1]])
100 Ptref = Config.ListObj[ObjID].PtCoor[[3,0][index-2]]
101 NewPt = ExtrapPoint (Ptref,Vref1,Vref2,NDelta)
102 CompositeBoxF (Ptref, Pt2, Pt3, NewPt, recursive=1, groups = GroupNames[0:2] + [WestGR,EastGR])
104 if Config.Count >= Config.Criterion :
106 if flag == 0 and Config.Count < Config.Criterion:
107 print "Creating NonOrtho object with the points:", Pt1,Pt2,Pt3,Pt4
108 MacObject('NonOrtho',[Pt1,Pt2,Pt3,Pt4] ,['auto'], groups = GroupNames)
111 if Config.DirIndex > 1 : Config.RefPts = [Pt1, Pt4]
112 elif Config.DirIndex==0 : Config.RefPts = [Pt4, Pt3]
113 else : Config.RefPts = [Pt1, Pt2]
115 def FindCommonSide (Int1, Int2) :
116 if max(Int1[0],Int2[0])<min(Int1[1],Int2[1]): return [max(Int1[0],Int2[0]), min(Int1[1],Int2[1])]
118 print "Can not find interval intersection, returning [0,0]..."
121 def IntLen (Interval) :
122 return float(abs(Interval[1]-Interval[0]))
124 def RemoveLastObj() :
125 Config.ListObj = Config.ListObj[:-1]
126 Config.Connections = Config.Connections[:-1]
128 def NormalizeVector(V):
129 Magnitude = math.sqrt(GenFunctions.DotProd(V,V))
130 return [ V[i]/Magnitude for i in range(len(V))]
132 def GetCriterion (ObjListIDs):
133 return max(Config.Criterion, max(len(ObjListIDs[0]),len(ObjListIDs[1]))*max(len(ObjListIDs[2]),len(ObjListIDs[3])))
135 def SortObjLists (List,X0,Y0,DX,DY) :
137 This function sorts the list of neighbouring objects on each side, according to their intersection
138 with the object being created. From South to North and from East to West
141 # First find the directions where no neighbour exists
142 # Important : Here we assume that exactly two directions have no neighbours !!!
143 # Should we change this to allow a more general case ????
144 dummy = IndexMultiOcc(List,(-1,))
146 # dummy[0] is either 0, meaning there is no neighbour on X- (West)
147 # or 1, meaning there is no neighbour on X+ (East)
148 # Similarly dummy[1] can be either 2 or 3 (South and North respectively)
149 # In order to get back to the formalism of groups (SNWE)
150 # => we do the following to define Sense of no neighbours and then the Direction list
151 # is calculated as to include uniquely the directions where we DO have neighbours
153 # This adds a second direction where neighbours are not regarded, it is either 0 or 2
154 dummy.append(2*(dummy[0]+2<4))
155 print("Careful, you have neighbours on 3 or more sides of the box, we will not check if on two parallel sides the boxes are compatible !!!")
156 if len(dummy) == 2 or len(dummy) == 1 :
157 # Sense contains : Vertical then Horizontal
158 Sense = [dummy[1]%2,dummy[0]]
159 DirList = [[1,0][dummy[0]],[3,2][dummy[1]%2]]
160 for index,Direction in enumerate(DirList) :
161 ObjList = List[Direction]
163 ToLook0 = [2,2,0,0][Direction]
164 ToLook1 = [3,2,1,0][Direction]
165 for index1,ObjID in enumerate(ObjList) :
166 RankMin.append([-1.,1.][Sense[index]] * FindCommonSide(Config.ListObj[ObjID].DirBoundaries(ToLook1),[X0-DX/2.,X0+DX/2.,Y0-DY/2.,Y0+DY/2.][ToLook0:ToLook0+2])[Sense[index]])
167 Output[Direction] = SortList(ObjList,RankMin)
169 elif len(dummy) == 3 :
170 # We find the direction where we do have neighbours and then we sort the object list along it
172 Direction = [ i not in dummy for i in range(4) ].index(True)
173 ObjList = List[Direction]
175 ToLook0 = [2,2,0,0][Direction]
176 ToLook1 = [3,2,1,0][Direction]
177 for index1,ObjID in enumerate(ObjList) :
178 RankMin.append([-1.,1.][Sense] * FindCommonSide(Config.ListObj[ObjID].DirBoundaries(ToLook1),[X0-DX/2.,X0+DX/2.,Y0-DY/2.,Y0+DY/2.][ToLook0:ToLook0+2])[Sense])
179 Output[Direction] = SortList(ObjList,RankMin)
181 print ("Error : the composite box being created has no neighbours, how on earth do you want us to inherit its mesh parameters!!!")
186 def IndexMultiOcc (Array,Element) :
188 This functions returns the occurrences indices of Element in Array.
189 As opposed to Array.index(Element) method, this allows determining
190 multiple entries rather than just the first one!
193 try : Array.index(Element)
194 except ValueError : print "No more occurrences"
195 else : Output.append(Array.index(Element))
197 if not(Output == []) and len(Array) > 1 :
198 for index, ArrElem in enumerate(Array[Output[0]+1:]) :
199 if ArrElem == Element : Output.append(index+Output[0]+1)
203 def SortList (ValList, CritList):
205 SortedCritList = copy.copy(CritList)
206 SortedCritList.sort()
207 for i in range(0,len(ValList)):
208 index = CritList.index(SortedCritList[i])
209 Output.append(ValList[index])
212 def ExtrapPoint (Ptref,Vref1,Vref2,Delta):
214 This function allows determining the absolute coordinates of an extrapolation point
215 as shown in the following :
218 ExtrapPoint x---Vref2->--------o
221 Ptref x---------------------+
224 Delta = (delta_loc - delta_glob) * Nseg
227 X = Ptref[0] + Vref1[0] + Delta*Vref2[0]
228 Y = Ptref[1] + Vref1[1] + Delta*Vref2[1]