40
41
42
43 USE elbufdef_mod
46
47
48
49#include "implicit_f.inc"
50
51
52
53#include "mvsiz_p.inc"
54#include "param_c.inc"
55
56
57
58#include "drape_c.inc"
59
60
61
62 INTEGER JFT,JLT,NPT,NEL,IGTYP,ISUBSTACK,NLAY,NFT,IDRAPE,NUMEL_DRAPE
63 INTEGER MAT(*), PID(*), (*), IGEO(NPROPGI,*)
64 my_real geo(npropg,*),posly(mvsiz,*),thk(*)
65 TYPE(ELBUF_STRUCT_), TARGET :: ELBUF_STR
66 TYPE (STACK_PLY) :: STACK
67 TYPE (DRAPE_), DIMENSION(NUMELC_DRAPE + NUMELTG_DRAPE), TARGET :: DRAPE
68 INTEGER , DIMENSION(NUMEL_DRAPE) :: INDX
69
70
71
72 INTEGER I, J, N, IADR, IPTHK, IPMAT, IPPOS ,IPPID, IPID,
73 . IPANG,MAT_LY(MVSIZ),IT,,,NPTT,MAX_NPTT,IPT,JMLY,IINT,
74 . IPID_LY,IPT_ALL,MAT_LAY,NSLICE,IPOS,IE,IP
75 parameter(max_nptt = 100)
77 . thk_it(max_nptt*nlay,mvsiz),zshift,thk_nptt,
78 . thkl,pos_nptt,pos_0,thickt,thinning,thk_ly(mvsiz),
79 . thkly,ratio_thkly(mvsiz,npt)
80
81 TYPE (DRAPE_PLY_) , POINTER :: DRAPE_PLY
82
84 . a_gauss(9,9),w_gauss(9,9)
85
86 DATA a_gauss /
87 1 0. ,0. ,0. ,
88 1 0. ,0. ,0. ,
89 1 0. ,0. ,0. ,
90 2 -.577350269189626,0.577350269189626,0. ,
91 2 0. ,0. ,0. ,
92 2 0. ,0. ,0. ,
93 3 -.774596669241483,0. ,0.774596669241483,
94 3 0. ,0. ,0. ,
95 3 0. ,0. ,0. ,
96 4 -.861136311594053,-.339981043584856,0.339981043584856,
97 4 0.861136311594053,0. ,0. ,
98 4 0. ,0. ,0. ,
99 5 -.906179845938664,-.538469310105683,0. ,
100 5 0.538469310105683,0.906179845938664,0. ,
101 5 0. ,0. ,0. ,
102 6 -.932469514203152,-.661209386466265,-.238619186083197,
103 6 0.238619186083197,0.661209386466265,0.932469514203152,
104 6 0. ,0. ,0. ,
105 7 -.949107912342759,-.741531185599394,-.405845151377397,
106 7 0. ,0.405845151377397,0.741531185599394,
107 7 0.949107912342759,0. ,0. ,
108 8 -.960289856497536,-.796666477413627,-.525532409916329,
109 8 -.183434642495650,0.183434642495650,0.525532409916329,
110 8 0.796666477413627,0.960289856497536,0. ,
111 9 -.968160239507626,-.836031107326636,-.613371432700590,
112 9 -.324253423403809,0. ,0.324253423403809,
113 9 0.613371432700590,0.836031107326636,0.968160239507626/
114 DATA w_gauss /
115 1 2. ,0. ,0. ,
116 1 0. ,0. ,0. ,
117 1 0. ,0. ,0. ,
118 2 1. ,1. ,0. ,
119 2 0. ,0. ,0. ,
120 2 0. ,0. ,0. ,
121 3 0.555555555555556,0.888888888888889,0.555555555555556,
122 3 0. ,0. ,0. ,
123 3 0. ,0. ,0. ,
124 4 0.347854845137454,0.652145154862546,0.652145154862546,
125 4 0.347854845137454,0. ,0. ,
126 4 0. ,0. ,0. ,
127 5 0.236926885056189,0.478628670499366,0.568888888888889,
128 5 0.478628670499366,0.236926885056189,0. ,
129 5 0. ,0. ,0. ,
130 6 0.171324492379170,0.360761573048139,0.467913934572691,
131 6 0.467913934572691,0.360761573048139,0.171324492379170,
132 6 0. ,0. ,0. ,
133 7 0.129484966168870,0.279705391489277,0.381830050505119,
134 7 0.417959183673469,0.381830050505119,0.279705391489277,
135 7 0.129484966168870,0. ,0. ,
136 8 0.101228536290376,0.222381034453374,0.313706645877887,
137 8 0.362683783378362,0.362683783378362,0.313706645877887,
138 8 0.222381034453374,0.101228536290376,0. ,
139 9 0.081274388361574,0.180648160694857,0.260610696402935,
140 9 0.312347077040003,0.330239355001260,0.312347077040003,
141 9 0.260610696402935,0.180648160694857,0.081274388361574/
142
143 ipthk = 300
144 ippos = 400
145 ipmat = 100
146
147
148 SELECT CASE (igtyp)
149
150 CASE (1,9,10,11,16)
151 DO n=1,npt
152 iadr = (n-1)*jlt
153 pos_0 = geo(ippos+n,pid(1))
154 DO i = jft,jlt
155 j = iadr+i
156
157 posly(i,n) = pos_0
158 matly(j) = mat(1)
159 ENDDO
160 ENDDO
161
162 CASE (17)
163
164 ippid = 2
165 ipmat = ippid + npt
166 ipang = 1
167 ipthk = ipang + npt
168 ippos = ipthk + npt
169 ipos = igeo(99,pid(1))
170 zshift = geo(199,pid(1))
171 thickt = stack%GEO(1,isubstack)
172 IF(ipos == 2 ) zshift = zshift /
max(thickt,em20)
173 IF(idrape == 0 ) THEN
174 DO n=1,npt
175 iadr = (n-1)*jlt
176 DO i=jft,jlt
177 j = iadr+i
178 matly(j) = stack%IGEO(ipmat + n ,isubstack)
179 posly(i,n) = stack%GEO (ippos + n ,isubstack)
180 ENDDO
181 ENDDO
182 ELSE
183 DO n=1,npt
184 iadr
185 DO i=jft,jlt
186 j = iadr+i
187 ie = indx(nft + i)
188 IF(ie == 0) THEN
189 matly(j) = stack%IGEO(ipmat + n ,isubstack)
190 posly(i,n) = stack%GEO (ippos + n ,isubstack)
191 thickt = stack%GEO(1,isubstack)
192 thkly = stack%GEO (ipthk + n ,isubstack)*thickt
193 ratio_thkly(i,n) = thkly/thk(i)
194 IF (n == 1) THEN
195 posly(i,n) = zshift + half*ratio_thkly(i,n)
196 ELSE
197 posly(i,n) = posly(i,n-1)
198 . + half*(ratio_thkly(i,n)+ratio_thkly(i,n-1))
199 ENDIF
200 ELSE
201 ip= drape(ie)%INDX_PLY(n)
202 IF(ip > 0) THEN
203 drape_ply => drape(ie)%DRAPE_PLY(ip)
204 nslice = drape_ply%NSLICE
205 thinning = drape_ply%RDRAPE(1,1)
206 thickt = stack%GEO(1,isubstack)
207 thkly = stack%GEO(ipthk + n,isubstack)*thickt
208 thkly = thkly*thinning
209 thkly = thkly/thk(i)
210 ratio_thkly(i,n) = thkly
211 IF (n == 1) THEN
212 posly(i,n) = zshift + half*ratio_thkly(i,n)
213 ELSE
214 posly(i,n) = posly(i,n-1)
215 . + half*(ratio_thkly(i,n)+ratio_thkly(i,n-1))
216 ENDIF
217 ELSE
218 thickt = stack%GEO(1,isubstack)
219 thkly = stack%GEO(ipthk + n,isubstack)*thickt
220 thkly
221 ratio_thkly(i,n) = thkly
222 IF (n == 1) THEN
223 posly(i,n) = zshift + half*ratio_thkly(i,n)
224 ELSE
225 posly(i,n) = posly(i,n-1)
226 . + half*(ratio_thkly(i,n)+ratio_thkly(i,n-1))
227 ENDIF
228 ENDIF
229 ENDIF
230 ENDDO
231 ENDDO
232 ENDIF
233
234 CASE (51, 52)
235
236 ipt_all = 0
237
238 ipang = 1
239 ippid = 2
240 ipmat = ippid + nlay
241 ipthk = ipang + nlay
242 ippos = ipthk + nlay
243 ipos = igeo(99,pid(1))
244 zshift = geo(199,pid(1))
245 thickt = stack%GEO(1,isubstack)
246 IF(ipos == 2 )zshift = zshift /
max(thickt,em20)
247 IF(idrape == 0) THEN
248 DO ilay=1,nlay
249 nptt = elbuf_str%BUFLY(ilay)%NPTT
250 ipid_ly = stack%IGEO(ippid + ilay,isubstack)
251 ipid = stack%IGEO(ippid,isubstack)
252 iint = igeo(47,ipid)
253 mat_ly = elbuf_str%BUFLY(ilay)%IMAT
254
255 IF(iint == 1) THEN
256 DO i=jft,jlt
257 mat_ly(i) = stack%IGEO(ipmat + ilay,isubstack)
258 thk_ly(i) = stack%GEO(ipthk + ilay,isubstack)
259 ratio_thkly(i,ilay) = thk_ly
260 jmly = (ilay-1)*jlt + i
261 DO it=1,nptt
262 ipt = ipt_all + it
263 thk_it(ipt,i) = thk_ly(i)/nptt
264 IF (ipt == 1) THEN
265 posly(i,ipt) = zshift + half*thk_it(ipt,i)
266 ELSE
267 posly(i,ipt) = posly(i,ipt - 1)
268 . + half*(thk_it(ipt,i
269 ENDIF
270 ENDDO
271 matly(jmly) = mat_ly(i)
272 ENDDO
273 ELSEIF(iint == 2) THEN
274 DO i=jft,jlt
275 mat_ly(i) = stack%IGEO(ipmat + ilay,isubstack)
276 thk_ly(i) = stack%GEO(ipthk + ilay,isubstack)
277 mat_ly(i) = stack%IGEO(ipmat + ilay,isubstack)
278 ratio_thkly(i,ilay) = thk_ly(i)
279 jmly = (ilay-1)*jlt + i
280 DO it=1,nptt
281 ipt = ipt_all + it
282 thk_it(ipt,i) = half*thk_ly(i)*w_gauss(it,nptt)
283 IF (ipt == 1) THEN
284 posly(i,ipt) = zshift
285 ELSE
286 posly(i,ipt) = posly(i,ipt - 1)
287 . + half*(thk_it(ipt,i) + thk_it(ipt-1,i))
288 ENDIF
289 ENDDO
290 matly(jmly) = mat_ly(i)
291 ENDDO
292 ENDIF
293 ipt_all = ipt_all + nptt
294 ENDDO
295 ELSE
296 DO ilay=1,nlay
297 nptt = elbuf_str%BUFLY(ilay)%NPTT
298 ipid_ly = stack%IGEO(ippid + ilay,isubstack)
299 ipid = stack%IGEO(ippid,isubstack)
300 iint = igeo(47,ipid)
301 mat_ly = elbuf_str%BUFLY(ilay)%IMAT
302
303 IF(iint == 1) THEN
304 DO i=jft,jlt
305 ie = indx(nft + i)
306 mat_ly(i) = stack%IGEO(ipmat + ilay,isubstack)
307 IF(ie == 0) THEN
308 thk_ly(i) = stack%GEO(ipthk + ilay,isubstack)
309 ratio_thkly(i,ilay) = thk_ly(i)
310 jmly = (ilay-1)*jlt + i
311 DO it=1,nptt
312 ipt = ipt_all + it
313 thk_it(ipt,i) = thk_ly(i)/nptt
314 IF (ipt == 1) THEN
315 posly(i,ipt) = zshift + half*thk_it(ipt,i)
316 ELSE
317 posly(i,ipt) = posly(i,ipt - 1)
318 . + half*(thk_it(ipt,i) + thk_it(ipt-1,i))
319 ENDIF
320 matly(jmly) = mat_ly(i)
321 ENDDO
322 ELSE
323 ip = drape(ie)%INDX_PLY(ilay)
324 IF(ip > 0 ) THEN
325 drape_ply => drape(ie)%DRAPE_PLY(ip)
326 nslice = drape_ply%NSLICE
327 thickt = stack%GEO(1,isubstack)
328 thk_ly(i) = stack%GEO(ipthk + ilay,isubstack)*thickt
329 jmly = (ilay-1)*jlt + i
330 DO it = 1,nptt
331 ipt = ipt_all + it
332 j = (ipt-1)*jlt + i
333 thinning = drape_ply%RDRAPE(it,1)
334 thk_it(ipt,i) = thk_ly(i)*thinning/nptt
335 thk_it(ipt,i) = thk_it(ipt,i)/thk(i)
336
337 IF (ipt == 1 ) THEN
338 posly(i,ipt) = zshift + half*thk_it(ipt,i)
339 ELSE
340 posly(i,ipt) = posly(i,ipt - 1)
341 . + half*(thk_it(ipt,i) + thk_it(ipt-1,i
342 ENDIF
343 matly(jmly) = mat_ly(i
344 ENDDO
345 ELSE
346 thickt = stack%GEO(1,isubstack)
347 thk_ly(i) = stack%GEO(ipthk + ilay,isubstack)*thickt
348 jmly = (ilay-1)*jlt + i
349 DO it = 1,nptt
350 ipt = ipt_all + it
351 j = (ipt-1)*jlt + i
352 thk_it(ipt,i) = thk_ly(i)/nptt
353 thk_it(ipt,i) = thk_it(ipt,i)/thk(i)
354
355 IF (ipt == 1 ) THEN
356 posly(i,ipt) = zshift + half*thk_it(ipt,i)
357 ELSE
358 posly(i,ipt) = posly(i,ipt - 1)
359 . + half*(thk_it(ipt,i) + thk_it(ipt-1,i))
360 ENDIF
361 matly(jmly) = mat_ly(i)
362 ENDDO
363 ENDIF
364 ENDIF
365 ENDDO
366 ELSEIF(iint == 2) THEN
367 DO i=jft,jlt
368 ie = indx(nft + i)
369 mat_ly(i) = stack%IGEO(ipmat + ilay,isubstack)
370 IF(ie == 0) THEN
371 thk_ly(i) = stack%GEO(ipthk + ilay,isubstack)
372 mat_ly(i) = stack%IGEO(ipmat + ilay,isubstack)
373 ratio_thkly(i,ilay) = thk_ly(i)
374 jmly = (ilay-1)*jlt + i
375 DO it=1,nptt
376 ipt = ipt_all + it
377 thk_it(ipt,i) = half*thk_ly(i)*w_gauss(it,nptt)
378 IF (ipt == 1) THEN
379 posly(i,ipt) = zshift + half*thk_it(ipt,i)
380 ELSE
381 posly(i,ipt) = posly(i,ipt - 1)
382 . + half*(thk_it(ipt,i) + thk_it(ipt-1,i))
383 ENDIF
384 matly(jmly) = mat_ly(i)
385 ENDDO
386 ELSE
387
388 ip = drape(ie)%INDX_PLY(ilay)
389 IF(ip > 0) THEN
390 drape_ply => drape(ie)%DRAPE_PLY(ip)
391 nslice = drape_ply%NSLICE
392 thickt = stack%GEO(1,isubstack)
393 thk_ly(i) = stack%GEO(ipthk + ilay,isubstack)*thickt
394 jmly = (ilay-1)*jlt + i
395 DO it = 1,nptt
396 ipt = ipt_all + it
397 j = (ipt-1)*jlt + i
398 thinning = drape_ply%RDRAPE(it,1)
399 thk_it(ipt,i) = half*thk_ly(i)*w_gauss(it,nptt)*thinning
400 thk_it(ipt,i) = thk_it(ipt,i)/thk(i)
401 IF (ipt == 1 ) THEN
402 posly(i,ipt) = zshift + half*thk_it(ipt,i)
403 ELSE
404 posly(i,ipt) = posly(i,ipt - 1)
405 . + half*(thk_it(ipt,i) + thk_it(ipt-1,i))
406 ENDIF
407 matly(jmly) = mat_ly(i)
408 ENDDO
409 ELSE
410 thickt = stack%GEO(1,isubstack)
411 thk_ly(i) = stack%GEO(ipthk + ilay,isubstack)*thickt ! initial thkly
412 jmly = (ilay-1)*jlt + i
413 DO it = 1,nptt
414 ipt = ipt_all + it
415 j = (ipt-1)*jlt + i
416 thk_it(ipt,i) = half*thk_ly(i)*w_gauss(it,nptt)
417 thk_it(ipt,i) = thk_it(ipt,i)/thk(i)
418 IF (ipt == 1 ) THEN
419 posly(i,ipt) = zshift + half*thk_it(ipt,i)
420 ELSE
421 posly(i,ipt) = posly(i,ipt - 1)
422 . + half*(thk_it(ipt,i) + thk_it(ipt-1,i))
423 ENDIF
424 matly(jmly) = mat_ly(i)
425 ENDDO
426 ENDIF
427 ENDIF
428 ENDDO
429 ENDIF
430 ipt_all = ipt_all + nptt
431 ENDDO
432 ENDIF
433
434 CASE DEFAULT
435
436 DO n=1,npt
437 iadr = (n-1)*jlt
438 pos_0 = geo(ippos+n,pid(1))
439 thk_nptt = geo(ipthk+n,pid(1))
440 DO i = jft,jlt
441 j = iadr+i
442
443 posly(i,n) = pos_0
444 matly(j) = mat(1)
445 ENDDO
446 ENDDO
447
448 END SELECT
449
450 RETURN