v0.14.0
Textile_RVE_TransIso.jou
Go to the documentation of this file.
1 #!python
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.
4 
5 cubit.cmd('new')
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;
12 Vgap=0.012
13 
14 s_weft=hgap_weft+W_weft; s_wraf=hgap_wrap+W_wrap
15 T=(H_weft/2+H_wrap/2)/2+Vgap
16 
17 L_RVE=2*s_weft; W_RVE=2*s_wraf; H_RVE=0.3;
18 
19 #=============================================================
20 #Create vartices for splines to create splines
21 #=============================================================
22 
23 coordx=-22.5*s_weft; coordy=T; coordz=0;
24 for i in range(0, 48):
25  if i % 2 == 0:
26  str1='create vertex ' + str(coordx) +' '+str(coordy)+' '+str(coordz)+' color' ; coordx=coordx+s_weft;
27  cubit.cmd(str1)
28  if i % 2 == 1:
29  str1='create vertex ' + str(coordx) +' '+str(-coordy)+' '+str(coordz)+' color' ; coordx=coordx+s_weft;
30  cubit.cmd(str1)
31 
32 coordx=-22.5*s_weft; coordy=-T; coordz=s_wraf;
33 for i in range(0, 48):
34  if i % 2 == 0:
35  str1='create vertex ' + str(coordx) +' '+str(coordy)+' '+str(coordz)+' color' ; coordx=coordx+s_weft;
36  cubit.cmd(str1)
37  if i % 2 == 1:
38  str1='create vertex ' + str(coordx) +' '+str(-coordy)+' '+str(coordz)+' color' ; coordx=coordx+s_weft;
39  cubit.cmd(str1)
40 
41 coordx=-0.5*s_weft; coordy=-T; coordz=-22*s_wraf;
42 for i in range(0, 48):
43  if i % 2 == 0:
44  str1='create vertex ' + str(coordx) +' '+str(coordy)+' '+str(coordz)+' color' ; coordz=coordz+s_wraf;
45  cubit.cmd(str1)
46  if i % 2 == 1:
47  str1='create vertex ' + str(coordx) +' '+str(-coordy)+' '+str(coordz)+' color' ; coordz=coordz+s_wraf;
48  cubit.cmd(str1)
49 
50 coordx=0.5*s_weft; coordy=T; coordz=-22*s_wraf;
51 for i in range(0, 48):
52  if i % 2 == 0:
53  str1='create vertex ' + str(coordx) +' '+str(coordy)+' '+str(coordz)+' color' ; coordz=coordz+s_wraf;
54  cubit.cmd(str1)
55  if i % 2 == 1:
56  str1='create vertex ' + str(coordx) +' '+str(-coordy)+' '+str(coordz)+' color' ; coordz=coordz+s_wraf;
57  cubit.cmd(str1)
58 
59 #=============================================================
60 #Joint vertices to create splines
61 #=============================================================
62 
63 Tvertices=48;
64 vertices=range(1,Tvertices+1); count1=0;
65 for i in range(0, 4):
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;
70 
71 #=============================================================
72 #Creating fibers with elliptical cross sections
73 #=============================================================
74 
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')
80 
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')
86 
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')
92 
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')
98 
99 #=============================================================
100 #Deletion of remaining free vertices, curves and bodies
101 #=============================================================
102 
103 cubit.cmd('delete vertex all')
104 cubit.cmd('delete curve all')
105 cubit.cmd('delete body 7 3 5 1')
106 
107 #=============================================================
108 #Create matrix and imprint and merge fibres and matrix to create common surfaces
109 #=============================================================
110 
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)
113 
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')
120 
121 #=============================================================
122 #Meshing negative x, y and z face of the RVE and coping the same to its positive faces
123 #=============================================================
124 
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')
131 
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')
137 
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')
143 
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')
148 
149 #=============================================================
150 #Mesh volume
151 #=============================================================
152 
153 cubit.cmd('volume all scheme Tetmesh')
154 cubit.cmd('mesh volume all')
155 
156 #=============================================================
157 #Defining blocks for elastic, transversely-isotropic and potential flow problems
158 #=============================================================
159 
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)
166 
167 #=============================================================
168 #Material properties for matrix part
169 #=============================================================
170 
171 cubit.cmd('block 1 attribute count 2')
172 Em=3.5e3; Enu=0.3; #gig to mega as we used dimension in mm
173 
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)
177 
178 #=============================================================
179 #Material properties for fibres
180 #=============================================================
181 
182 #to use as isotropic
183 #cubit.cmd('block 2 attribute count 5')
184 #Ep=Em; Ez=Em; nup=Enu; nupz=Enu; Gzp=Em/(2*(1+Enu));
185 
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;
190 
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)
194 
195 #=============================================================
196 #Material properties for interface between fibres and matrix
197 #=============================================================
198 
199 alpha_interf=500
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)
208 
209 
210 #=============================================================
211 #Defining interfaces
212 #=============================================================
213 
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)
218 
219 #=============================================================
220 #Defining pressures for potential flow problem
221 #=============================================================
222 
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;
230 
231 #=============================================================
232 #Defining surfaces for dispacement, traction and periodic boundary conditions
233 #=============================================================
234 
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
238 
239 
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')
249 
250 #=============================================================
251 #Definign zero proessrues for potential flow problem (This should be of the same order as PotentialFlow blocks )
252 #=============================================================
253 
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)
258 
259 #=============================================================
260 #Saving input RVE file
261 #=============================================================
262 
263 cubit.cmd('save as "/Users/zahur/Documents/moFEM/mofem-cephas/mofem_v0.2/users_modules/homogenisation/meshes/Textile_RVE_TransIso.cub" overwrite')