OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
stackgroup_drape.F
Go to the documentation of this file.
1Copyright> OpenRadioss
2Copyright> Copyright (C) 1986-2025 Altair Engineering Inc.
3Copyright>
4Copyright> This program is free software: you can redistribute it and/or modify
5Copyright> it under the terms of the GNU Affero General Public License as published by
6Copyright> the Free Software Foundation, either version 3 of the License, or
7Copyright> (at your option) any later version.
8Copyright>
9Copyright> This program is distributed in the hope that it will be useful,
10Copyright> but WITHOUT ANY WARRANTY; without even the implied warranty of
11Copyright> MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12Copyright> GNU Affero General Public License for more details.
13Copyright>
14Copyright> You should have received a copy of the GNU Affero General Public License
15Copyright> along with this program. If not, see <https://www.gnu.org/licenses/>.
16Copyright>
17Copyright>
18Copyright> Commercial Alternative: Altair Radioss Software
19Copyright>
20Copyright> As an alternative to this open-source version, Altair also offers Altair Radioss
21Copyright> software under a commercial license. Contact Altair to discuss further if the
22Copyright> commercial version may interest you: https://www.altair.com/radioss/.
23!||====================================================================
24!|| stackgroup_drape ../starter/source/stack/stackgroup_drape.F
25!||--- called by ------------------------------------------------------
26!|| lectur ../starter/source/starter/lectur.F
27!||--- calls -----------------------------------------------------
28!||--- uses -----------------------------------------------------
29!|| drape_mod ../starter/share/modules1/drape_mod.F
30!|| message_mod ../starter/share/message_module/message_mod.F
31!|| stack_mod ../starter/share/modules1/stack_mod.F
32!|| submodel_mod ../starter/share/modules1/submodel_mod.F
33!||====================================================================
34 SUBROUTINE stackgroup_drape(DRAPE ,DRAPEG , IWORK_T, IWORKSH ,
35 . IGRSH3N ,IGRSH4N ,IXC ,IXTG ,
36 . IGEO ,GEO ,THK ,STACK ,
37 . IGEO_STACK ,GEO_STACK , STACK_INFO ,
38 . NUMGEO_STACK,NPROP_STACK , PLY_INFO)
39C-----------------------------------------------
40C M o d u l e s
41C-----------------------------------------------
42 USE submodel_mod
43 USE stack_mod
44 USE message_mod
45 USE groupdef_mod
46 USE drape_mod
47 USE names_and_titles_mod , ONLY : nchartitle
48C-----------------------------------------------
49C I m p l i c i t T y p e s
50C-----------------------------------------------
51#include "implicit_f.inc"
52C-----------------------------------------------
53C C o m m o n B l o c k s
54C-----------------------------------------------
55#include "com01_c.inc"
56#include "scr17_c.inc"
57#include "scr03_c.inc"
58#include "com04_c.inc"
59#include "units_c.inc"
60#include "warn_c.inc"
61#include "param_c.inc"
62#include "remesh_c.inc"
63#include "sphcom.inc"
64#include "drape_c.inc"
65C-----------------------------------------------
66C D u m m y A r g u m e n t s
67C-----------------------------------------------
68 INTEGER IXC(NIXC,NUMELC),
69 . IXTG(NIXTG,NUMELTG),IGEO(NPROPGI,NUMGEO),IWORKSH(3,NUMELC+NUMELTG),
70 . IGEO_STACK(NPROPGI,NUMSTACK + NUMPLY),NUMGEO_STACK(NUMGEO+NUMSTACK),
71 . NPROP_STACK
72 INTEGER , INTENT(INOUT) :: PLY_INFO(2,NUMPLY)
73 my_real
74 . GEO(NPROPG,NUMGEO),THK(NUMELC+NUMELTG),GEO_STACK(NPROPG,NUMSTACK + NUMPLY)
75C-----------------------------------------------
76 TYPE (GROUP_) , DIMENSION(NGRSH3N) :: IGRSH3N
77 TYPE (GROUP_) , DIMENSION(NGRSHEL) :: IGRSH4N
78 TYPE (DRAPE_) , DIMENSION(NUMELC + NUMELTG) , TARGET :: DRAPE
79 TYPE (DRAPEG_) , TARGET :: DRAPEG
80 TYPE(drape_work_) , DIMENSION(NUMELC + NUMELTG) , TARGET :: IWORK_T
81C-----------------------------------------------
82C L o c a l V a r i a b l e s
83C-----------------------------------------------
84 INTEGER I,J,II,IGTYP,ID,JD,IDPLY,NEL,
85 . IAD,ITY,IDSHEL,PID,IS,IDS,NSH,MODE,JJ,NGEO_STACK,
86 . IGRTYP,N1,IPMAT,IPANG,IPTHK,IIGEO,NSS,IPPOS,NPT,IIS,NP,
87 . JJPID,JSTACK,JPID,ITG,IPMAT_IPLY,ISH3N,J4N,J3N,IPOS,
88 . mat_ly,nlay,nptt,ipidl,it,ilay,ipthk_nptt,ippos_nptt,
89 . iint,ipid_ly,ipdir ,ns_stack0 ,npt_stack0,is0,js,pids,ip,
90 . ii1,ii2,jj1,jj2,nslice,ie_drp,npt_lay,ipnpt_lay,
91 . ibit, nkey, ikey,irest,n_ply,nbit,nply,ns_total,ns_first,
92 . ns_sub, nstack,i1,i2,ies
93 INTEGER , DIMENSION(:), ALLOCATABLE :: WORK,INDX_SH,PID_SH,
94 . NFIRST,NLAST,ISUBSTACK,IPTPLY,
95 . NFIRST1,NLAST1,NFIRST2,NLAST2,
96 . INDX1,INDX_TOTAL
97 INTEGER :: NBFI,IPPID,NGL,IPID_1,NUMS,IPWEIGHT,IPTHKLY,NSHQ4,NSHT3
98 my_real
99 . THICKT,ZSHIFT,TMIN,TMAX,DT,THK_LY,POS_LY,THK_IT(100),
100 . POS_IT(100),POS_NPTT,THK_NPTT,POS_0,THINNING,POS
101
102 INTEGER, DIMENSION(:,:) ,ALLOCATABLE :: ITRI,ACTIV_PLY
103 INTEGER, DIMENSION(:) ,ALLOCATABLE ::INDX,IDSTACK,INDX_SUB
104 TYPE (STACK_PLY) :: STACK, IWORKS
105 TYPE (STACK_INFO_ ) , DIMENSION (1:NPROP_STACK) :: STACK_INFO
106 TYPE (DRAPE_PLY_), POINTER :: DRAPE_PLY
107 CHARACTER(LEN=NCHARTITLE) :: TITR,TITR1
108C-----------------------------------------------
109 my_real
110 . a_gauss(9,9),w_gauss(9,9)
111C-----------------------------------------------
112 DATA a_gauss /
113 1 0. ,0. ,0. ,
114 1 0. ,0. ,0. ,
115 1 0. ,0. ,0. ,
116 2 -.577350269189626,0.577350269189626,0. ,
117 2 0. ,0. ,0. ,
118 2 0. ,0. ,0. ,
119 3 -.774596669241483,0. ,0.774596669241483,
120 3 0. ,0. ,0. ,
121 3 0. ,0. ,0. ,
122 4 -.861136311594053,-.339981043584856,0.339981043584856,
123 4 0.861136311594053,0. ,0. ,
124 4 0. ,0. ,0. ,
125 5 -.906179845938664,-.538469310105683,0. ,
126 5 0.538469310105683,0.906179845938664,0. ,
127 5 0. ,0. ,0. ,
128 6 -.932469514203152,-.661209386466265,-.238619186083197,
129 6 0.238619186083197,0.661209386466265,0.932469514203152,
130 6 0. ,0. ,0. ,
131 7 -.949107912342759,-.741531185599394,-.405845151377397,
132 7 0. ,0.405845151377397,0.741531185599394,
133 7 0.949107912342759,0. ,0. ,
134 8 -.960289856497536,-.796666477413627,-.525532409916329,
135 8 -.183434642495650,0.183434642495650,0.525532409916329,
136 8 0.796666477413627,0.960289856497536,0. ,
137 9 -.968160239507626,-.836031107326636,-.613371432700590,
138 9 -.324253423403809,0. ,0.324253423403809,
139 9 0.613371432700590,0.836031107326636,0.968160239507626/
140 DATA w_gauss /
141 1 2. ,0. ,0. ,
142 1 0. ,0. ,0. ,
143 1 0. ,0. ,0. ,
144 2 1. ,1. ,0. ,
145 2 0. ,0. ,0. ,
146 2 0. ,0. ,0. ,
147 3 0.555555555555556,0.888888888888889,0.555555555555556,
148 3 0. ,0. ,0. ,
149 3 0. ,0. ,0. ,
150 4 0.347854845137454,0.652145154862546,0.652145154862546,
151 4 0.347854845137454,0. ,0. ,
152 4 0. ,0. ,0. ,
153 5 0.236926885056189,0.478628670499366,0.568888888888889,
154 5 0.478628670499366,0.236926885056189,0. ,
155 5 0. ,0. ,0. ,
156 6 0.171324492379170,0.360761573048139,0.467913934572691,
157 6 0.467913934572691,0.360761573048139,0.171324492379170,
158 6 0. ,0. ,0. ,
159 7 0.129484966168870,0.279705391489277,0.381830050505119,
160 7 0.417959183673469,0.381830050505119,0.279705391489277,
161 7 0.129484966168870,0. ,0. ,
162 8 0.101228536290376,0.222381034453374,0.313706645877887,
163 8 0.362683783378362,0.362683783378362,0.313706645877887,
164 8 0.222381034453374,0.101228536290376,0. ,
165 9 0.081274388361574,0.180648160694857,0.260610696402935,
166 9 0.312347077040003,0.330239355001260,0.312347077040003,
167 9 0.260610696402935,0.180648160694857,0.081274388361574/
168C----------------------------f-------------------
169c=======================================================================
170c define temporary work structure
171c=======================================================================
172C=======================================================================
173C For Shell
174C-----------------------------------------------
175 ns_stack = 0
176 npt_stack = 0
177C
178 ALLOCATE (indx_sh(numelc+numeltg),pid_sh(numelc+numeltg),
179 . isubstack(numgeo+numstack),
180 . iptply(numgeo+numply), work(70000) )
181
182 work = 0
183 indx_sh = 0
184 pid_sh = 0
185 isubstack = 0
186 iptply = 0
187 indx_sh = 0
188 pid_sh = 0
189 !! Type51 and TYPE17
190 IF(ipart_stack > 0) THEN
191 !
192 nply = 0
193 DO i = 1, numgeo
194 igtyp=igeo(11,i)
195 nstack = igeo(42,i) ! number of stack where ply is attached
196 IF (igtyp == 19 .AND. nstack > 0) nply = nply+1
197 ENDDO
198 !!
199 nsh = 0
200 DO i=1,numelc
201 pid = ixc(6,i)
202 igtyp = igeo(11,pid)
203 IF(igtyp == 17 .OR. igtyp == 51)THEN
204 nsh = nsh + 1
205 indx_sh(nsh) = i
206 pid_sh(nsh) = pid
207 ENDIF
208 ENDDO
209C
210 DO i=1,numeltg
211 pid = ixtg(5,i)
212 igtyp = igeo(11,pid)
213 IF(igtyp == 17 .OR. igtyp == 51)THEN
214 nsh = nsh + 1
215 indx_sh(nsh) = i + numelc
216 pid_sh(nsh) = pid
217 ENDIF
218 ENDDO
219C #####################################################"
220 nbit = bit_size(nply)
221 irest = mod(nply,nbit)
222 nkey = nply / nbit
223 IF(irest > 0) nkey = nkey + 1
224 ALLOCATE( activ_ply(numelc+numeltg,nkey))
225 IF(numelc + numeltg > 0)activ_ply = 0
226 !!
227 DO i=1,nsh
228 ii = indx_sh(i)
229 pid = pid_sh(i)
230 igtyp = igeo(11,pid)
231 npt = iworksh(1,ii)
232 ie_drp = drapeg%INDX(ii)
233 IF(igtyp == 17) THEN
234 IF(npt > 0) THEN
235 ALLOCATE(iwork_t(ii)%NPT_PLY(npt))
236 iwork_t(ii)%NPT_PLY = 0
237 DO j=1,npt
238 n_ply = iwork_t(ii)%PLYNUM(j)
239 ikey = n_ply / nbit ! 32 or 64 bits
240 IF(mod(n_ply,nbit) > 0 ) ikey = ikey + 1
241 ikey = min(ikey, nkey)
242 ibit = n_ply - (ikey - 1)*nbit ! Bit corresponding to PLy
243 idply = iwork_t(ii)%PLYID(j)
244 iwork_t(ii)%NPT_PLY(j) = 1
245 activ_ply(ii,ikey) = ibset(activ_ply(ii,ikey),ibit)
246 ENDDO
247 ENDIF
248 ELSEIF(igtyp == 51) THEN
249 ALLOCATE(iwork_t(ii)%NPT_PLY(npt))
250 iwork_t(ii)%NPT_PLY = 0
251 IF(ie_drp > 0 .AND. npt > 0) THEN
252 DO j=1,npt
253 ip = drape(ie_drp)%INDX_PLY(j)
254 n_ply = iwork_t(ii)%PLYNUM(j)
255 ikey = n_ply / nbit ! 32 or 64 bits
256 IF(mod(n_ply,nbit) > 0 ) ikey = ikey + 1
257 ikey = min(ikey, nkey)
258 ibit = n_ply - (ikey - 1)*nbit ! Bit corresponding to PLy
259 activ_ply(ii,ikey) = ibset(activ_ply(ii,ikey),ibit)
260 IF(ip > 0 ) THEN
261 drape_ply => drape(ie_drp)%DRAPE_PLY(ip)
262 nslice = drape_ply%NSLICE
263 idply = iwork_t(ii)%PLYID(j)
264 iwork_t(ii)%NPT_PLY(j) = nslice
265 igeo(44,idply) = max(igeo(4,idply),nslice)
266 ELSE
267 idply = iwork_t(ii)%PLYID(j)
268 npt_lay = igeo(4,idply)
269 iwork_t(ii)%NPT_PLY(j) = npt_lay
270 ENDIF
271 ENDDO
272 ELSE
273 IF(npt > 0) THEN
274 DO j=1,npt
275 idply = iwork_t(ii)%PLYID(j)
276 npt_lay = igeo(4,idply)
277 iwork_t(ii)%NPT_PLY(j) = npt_lay
278 n_ply = iwork_t(ii)%PLYNUM(j)
279 ikey = n_ply / nbit ! 32 or 64 bits
280 IF(mod(n_ply,nbit) > 0 ) ikey = ikey + 1
281 ikey = min(ikey, nkey)
282 ibit = n_ply - (ikey - 1)*nbit ! Bit corresponding to PLy
283 activ_ply(ii,ikey) = ibset(activ_ply(ii,ikey),ibit)
284 ENDDO
285 ENDIF
286 ENDIF
287 ENDIF ! IGTYP
288 ENDDO
289!!###########################################################"
290 ALLOCATE(indx(2*nsh),itri(2 + nkey,nsh),indx_total(nsh))
291 ALLOCATE (nfirst(nsh) ,nlast(nsh),
292 . nfirst1(nsh) ,nlast1(nsh))
293 IF(nsh > 0) THEN
294 indx = 0
295 itri = 0
296 indx_total = 0
297 nlast = 0
298 nfirst = 0
299 nfirst1 = 0
300 nlast1 = 0
301 ENDIF
302 DO i = 1,nsh
303 indx(i) = i
304 ii = indx_sh(i)
305 itri(1,i) = pid_sh(i)
306 itri(2,i) = iworksh(1,ii)
307 DO j=1,nkey
308 itri(2+j,i) = activ_ply(ii,j)
309 ENDDO
310 ENDDO
311C
312 mode = 0
313C
314 nkey = nkey + 2
315 CALL my_orders(mode, work, itri, indx, nsh , nkey)
316 ns_first = 1
317 nfirst1(1) = 1
318 nlast1(1) = 1
319 DO i=2,nsh
320 DO ikey = 1, nkey
321 ii = itri(ikey,indx(i))
322 jj = itri(ikey,indx(i-1))
323 IF(ii /= jj) THEN
324 ns_first = ns_first + 1
325 nfirst1(ns_first) = i
326 nlast1(ns_first) = nfirst1(ns_first)
327 EXIT
328 ELSEIF(ikey == nkey) THEN
329 nlast1(ns_first) = nlast1(ns_first) + 1
330 ENDIF
331 ENDDO
332 ENDDO
333 DEALLOCATE (itri)
334 !! sorting depending on number of slice = int point by ply for Type51
335 ns_total = 0
336 ns_sub = 0
337 !!
338 DO is =1,ns_first
339 id = nfirst1(is)
340 i = indx(id)
341 ii = indx_sh(i)
342 pid = pid_sh(i)
343 igtyp = igeo(11,pid)
344 nkey = iworksh(1,ii)
345 IF(igtyp == 51) THEN
346 nsh = nlast1(is) - nfirst1(is) + 1
347 ALLOCATE(indx1(2*nsh),itri(nkey,nsh),indx_sub(nsh))
348 indx1 = 0
349 itri = 0
350 i1 = 0
351 DO i= nfirst1(is), nlast1(is)
352 id = indx(i)
353 ii = indx_sh(id)
354 i1 = i1 + 1
355 indx1(i1) = i1
356 indx_sub(i1) = id
357 DO j=1,nkey
358 itri(j,i1) = iwork_t(ii)%NPT_PLY(j)
359 ENDDO
360 ENDDO
361 mode = 0
362 CALL my_orders(mode, work, itri, indx1, nsh , nkey)
363 ALLOCATE (nfirst2(nsh) ,nlast2(nsh))
364 ns_sub = 1
365 nfirst2(1) = 1
366 nlast2(1) = 1
367 DO i=2,nsh
368 DO ikey = 1, nkey
369 ii = itri(ikey,indx1(i))
370 jj = itri(ikey,indx1(i-1))
371 IF(ii /= jj) THEN
372 ns_sub = ns_sub + 1
373 nfirst2(ns_sub) = i
374 nlast2(ns_sub) = nfirst2(ns_sub)
375 EXIT
376 ELSEIF(ikey == nkey) THEN
377 nlast2(ns_sub) = nlast2(ns_sub) + 1
378 ENDIF
379 ENDDO
380 ENDDO
381 !!
382 DO iis = 1,ns_sub
383 ns_total = ns_total + 1
384 nfirst(ns_total) = nfirst1(is) + nfirst2(iis) - 1
385 nlast(ns_total ) = nfirst1(is) + nlast2(iis) - 1
386 DO i = nfirst2(iis),nlast2(iis)
387 i2 = nfirst1(is) + i - 1
388 indx_total(i2) = indx_sub(indx1(i)) !
389 ENDDO
390 ENDDO
391 DEALLOCATE(indx1,nfirst2,nlast2, itri,indx_sub)
392 ELSE
393 ns_total = ns_total + 1
394 nfirst(ns_total) = nfirst1(is)
395 nlast(ns_total ) = nlast1(is)
396 DO i= nfirst1(is), nlast1(is)
397 indx_total(i) = indx(i)
398 ENDDO
399 ENDIF
400 END DO
401 DEALLOCATE(nfirst1,nlast1)
402C substack
403C
404 npt_stack = 0
405 ns_stack = ns_total
406C
407 DO is = 1,ns_total
408 id = nfirst(is)
409 i = indx_total(id)
410 ii = indx_sh(i)
411 npt = iworksh(1,ii)
412 npt_stack = max(npt_stack,npt)
413 ENDDO
414C allocation
415 ALLOCATE(iworks%IGEO(4*npt_stack+2,ns_stack))
416 ALLOCATE(iworks%GEO(6*npt_stack+1,ns_stack))
417C
418 iworks%IGEO = 0
419 iworks%GEO = zero
420C
421 DO is = 1,ns_total
422 ngeo_stack = numgeo + is
423 id = nfirst(is)
424C
425 i = indx_total(id)
426 ii = indx_sh(i)
427 pid = pid_sh(i)
428!! MAIN_PID(1,IS) = PID_SH(II) ! pid
429!! MAIN_PID(2,IS) = IWORKSH(1,II) ! npt
430 npt = iworksh(1,ii)
431!C ISUBSTACK(PID) = ISUBSTACK(PID) + 1
432 ies = ii
433C
434 DO i= nfirst(is) , nlast(is)
435 id = indx_total(i)
436 ii = indx_sh(id)
437 iworksh(2,ii) = ngeo_stack
438 iworksh(3,ii) = is ! compter for all stack
439!! IWORKSH(3,II) = ISUBSTACK(PID) ! computer by stack ! old conception
440 ENDDO
441C Geo and Igeo organisation
442!! I = NFIRST(IS)
443!! II = INDX_Q4(I)
444!! PID = PID_Q4(II)
445!! NSS = ISUBSTACK(PID) ! number of substack in each Pid
446C geometric property for each sub stack
447!!! IGEO( ,PID) = NGEO_STACK
448 n1 = int(geo(6,pid))
449 np = 0
450 nums = numgeo_stack(pid)
451 DO 700 j = 1,n1
452!! JPID = IGEO(100 + J, PID)
453 jpid = stack_info(nums)%PID(j)
454 IF(np <= npt) THEN
455 DO jj = 1,npt
456 jjpid = iwork_t(ies)%PLYID(jj)
457 IF(jjpid == jpid) THEN
458 np = np + 1
459 iptply(np) = j
460 GOTO 700
461 ENDIF
462 ENDDO
463 ENDIF
464 700 CONTINUE
465C geometric property
466 iworks%IGEO(1,is) = npt
467 iworks%IGEO(2,is) = pid
468 ippid = 2
469 ipmat = ippid + npt
470 ipmat_iply = ipmat + npt
471 ipnpt_lay = ipmat_iply + npt
472 ipang = 1
473 ipthk = ipang + npt
474 ippos = ipthk + npt
475 ipdir = ippos + npt
476 ipthkly = ipdir + npt
477 ipweight = ipthkly + npt
478 nums= numgeo_stack(pid)
479 DO j=1,npt
480 jstack = iptply(j)
481 iworks%IGEO(ippid + j,is) = stack_info(nums)%PID(jstack)
482 iworks%IGEO(ipmat + j,is) = stack_info(nums)%MID(jstack)
483 iworks%IGEO(ipmat_iply + j,is) = stack_info(nums)%MID_IP(jstack)
484 iworks%IGEO(ipnpt_lay + j,is) = iwork_t(ies)%NPT_PLY(j) ! npt per laye
485 iworks%GEO(ipang + j,is) = stack_info(nums)%ANG(jstack)
486 iworks%GEO(ipthk + j,is) = stack_info(nums)%THK(jstack)
487 iworks%GEO(ippos + j,is) = stack_info(nums)%POS(jstack)
488 iworks%GEO(ipdir + j,is) = stack_info(nums)%DIR(jstack)
489 iworks%GEO(ipthkly + j,is) = stack_info(nums)%THKLY(jstack)
490 iworks%GEO(ipweight + j,is) = stack_info(nums)%WEIGHT(jstack)
491 ENDDO
492
493C
494C to be clarified later IPOS > 0
495 ipos = igeo(99,pid)
496 zshift = geo(199,pid)
497 IF(ipos == 1)THEN
498 tmin = ep20
499 tmax = -ep20
500 DO j=1,npt
501 dt = half*iworks%GEO(ipthk + j ,is)
502 tmin = min(tmin,iworks%GEO(ippos + j ,is)-dt)
503 tmax = max(tmax,iworks%GEO(ippos + j ,is)+dt)
504 ENDDO
505 thickt = tmax - tmin
506 DO j=1,npt
507 iworks%GEO(ipthk+j,is)=iworks%GEO(ipthk+j,is)/max(thickt,em20)
508 iworks%GEO(ippos+j,is)=iworks%GEO(ippos+j,is)/max(thickt,em20)
509 ENDDO
510
511 ELSE
512 thickt = zero
513 DO j=1,npt
514 thickt = thickt + iworks%GEO(ipthk+j,is)
515 ENDDO
516 DO j=1,npt
517 iworks%GEO(ipthk+j,is) =
518 . iworks%GEO(ipthk+j,is)/max(thickt,em20)
519 ENDDO
520C
521 IF(ipos == 2 ) zshift = zshift /max(thickt,em20)
522C--- calcul automatique de position des couches
523 iworks%GEO(ippos+1,is) = zshift + half*iworks%GEO(ipthk+1,is)
524 DO j=2,npt
525 iworks%GEO(ippos+j,is) = iworks%GEO(ippos+j-1,is)
526 . + half*(iworks%GEO(ipthk+j,is)+iworks%GEO(ipthk+j-1,is))
527 ENDDO
528C
529 ENDIF ! IPOS =0,1,2,3,4
530C calcul du thk by shell
531 iworks%GEO(1,is) = thickt
532C============================================================================
533C============================================================================
534 ENDDO ! DO IS = 1,NS_TOTAL
535C
536 DEALLOCATE(indx,nfirst,nlast,indx_total,activ_ply)
537 ENDIF
538C
539C pccommp part
540C
541C---
542 ns_stack0 = ns_stack
543 npt_stack0 = npt_stack
544C
545 IF(ipart_pcompp > 0) THEN
546!!!------------------------------------------
547 nply = 0
548 nstack = 0
549 DO i = 1, numply
550 ids = igeo_stack(42,numstack + i)
551 IF (ids > 0) nply = nply+1
552 ENDDO
553 !
554 nsh = 0
555 DO i=1,numelc
556 pid = ixc(6,i)
557 igtyp = igeo(11,pid)
558 IF(igtyp == 52)THEN
559 nsh = nsh + 1
560 indx_sh(nsh) = i
561 pid_sh(nsh) = pid
562 ENDIF
563 ENDDO
564C
565 DO i=1,numeltg
566 pid = ixtg(5,i)
567 igtyp = igeo(11,pid)
568 IF(igtyp == 52)THEN
569 nsh = nsh + 1
570 indx_sh(nsh) = i + numelc
571 pid_sh(nsh) = pid
572 ENDIF
573 ENDDO
574C #####################################################"
575 nbit = bit_size(nply)
576 irest = mod(nply,nbit)
577 nkey = nply / nbit
578 IF(irest > 0) nkey = nkey + 1
579 ALLOCATE( activ_ply(numelc+numeltg,nkey))
580 IF(numelc + numeltg > 0)activ_ply = 0
581 DO i=1,nsh
582 ii = indx_sh(i)
583 npt = iworksh(1,ii)
584 ie_drp = drapeg%INDX(ii)
585 ALLOCATE(iwork_t(ii)%NPT_PLY(npt))
586 iwork_t(ii)%NPT_PLY = 0
587 IF(ie_drp > 0 .AND. npt > 0) THEN
588 DO j=1,npt
589 ip = drape(ie_drp)%INDX_PLY(j)
590 n_ply = iwork_t(ii)%PLYNUM(j)
591 ikey = n_ply / nbit ! 32 or 64 bits
592 IF(mod(n_ply,nbit) > 0 ) ikey = ikey + 1
593 ikey = min(ikey, nkey)
594 ibit = n_ply - (ikey - 1)*nbit ! Bit corresponding to PLy
595 activ_ply(ii,ikey) = ibset(activ_ply(ii,ikey),ibit)
596 IF(ip > 0 ) THEN
597 drape_ply => drape(ie_drp)%DRAPE_PLY(ip)
598 nslice = drape_ply%NSLICE
599 idply = iwork_t(ii)%PLYID(j)
600 iwork_t(ii)%NPT_PLY(j) = nslice
601 ply_info(2,idply - numstack) = max(ply_info(2,idply - numstack),nslice)
602 ELSE
603 idply = iwork_t(ii)%PLYID(j)
604 npt_lay = igeo_stack(4,idply)
605 iwork_t(ii)%NPT_PLY(j) = npt_lay
606 ENDIF
607 ENDDO
608 ELSE
609 IF(npt > 0) THEN
610 DO j=1,npt
611 n_ply = iwork_t(ii)%PLYNUM(j)
612 ikey = n_ply / nbit ! 32 or 64 bits
613 IF(mod(n_ply,nbit) > 0 ) ikey = ikey + 1
614 ikey = min(ikey, nkey)
615 ibit = n_ply - (ikey - 1)*nbit ! Bit corresponding to PLy
616 activ_ply(ii,ikey) = ibset(activ_ply(ii,ikey),ibit)
617 idply = iwork_t(ii)%PLYID(j)
618 npt_lay = igeo_stack(4,idply)
619 iwork_t(ii)%NPT_PLY(j) = npt_lay
620 ENDDO
621 ENDIF
622 ENDIF
623 ENDDO
624!!###########################################################"
625 ALLOCATE(indx(2*nsh),itri(2 + nkey,nsh),indx_total(nsh))
626 ALLOCATE (nfirst(nsh) ,nlast(nsh),
627 . nfirst1(nsh) ,nlast1(nsh))
628 IF(nsh > 0) THEN
629 indx = 0
630 itri = 0
631 indx_total = 0
632 nlast = 0
633 nfirst = 0
634 nfirst1 = 0
635 nlast1 = 0
636 ENDIF
637 !! sorting depending on number of number & type of ply
638 DO i = 1,nsh
639 indx(i) = i
640 ii = indx_sh(i)
641 itri(1,i) = pid_sh(i)
642 itri(2,i) = iworksh(1,ii)
643 DO j=1,nkey
644 itri(2+j,i) = activ_ply(ii,j)
645 ENDDO
646 ENDDO
647C
648 mode = 0
649C
650 nkey = nkey + 2
651C
652 CALL my_orders(mode, work, itri, indx, nsh , nkey)
653 ns_first = 1
654 nfirst1(1) = 1
655 nlast1(1) = 1
656 DO i=2,nsh
657 DO ikey = 1, nkey
658 ii = itri(ikey,indx(i))
659 jj = itri(ikey,indx(i-1))
660 IF(ii /= jj) THEN
661 ns_first = ns_first + 1
662 nfirst1(ns_first) = i
663 nlast1(ns_first) = nfirst1(ns_first)
664 EXIT
665 ELSEIF(ikey == nkey) THEN
666 nlast1(ns_first) = nlast1(ns_first) + 1
667 ENDIF
668 ENDDO
669 ENDDO
670 DEALLOCATE (itri)
671 !! sorting depending on number of slice by ply = int point by ply for Type51
672 ns_total = 0
673 ns_sub = 0
674 !!
675 DO is =1,ns_first
676 id = nfirst1(is)
677 i = indx(id)
678 ii = indx_sh(i)
679 pid = pid_sh(i)
680 igtyp = igeo(11,pid)
681 nkey = iworksh(1,ii)
682 nsh = nlast1(is) - nfirst1(is) + 1
683 ALLOCATE(indx1(2*nsh),itri(nkey,nsh),indx_sub(nsh))
684 indx1 = 0
685 itri = 0
686 i1 = 0
687 DO i= nfirst1(is), nlast1(is)
688 id = indx(i)
689 ii = indx_sh(id)
690 i1 = i1 + 1
691 indx1(i1) = i1
692 indx_sub(i1) = id
693 DO j=1,nkey
694 itri(j,i1) = iwork_t(ii)%NPT_PLY(j)
695 ENDDO
696 ENDDO
697 mode = 0
698 CALL my_orders(mode, work, itri, indx1, nsh , nkey)
699 ALLOCATE (nfirst2(nsh) ,nlast2(nsh))
700 ns_sub = 1
701 nfirst2(1) = 1
702 nlast2(1) = 1
703 DO i=2,nsh
704 DO ikey = 1, nkey
705 ii = itri(ikey,indx1(i))
706 jj = itri(ikey,indx1(i-1))
707 IF(ii /= jj) THEN
708 ns_sub = ns_sub + 1
709 nfirst2(ns_sub) = i
710 nlast2(ns_sub) = nfirst2(ns_sub)
711 EXIT
712 ELSEIF(ikey == nkey) THEN
713 nlast2(ns_sub) = nlast2(ns_sub) + 1
714 ENDIF
715 ENDDO
716 ENDDO
717 !!
718 DO iis = 1,ns_sub
719 ns_total = ns_total + 1
720 nfirst(ns_total) = nfirst1(is) + nfirst2(iis) - 1
721 nlast(ns_total ) = nfirst1(is) + nlast2(iis) - 1
722 DO i = nfirst2(iis),nlast2(iis)
723 i2 = nfirst1(is) + i - 1
724 indx_total(i2) = indx_sub(indx1(i)) !
725 ENDDO
726 ENDDO
727 DEALLOCATE(indx1,nfirst2,nlast2, itri,indx_sub)
728 END DO
729 DEALLOCATE(nfirst1,nlast1)
730C
731C sous stack
732C
733 ALLOCATE(idstack(ns_total))
734 idstack = 0
735 ns_stack = ns_stack + ns_total
736 DO is = 1,ns_total
737 id = nfirst(is)
738 i = indx_total(id)
739 ii = indx_sh(i)
740 npt = iworksh(1,ii)
741 npt_stack = max(npt_stack,npt)
742C
743 ids = iwork_t(ii)%IDSTACK
744 idstack(is) = ids
745 ENDDO
746C
747C allocation
748C
749 ALLOCATE(stack%IGEO(4*npt_stack+2,ns_stack))
750 ALLOCATE(stack%GEO(6*npt_stack+1,ns_stack))
751 ALLOCATE(stack%PM(20,ns_stack))
752C
753 stack%IGEO = 0
754 stack%GEO = zero
755 stack%PM = zero
756C
757 DO is = 1,ns_total
758C a changer
759 ngeo_stack = numgeo + numstack + numply + is !!!!!!! limit id i will change it
760 id = nfirst(is)
761C
762 i = indx_total(id)
763 ii = indx_sh(i)
764 pid = pid_sh(i)
765!! MAIN_PID(1,IS) = PID_SH(II) ! pid
766!! MAIN_PID(2,IS) = IWORKSH(1,II) ! npt
767 npt = iworksh(1,ii)
768 ies = ii
769C
770 DO i= nfirst(is) , nlast(is)
771 id = indx_total(i)
772 ii = indx_sh(id)
773 iworksh(2,ii) = ngeo_stack
774 iworksh(3,ii) = ns_stack0 + is ! compter for all stack
775 ENDDO
776C Cp igeo_stack and Geo_stack dans IGEO, GEO --ppccomp
777 igtyp = igeo(11,pid)
778 DO j=2,npropgi - ltitr
779 igeo(j,pid) = igeo_stack(j,idstack(is))
780 ENDDO
781 igeo(11,pid) = igtyp
782!
783 DO j=1,npropg
784 geo(j,pid) = geo_stack(j,idstack(is))
785 ENDDO
786C
787 n1 = int(geo(6,pid))
788 np = 0
789 nums = numgeo_stack(numgeo + idstack(is))
790 DO 777 j = 1,n1
791 jpid = stack_info(nums)%PID(j)
792 IF(np <= npt) THEN
793 DO jj = 1,npt
794 jjpid = iwork_t(ies)%PLYID(jj)
795 IF(jjpid == jpid) THEN
796 np = np + 1
797 iptply(np) = j
798 GOTO 777
799 ENDIF
800 ENDDO
801 ENDIF
802 777 CONTINUE
803C geometric property
804C
805 iis = ns_stack0 + is
806 stack%IGEO(1,iis) = npt
807 stack%IGEO(2,iis) = pid
808 ippid = 2
809 ipmat = ippid + npt
810 ipmat_iply = ipmat + npt
811 ipnpt_lay = ipmat_iply + npt
812C
813 ipang = 1
814 ipthk = ipang + npt
815 ippos = ipthk + npt
816 ipdir = ippos + npt
817 ipthkly = ipdir + npt
818 ipweight =ipthkly + npt
819C stack id
820 pids = idstack(is)
821 nums = numgeo_stack(numgeo + pids)
822 DO j=1,npt
823 js = iptply(j)
824 stack%IGEO(ippid+j ,iis) = stack_info(nums)%PID(js)
825 stack%IGEO(ipmat + j ,iis) = stack_info(nums)%MID(js)
826 stack%IGEO(ipmat_iply+j ,iis) = stack_info(nums)%MID_IP(js)
827 stack%IGEO(ipnpt_lay + j,iis) = iwork_t(ies)%NPT_PLY(j) ! npt per laye
828 stack%GEO(ipang + j ,iis) = stack_info(nums)%ANG(js)
829 stack%GEO(ipthk + j ,iis) = stack_info(nums)%THK(js)
830 stack%GEO(ippos + j ,iis) = stack_info(nums)%POS(js)
831 stack%GEO(ipdir + j ,iis) = stack_info(nums)%DIR(js)
832 stack%GEO(ipthkly + j ,iis) = stack_info(nums)%THKLY(js)
833 stack%GEO(ipweight + j ,iis) = stack_info(nums)%WEIGHT(js)
834 ENDDO
835C
836C to be clarified later IPOS > 0
837 ipos = igeo(99,pid)
838 zshift = geo(199,pid)
839 IF(ipos == 1)THEN
840 tmin = ep20
841 tmax = -ep20
842 DO j=1,npt
843 dt = half*stack%GEO(ipthk + j ,iis)
844 tmin = min(tmin,stack%GEO(ippos + j ,iis)-dt)
845 tmax = max(tmax,stack%GEO(ippos + j ,iis)+dt)
846 ENDDO
847 thickt = tmax - tmin
848 DO j=1,npt
849 stack%GEO(ipthk+j,iis)=
850 . stack%GEO(ipthk+j,iis)/max(thickt,em20)
851 stack%GEO(ippos+j,iis)=
852 . stack%GEO(ippos+j,iis)/max(thickt,em20)
853 ENDDO
854
855 ELSE
856 thickt = zero
857 DO j=1,npt
858 thickt = thickt + stack%GEO(ipthk+j,iis)
859 ENDDO
860 DO j=1,npt
861 stack%GEO(ipthk+j,iis) =
862 . stack%GEO(ipthk+j,iis)/max(thickt,em20)
863 ENDDO
864C
865 IF(ipos == 2 )zshift = zshift /max(thickt,em20)
866C--- calcul automatique de position des couches
867 stack%GEO(ippos+1,iis) = zshift +
868 . half*stack%GEO(ipthk+1,iis)
869 DO j=2,npt
870 stack%GEO(ippos+j,iis) =
871 . stack%GEO(ippos+j-1,iis) +
872 . half*(stack%GEO(ipthk+j,iis)+
873 . stack%GEO(ipthk+j-1,iis))
874 ENDDO
875C
876 ENDIF ! IPOS =0,1,2,3,4
877C calcul du thk by shell
878 stack%GEO(1,iis) = thickt
879C============================================================================
880C---
881 ippid = 2
882 DO ilay=1,npt
883 pids = stack%IGEO(ippid + ilay ,iis)
884 nptt = igeo_stack(4,pids)
885 igeo(4,pid) = max(igeo(4,pid),nptt)
886 ENDDO
887 ENDDO ! DO IS = 1,NS
888C---
889 DEALLOCATE(indx,nfirst,nlast,indx_total,idstack,activ_ply)
890 ENDIF
891C
892 DO i=1,numelc + numeltg
893 npt = iworksh(1,i)
894 IF(npt > 0) THEN
895 DEALLOCATE(iwork_t(i)%PLYID)
896 DEALLOCATE(iwork_t(i)%NPT_PLY)
897 ENDIF
898 ENDDO
899 IF(ipart_stack > 0) THEN
900 IF(ipart_pcompp == 0) THEN
901 ALLOCATE(stack%IGEO(4*npt_stack0+2,ns_stack0))
902 ALLOCATE(stack%GEO(6*npt_stack0+1,ns_stack0))
903 ALLOCATE(stack%PM(20,ns_stack0))
904 stack%IGEO = 0
905 stack%GEO = zero
906 stack%PM = zero
907 ENDIF
908 DO is = 1, ns_stack0
909 DO j = 1, 4*npt_stack0 + 2
910 stack%IGEO(j, is ) = iworks%IGEO(j,is)
911 ENDDO
912 DO j = 1, 6*npt_stack0+1
913 stack%GEO(j, is ) = iworks%GEO(j,is)
914
915 ENDDO
916 ENDDO
917 DEALLOCATE(iworks%IGEO, iworks%GEO)
918 ENDIF
919C --- update of sub-stack
920 IF(ns_stack > 0) THEN
921 DO is = 1,ns_stack
922 npt = stack%IGEO(1,is)
923 pid = stack%IGEO(2,is)
924 thickt = stack%GEO(1,is)
925 id = igeo(1,pid)
926 igtyp = igeo(11,pid)
927C
928 WRITE(iout,1000)id, is
929 WRITE(iout,1100) thickt,npt
930!! IPANG = 1
931!! IPTHK = IPANG + NPT
932!! IPPOS = IPTHK + NPT
933 ippos = 1 + 2*npt
934 ippid = 2
935 IF(igtyp == 52) THEN
936 DO j = 1,npt
937 pid = stack%IGEO(ippid + j ,is)
938 pos = stack%GEO( ippos + j ,is)
939 pos = pos*thickt
940 id = igeo_stack(1,pid)
941 WRITE(iout,2000)j, id , pos
942 ENDDO
943 ELSE
944 DO j = 1,npt
945 pid = stack%IGEO(ippid + j ,is)
946 pos = stack%GEO( ippos + j ,is)
947 pos = pos*thickt
948 id = igeo(1,pid)
949 WRITE(iout,2000)j, id , pos
950 ENDDO
951 ENDIF
952 ENDDO
953 ENDIF
954C for restart
955 IF(ipart_pcompp > 0 .AND. ipart_stack == 0) ipart_stack = 1
956
957 DEALLOCATE (indx_sh,pid_sh,isubstack,
958 . iptply, work )
959C--------
960 RETURN
961 1000 FORMAT(//,
962 & 5x,'COMPOSITE STACK SHELL PROPERTY SET ',
963 & 'WITH VARIABLE THICKNESSES AND MATERIALS'//,
964 & 7x,'PROPERTY SET NUMBER . . . . . . . . . . ..=',i10/,
965 & 7x,'SUB PROPERTY SET NUMBER . . . . . . . . . .=',i10/)
966 1100 FORMAT(
967 & 8x,'SHELL THICKNESS . . . . . . . . . . . .=',1pg20.13/
968 & 8x,'NUMBER OF PLIES. . . . . . . . . . . . =',i10/)
969 2000 FORMAT(
970 & 8x,' PLY ',i3/,
971 & 8x,' PLY PID NUMBER . . . . . . . . .=',i10/
972 & 8x,' POSITION. . . . . . . . . . . . .=',1pg20.13/)
973
974 END
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21
my_orders
void my_orders(int *mode, int *iwork, int *data, int *index, int *n, int *irecl)
Definition my_orders.c:82
drape_mod
Definition drape_mod.F:54
groupdef_mod
Definition groupdef_mod.F:662
message_mod
Definition message_mod.F:1249
names_and_titles_mod
Definition names_and_titles_mod.F:997
names_and_titles_mod::nchartitle
integer, parameter nchartitle
Definition names_and_titles_mod.F:1003
stack_mod
Definition stack_mod.F:72
submodel_mod
Definition submodel_mod.F:712
stackgroup_drape
subroutine stackgroup_drape(drape, drapeg, iwork_t, iworksh, igrsh3n, igrsh4n, ixc, ixtg, igeo, geo, thk, stack, igeo_stack, geo_stack, stack_info, numgeo_stack, nprop_stack, ply_info)
Definition stackgroup_drape.F:39
drape_mod::drape_work_
Definition drape_mod.F:85