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