2 #This scrip will create RVE including fibres and matrix, which can be used for all
3 #three type of boundary conditions including dispalcement, traction and periodic.
6 #=============================================================
7 #Geometry Input parametrs
8 #=============================================================
9 W_wrap=0.3; H_wrap=0.1514; hgap_wrap=0.09
10 W_weft=0.3; H_weft=0.0757; hgap_weft=1.2
11 L_RVE=3.0; W_RVE=0.8; H_RVE=0.3;
14 s_weft=hgap_weft+W_weft; s_wraf=hgap_wrap+W_wrap
15 T=(H_weft/2+H_wrap/2)/2+Vgap
17 L_RVE=2*s_weft; W_RVE=2*s_wraf; H_RVE=0.3;
19 #=============================================================
20 #Create vartices for splines to create splines
21 #=============================================================
23 coordx=-22.5*s_weft; coordy=T; coordz=0;
24 for i in range(0, 48):
26 str1='create vertex ' + str(coordx) +' '+str(coordy)+' '+str(coordz)+' color' ; coordx=coordx+s_weft;
29 str1='create vertex ' + str(coordx) +' '+str(-coordy)+' '+str(coordz)+' color' ; coordx=coordx+s_weft;
32 coordx=-22.5*s_weft; coordy=-T; coordz=s_wraf;
33 for i in range(0, 48):
35 str1='create vertex ' + str(coordx) +' '+str(coordy)+' '+str(coordz)+' color' ; coordx=coordx+s_weft;
38 str1='create vertex ' + str(coordx) +' '+str(-coordy)+' '+str(coordz)+' color' ; coordx=coordx+s_weft;
41 coordx=-0.5*s_weft; coordy=-T; coordz=-22*s_wraf;
42 for i in range(0, 48):
44 str1='create vertex ' + str(coordx) +' '+str(coordy)+' '+str(coordz)+' color' ; coordz=coordz+s_wraf;
47 str1='create vertex ' + str(coordx) +' '+str(-coordy)+' '+str(coordz)+' color' ; coordz=coordz+s_wraf;
50 coordx=0.5*s_weft; coordy=T; coordz=-22*s_wraf;
51 for i in range(0, 48):
53 str1='create vertex ' + str(coordx) +' '+str(coordy)+' '+str(coordz)+' color' ; coordz=coordz+s_wraf;
56 str1='create vertex ' + str(coordx) +' '+str(-coordy)+' '+str(coordz)+' color' ; coordz=coordz+s_wraf;
59 #=============================================================
60 #Joint vertices to create splines
61 #=============================================================
64 vertices=range(1,Tvertices+1); count1=0;
66 str1='create curve spline vertex '
67 for j in range(0, 48):
68 str1=str1+' '+ str(count1+vertices[j])
69 cubit.cmd(str1); count1=count1+48;
71 #=============================================================
72 #Creating fibers with elliptical cross sections
73 #=============================================================
75 cubit.cmd('create planar surface with plane normal to curve 1 distance 0 from vertex 1')
76 str1='create curve location vertex 1 direction curve 8 length '+ str(W_wrap/2); cubit.cmd(str1)
77 str1='create curve location vertex 1 direction curve 7 length '+str(H_wrap/2); cubit.cmd(str1)
78 cubit.cmd('create surface ellipse vertex 198 200 1')
79 cubit.cmd('sweep surface 2 along curve 1')
81 cubit.cmd('create planar surface with plane normal to curve 2 distance 0 from vertex 49')
82 str1='create curve location vertex 49 direction curve 17 length '+ str(W_wrap/2); cubit.cmd(str1)
83 str1='create curve location vertex 49 direction curve 16 length '+str(H_wrap/2); cubit.cmd(str1)
84 cubit.cmd('create surface ellipse vertex 208 210 49')
85 cubit.cmd('sweep surface 6 along curve 2')
87 cubit.cmd('create planar surface with plane normal to curve 3 distance 0 from vertex 97')
88 str1='create curve location vertex 97 direction curve 24 length '+ str(W_weft/2); cubit.cmd(str1)
89 str1='create curve location vertex 97 direction curve 23 length '+str(H_weft/2); cubit.cmd(str1)
90 cubit.cmd('create surface ellipse vertex 218 220 97')
91 cubit.cmd('sweep surface 10 along curve 3')
93 cubit.cmd('create planar surface with plane normal to curve 4 distance 0 from vertex 145')
94 str1='create curve location vertex 145 direction curve 33 length '+ str(W_weft/2); cubit.cmd(str1)
95 str1='create curve location vertex 145 direction curve 32 length '+str(H_weft/2); cubit.cmd(str1)
96 cubit.cmd('create surface ellipse vertex 228 230 145')
97 cubit.cmd('sweep surface 14 along curve 4')
99 #=============================================================
100 #Deletion of remaining free vertices, curves and bodies
101 #=============================================================
103 cubit.cmd('delete vertex all')
104 cubit.cmd('delete curve all')
105 cubit.cmd('delete body 7 3 5 1')
107 #=============================================================
108 #Create matrix and imprint and merge fibres and matrix to create common surfaces
109 #=============================================================
111 str1='brick x '+str(L_RVE)+' y '+str(H_RVE)+' z '+str(W_RVE); cubit.cmd(str1)
112 str1='move Volume '+str(9)+' x '+str(0)+' y '+str(0)+' z '+str(s_wraf/2); cubit.cmd(str1)
114 cubit.cmd('intersect volume all keep')
115 cubit.cmd('delete volume 6 2 8 4')
116 cubit.cmd('subtract volume 10 11 12 13 from volume 9 keep')
117 cubit.cmd('delete volume 9')
118 cubit.cmd('imprint volume 10 11 12 13 14')
119 cubit.cmd('merge volume 10 11 12 13 14')
121 #=============================================================
122 #Meshing negative x, y and z face of the RVE and coping the same to its positive faces
123 #=============================================================
125 cubit.cmd('group 1 add surface 42 27 24')
126 cubit.cmd('group "g2" add surface 44 25 28')
127 cubit.cmd('group "g3" add surface 39 30 33')
128 cubit.cmd('group "g4" add surface 31 40 34')
129 cubit.cmd('group "g5" add surface 43')
130 cubit.cmd('group "g6" add surface 41')
132 cubit.cmd('surface 42 27 24 size auto factor 6')
133 cubit.cmd('surface 42 27 24 scheme trimesh')
134 cubit.cmd('mesh surface 42 27 24')
135 cubit.cmd('surface 28 25 44 scheme mirror source surface 27 24 42 source vertex 263 target vertex 262 nosmoothing')
136 cubit.cmd('mesh surface 28 25 44')
138 cubit.cmd('surface 39 30 33 size auto factor 5')
139 cubit.cmd('surface 39 30 33 scheme trimesh')
140 cubit.cmd('mesh surface 39 30 33')
141 cubit.cmd('surface 31 34 40 scheme mirror source surface 30 33 39 source vertex 257 target vertex 262 nosmoothing')
142 cubit.cmd('mesh surface 31 34 40')
144 cubit.cmd('surface 43 scheme trimesh')
145 cubit.cmd('mesh surface 43')
146 cubit.cmd('surface 41 scheme mirror source surface 43 source vertex 259 target vertex 260 nosmoothing')
147 cubit.cmd('mesh surface 41')
149 #=============================================================
151 #=============================================================
153 cubit.cmd('volume all scheme Tetmesh')
154 cubit.cmd('mesh volume all')
156 #=============================================================
157 #Defining blocks for elastic, transversely-isotropic and potential flow problems
158 #=============================================================
160 vol=['14', '10,11,12,13','10', '11', '12', '13']
161 mat=['MAT_ELASTIC_1','MAT_ELASTIC_TRANSISO_1','PotentialFlow_1','PotentialFlow_2','PotentialFlow_3','PotentialFlow_4']
162 for i in range(0, 6):
163 cubit.cmd('set duplicate block elements on')
164 str1='block ' + str(i+1) +' volume '+vol[i]; cubit.cmd(str1)
165 str1='block ' + str(i+1) +' name "'+mat[i] + '"'; cubit.cmd(str1)
167 #=============================================================
168 #Material properties for matrix part
169 #=============================================================
171 cubit.cmd('block 1 attribute count 2')
172 Em=3.5e3; Enu=0.3; #gig to mega as we used dimension in mm
174 Elastic=[str(Em), str(Enu)]
175 for i in range(0, 2):
176 str1='block 1 attribute index ' + str(i+1) +' '+Elastic[i]; cubit.cmd(str1)
178 #=============================================================
179 #Material properties for fibres
180 #=============================================================
183 #cubit.cmd('block 2 attribute count 5')
184 #Ep=Em; Ez=Em; nup=Enu; nupz=Enu; Gzp=Em/(2*(1+Enu));
186 cubit.cmd('block 2 attribute count 5')
187 #Ep=40e3; Ez=230e3; nup=0.26; nupz=0.26; Gzp=24e3;
188 Ep=10*Em; Ez=20*Em; nup=0.26; nupz=0.26; Gzp=5*Em;
189 #Ep=2*Em; Ez=5*Em; nup=0.26; nupz=0.26; Gzp=2*Em;
191 TransIso=[str(Ep), str(Ez), str(nup), str(nupz), str(Gzp)]
192 for i in range(0, 5):
193 str1='block 2 attribute index ' + str(i+1) +' '+TransIso[i]; cubit.cmd(str1)
195 #=============================================================
196 #Material properties for interface between fibres and matrix
197 #=============================================================
200 cubit.cmd('set duplicate block elements on')
201 str1='block 7 surface 23 26 29 32'; cubit.cmd(str1)
202 str1='block 7 name "MAT_INTERF_1"'; cubit.cmd(str1)
203 cubit.cmd('block 7 attribute count 4')
204 str1='block 7 attribute index 1 '+str(alpha_interf); cubit.cmd(str1) #now we use 4 parameters for interface
205 str1='block 7 attribute index 2 '+str(0.0); cubit.cmd(str1)
206 str1='block 7 attribute index 3 '+str(0.0); cubit.cmd(str1)
207 str1='block 7 attribute index 4 '+str(0.0); cubit.cmd(str1)
210 #=============================================================
212 #=============================================================
214 Interface=['23', '26','29','32']
215 for i in range(0, 4):
216 str1='sideset ' + str(i+1) +' surface '+Interface[i]; cubit.cmd(str1)
217 str1='sideset ' + str(i+1) +' name "interface'+str(i+1); cubit.cmd(str1)
219 #=============================================================
220 #Defining pressures for potential flow problem
221 #=============================================================
223 Pres=['24', '25','27','28','30','31','33','34']; count=0; count1=len(Interface);
224 for i in range(0, 4):
225 str1='create pressure '+str(count+1)+' on surface '+str(Pres[count])+' magnitude 1'; cubit.cmd(str1)
226 str1='create pressure '+str(count+2)+' on surface '+str(Pres[count+1])+' magnitude -1'; cubit.cmd(str1)
227 str1='sideset '+str(count1+1)+' name "PressureIO_' + str(i+1) + '_1"'; cubit.cmd(str1)
228 str1='sideset '+str(count1+2)+' name "PressureIO_' + str(i+1) + '_2"'; cubit.cmd(str1)
229 count=count+2; count1=count1+2;
231 #=============================================================
232 #Defining surfaces for dispacement, traction and periodic boundary conditions
233 #=============================================================
235 cubit.cmd('sideset 101 surface 41 31 34 40 42 27 24') # all -ve boundary surfaces for periodic boundary conditions
236 cubit.cmd('sideset 102 surface 43 30 33 39 25 44 28') # all +ve boundary surfaces for periodic boundary conditions
237 cubit.cmd('sideset 103 surface 41 31 34 40 42 27 24 43 30 33 39 25 44 28') # all boundary surfaces
240 #To check for indiviual side for inserstion of prisms of any problem
241 #cubit.cmd('sideset 101 surface 42 27 24')
242 #cubit.cmd('sideset 102 surface 25 44 28')
243 #cubit.cmd('sideset 101 surface 41')
244 #cubit.cmd('sideset 102 surface 43')
245 #cubit.cmd('sideset 101 surface 31 34 40')
246 #cubit.cmd('sideset 102 surface 30 33 39')
247 #cubit.cmd('sideset 101 surface 40')
248 #cubit.cmd('sideset 102 surface 39')
250 #=============================================================
251 #Definign zero proessrues for potential flow problem (This should be of the same order as PotentialFlow blocks )
252 #=============================================================
254 zeroPressureNode=[140, 106, 279, 272]
255 for i in range(0, 4):
256 str1='nodeset ' + str(i+1) + ' node ' + str(zeroPressureNode[i]); cubit.cmd(str1)
257 str1='nodeset ' + str(i+1)+' name "ZeroPressure_' + str(i+1)+ '"'; cubit.cmd(str1)
259 #=============================================================
260 #Saving input RVE file
261 #=============================================================
263 cubit.cmd('save as "/Users/zahur/Documents/moFEM/mofem-cephas/mofem_v0.2/users_modules/homogenisation/meshes/Textile_RVE_TransIso.cub" overwrite')