OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
thcoq.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!|| thcoq ../engine/source/output/th/thcoq.F
25!||--- called by ------------------------------------------------------
26!|| hist2 ../engine/source/output/th/hist2.F
27!||--- uses -----------------------------------------------------
28!|| elbufdef_mod ../common_source/modules/mat_elem/elbufdef_mod.F90
29!|| matparam_def_mod ../common_source/modules/mat_elem/matparam_def_mod.F90
30!|| pinchtype_mod ../common_source/modules/pinchtype_mod.F
31!|| stack_mod ../engine/share/modules/stack_mod.F
32!||====================================================================
33 SUBROUTINE thcoq(ELBUF_TAB,MATPARAM_TAB,NTHGRP2 , ITHGRP ,
34 . IPARG ,ITHBUF ,WA ,
35 . IPM ,IGEO ,IXC ,IXTG ,PM ,
36 . RTHBUF ,THKE ,STACK)
37C-----------------------------------------------
38C M o d u l e s
39C-----------------------------------------------
40 USE elbufdef_mod
41 USE stack_mod
42 USE pinchtype_mod
43 USE matparam_def_mod
44C-----------------------------------------------
45C I m p l i c i t T y p e s
46C-----------------------------------------------
47#include "implicit_f.inc"
48C-----------------------------------------------
49C C o m m o n B l o c k s
50C-----------------------------------------------
51#include "mvsiz_p.inc"
52#include "com01_c.inc"
53#include "com04_c.inc"
54#include "task_c.inc"
55#include "param_c.inc"
56C-----------------------------------------------
57C D u m m y A r g u m e n t s
58C-----------------------------------------------
59 INTEGER IPARG(NPARG,*),ITHBUF(*),IXC(NIXC,*),
60 . IXTG(NIXTG,*),IPM(NPROPMI,*),IGEO(NPROPGI,*)
61 INTEGER, INTENT(in) :: NTHGRP2
62 INTEGER, DIMENSION(NITHGR,*), INTENT(in) :: ITHGRP
63 my_real
64 . wa(*),pm(npropm,*),rthbuf(*),thke(*)
65 TYPE (ELBUF_STRUCT_), DIMENSION(NGROUP), TARGET :: ELBUF_TAB
66 TYPE (MATPARAM_STRUCT_) ,DIMENSION(NUMMAT) ,INTENT(IN) :: MATPARAM_TAB
67C-----------------------------------------------
68C L o c a l V a r i a b l e s
69C-----------------------------------------------
70 INTEGER I,J,K,L,II,JJ,N, IH, NG, ITY, MTE, M2, M3, M5,M8,
71 . NPT,MPT,NPG,NPTR,NPTS,NPTT,NLAY,IP,IR,IS,IT,IL,IPT,
72 . LWA,NEL,NFT,I1,I2,I3,I4,IUV,IAA,IADR,N16,N16A,
73 . istrain,nu,nuvar,nuvarv,nuvard,igtyp,ihbe,nbd1,nbd2,nbd3,
74 . ifailure,iadd,isrot,ivisc,ipmat,ptmat,ishplyxfem,ipmat_iply,
75 . mat_iply,nbdelm,iwa,nv,ngl,iigeo,iadi,isubstack,ithk,npt_all,
76 . matly,kk(8),ipinch,ipg,imat,mat_orth, idrape
77 INTEGER PID(MVSIZ),MAT(MVSIZ)
78 INTEGER :: NITER,IAD,NN,IADV,NVAR,ITYP,IJK
79 my_real :: WWA(50000),FUNC(6),SIG(5),SIGG(5)
80 my_real ,DIMENSION(MVSIZ) :: dam1,dam2,wpla,dmax,wpmax,
81 . fail,fail1,fail2,fail3
82 my_real :: f1,f2,f3,f4,f5,f11,f22,f33,f44,f55,cp,sp,mm1,mm2,mm3,
83 . mm11,mm22,mm33,d1,d2,d11,d12,d22,val_ly_ip,val_ly_average
84 TYPE(g_bufel_) ,POINTER :: GBUF
85 TYPE(l_bufel_) ,POINTER :: LBUF
86 TYPE(BUF_LAY_) ,POINTER :: BUFLY
87 my_real ,DIMENSION(:), POINTER :: uvar,dir_a
88 my_real ,DIMENSION(:,:), ALLOCATABLE :: var
89 TYPE (STACK_PLY) :: STACK
90C-------------------------
91C ELEMENTS COQUES
92C=======================================================================
93 ijk = 0
94 DO niter=1,nthgrp2
95 ityp=ithgrp(2,niter)
96 nn =ithgrp(4,niter)
97 iad =ithgrp(5,niter)
98 nvar=ithgrp(6,niter)
99 iadv=ithgrp(7,niter)
100 ii=0
101 IF(ityp==3.OR.ityp==7)THEN
102! -------------------------------
103 ii=0
104 ih=iad
105C specifique spmd
106C decalage IH
107 DO WHILE((ithbuf(ih+nn)/=ispmd).AND.(ih<iad+nn))
108 ih = ih + 1
109 ENDDO
110 IF (ih>=iad+nn) GOTO 666
111C-------------------
112 DO ng=1,ngroup
113 ity=iparg(5,ng)
114 IF (ity == ityp) THEN
115 mte=iparg(1,ng)
116 nel=iparg(2,ng)
117 nft=iparg(3,ng)
118 npt = iparg(6,ng)
119 igtyp = iparg(38,ng)
120 istrain = iparg(44,ng)
121 ihbe = iparg(23,ng)
122 ifailure = iparg(43,ng)
123 ishplyxfem = iparg(50,ng)
124 isubstack = iparg(71,ng)
125 ithk =iparg(28,ng)
126 gbuf => elbuf_tab(ng)%GBUF
127 nptr = elbuf_tab(ng)%NPTR
128 npts = elbuf_tab(ng)%NPTS
129 nptt = elbuf_tab(ng)%NPTT
130 nlay = elbuf_tab(ng)%NLAY
131 idrape = elbuf_tab(ng)%IDRAPE
132 npg = nptr*npts
133cc NPT = NLAY*NPTT ! not compatible with PID51 (shell)
134 mpt = max(1,npt)
135!
136 DO i=1,8 ! length max of GBUF%G_STRA = 8
137 kk(i) = nel*(i-1)
138 ENDDO
139!
140C
141 IF (igtyp == 51 .OR. igtyp == 52) THEN
142 npt_all = 0
143 DO ipt=1,nlay
144 npt_all = npt_all + elbuf_tab(ng)%BUFLY(ipt)%NPTT
145 ENDDO
146 IF (nlay == 1) mpt = max(1,npt_all)
147 ENDIF
148C
149 ivisc = 0
150 nuvar = 0
151 nuvarv = 0
152 nuvard = 0
153c
154 IF (mte /= 13 .and. mte /= 0) THEN
155c
156 IF ((mte>=29.AND.mte<=31).OR.
157 . mte == 35.OR.mte == 36.OR.mte == 43.OR.
158 . mte == 44.OR.mte == 45.OR.mte == 48.OR.mte>=50) THEN
159 CONTINUE
160c
161 ELSEIF (mte == 25) THEN
162C
163 DO i=1,nel
164 dam1(i)=zero
165 dam2(i)=zero
166 wpla(i)=zero
167 fail(i)=zero
168 fail1(i)=zero
169 fail2(i)=zero
170 fail3(i)=zero
171 ENDDO
172c
173 IF (ity == 3) THEN
174 DO i=1,nel
175 mat(i)=ixc(1,nft+i)
176 pid(i)=ixc(6,nft+i)
177 ENDDO
178 ELSE
179 DO i=1,nel
180 mat(i)=ixtg(1,nft+i)
181 pid(i)=ixtg(5,nft+i)
182 ENDDO
183 ENDIF
184c---
185 IF (igtyp == 11) THEN
186 ipmat = 100
187 DO n=1,mpt
188 DO i=1,nel
189 matly = igeo(ipmat+n,pid(i))
190 IF (matparam_tab(matly)%IVISC > 0) THEN
191 ivisc = 1
192 nuvarv = max(nuvarv, matparam_tab(matly)%VISC%NUVAR)
193 END IF
194 ENDDO
195 ENDDO
196 ELSEIF (igtyp == 9 .OR. igtyp == 10) THEN
197 DO n=1,mpt
198 DO i=1,nel
199 matly=mat(i)
200 IF (matparam_tab(matly)%IVISC > 0) THEN
201 ivisc = 1
202 nuvarv = max(nuvarv, matparam_tab(matly)%VISC%NUVAR)
203 END IF
204 ENDDO
205 ENDDO
206 ELSEIF (igtyp == 17 .OR. igtyp == 51 .OR. igtyp == 52) THEN
207 ipmat = 2 + nlay
208 DO n=1,nlay
209 DO i=1,nel
210 matly = stack%IGEO(ipmat+n,isubstack)
211 IF (matparam_tab(matly)%IVISC > 0) THEN
212 ivisc = 1
213 nuvarv = max(nuvarv, matparam_tab(matly)%VISC%NUVAR)
214 END IF
215 ENDDO
216 ENDDO
217c
218 IF (ishplyxfem > 0) THEN
219 ipmat_iply = ipmat + mpt
220 DO j=1,npt -1
221 DO i=1,nel
222 mat_iply = stack%IGEO(ipmat_iply + j ,isubstack)
223 nuvard = max(nuvard, ipm(221,mat_iply))
224 ENDDO
225 ENDDO
226 ENDIF
227 ENDIF
228 ENDIF ! MTE
229c---------
230c
231c---------
232 DO i=1,nel
233 n=i+nft
234 k=ithbuf(ih)
235 ip=ithbuf(ih+nn)
236 iadr=ithbuf(ih+3*nn) ! Adress of cos sin related to skew
237C
238 IF (k == n) THEN
239 ih=ih+1
240C traitement specifique spmd
241C recherche du ii correct
242 ii = ((ih-1) - iad)*nvar
243 DO WHILE((ithbuf(ih+nn) /= ispmd) .AND. (ih < iad+nn))
244 ih = ih + 1
245 ENDDO
246C
247 IF (ih > iad+nn) GOTO 666
248C
249 m5=5*(i-1)
250 m8=8*(i-1)
251c
252 IF (iadr /= 0) THEN ! output with respect to a (non global) SKEW
253 cp=rthbuf(iadr)
254 sp=rthbuf(iadr+1)
255c
256 f11 = gbuf%FOR(kk(1)+i)
257 f22 = gbuf%FOR(kk(2)+i)
258 f33 = gbuf%FOR(kk(3)+i)
259 f44 = gbuf%FOR(kk(4)+i)
260 f55 = gbuf%FOR(kk(5)+i)
261c
262 mm11 = gbuf%MOM(kk(1)+i)
263 mm22 = gbuf%MOM(kk(2)+i)
264 mm33 = gbuf%MOM(kk(3)+i)
265c
266 f1 = cp*cp*f11
267 . + sp*sp*f22
268 . + two*cp*sp*f33
269c
270 f2 = sp*sp*f11
271 . + cp*cp*f22
272 . - two*cp*sp*f33
273c
274 f3 =-cp*sp*f11
275 . + cp*sp*f22
276 . + (cp*cp-sp*sp )*f33
277c
278 f4 =-sp*f55+cp*f44
279 f5 = cp*f55+sp*f44
280c
281 mm1 = cp*cp*mm11
282 . + sp*sp*mm22
283 . + two*cp*sp*mm33
284c
285 mm2 = sp*sp*mm11
286 . + cp*cp*mm22
287 . - two*cp*sp*mm33
288c
289 mm3 =-cp*sp*mm11
290 . + cp*sp*mm22
291 . + (cp*cp-sp*sp )*mm33
292 ELSE !output with respect to the global SKEW.
293 f1 = gbuf%FOR(kk(1)+i)
294 f2 = gbuf%FOR(kk(2)+i)
295 f3 = gbuf%FOR(kk(3)+i)
296 f4 = gbuf%FOR(kk(4)+i)
297 f5 = gbuf%FOR(kk(5)+i)
298c
299 mm1 = gbuf%MOM(kk(1)+i)
300 mm2 = gbuf%MOM(kk(2)+i)
301 mm3 = gbuf%MOM(kk(3)+i)
302 ENDIF
303 wwa(1) = f1
304 wwa(2) = f2
305 wwa(3) = f3
306 wwa(4) = f4
307 wwa(5) = f5
308 wwa(6) = mm1
309 wwa(7) = mm2
310 wwa(8) = mm3
311 wwa(9) = gbuf%EINT(i)
312 wwa(10)= gbuf%EINT(i+nel)
313 wwa(11)= gbuf%OFF(i)
314 IF (ithk > 0) THEN
315 wwa(12)= gbuf%THK(i)
316 ELSE
317 wwa(12)= thke(n)
318 ENDIF
319 wwa(13)=zero
320 wwa(14)=zero
321 wwa(15)=zero
322 wwa(16)=zero
323 wwa(17)=zero
324 wwa(18)=zero
325 wwa(19)=zero
326 wwa(20)=zero
327 wwa(21)=zero
328 wwa(22)=zero
329 IF (gbuf%G_EPSD == 0) THEN
330 wwa(23)=zero
331 ELSE
332 wwa(23)=gbuf%EPSD(i)
333 ENDIF
334 DO j = 24,50000
335 wwa(j)=zero
336 ENDDO
337c----------------------
338c Stress tensor
339c----------------------
340c---- mean stress over Gauss points in each layer
341 IF(idrape > 0 .AND. (igtyp == 51 .OR. igtyp == 52) )THEN
342 DO il = 1,nlay
343 bufly => elbuf_tab(ng)%BUFLY(il)
344 imat = bufly%IMAT
345 nptt = bufly%NPTT
346 k = 183 + (il-1)*5
347 sigg(1:5) = zero
348 DO it=1,nptt
349 dir_a => elbuf_tab(ng)%BUFLY(il)%LBUF_DIR(it)%DIRA
350 sig(1:5) = zero
351 DO ir=1,nptr
352 DO is=1,npts
353 lbuf => bufly%LBUF(ir,is,it)
354 DO j = 1,5
355 sig(j) = sig(j) + lbuf%SIG(kk(j) + i) / npg
356 ENDDO
357 ENDDO
358 ENDDO
359 d1 = dir_a(i)
360 d2 = dir_a(i+nel)
361 d11 = d1*d1
362 d22 = d2*d2
363 d12 = d1*d2
364 sigg(1) = sigg(1) + (d11*sig(1) + d22*sig(2) + two*d12 *sig(3)) /nptt
365 sigg(2) = sigg(2) + (d22*sig(1) + d11*sig(2) - two*d12 *sig(3)) /nptt
366 sigg(3) = sigg(3) + (d12*sig(2) + (d11-d22)*sig(3) -d12*sig(1)) /nptt
367 sigg(4) = sigg(4) + (d1 *sig(4) - d2 *sig(5)) / nptt
368 sigg(5) = sigg(5) + (d1 *sig(5) + d2 *sig(4)) / nptt
369 ENDDO
370 wwa(k + 1) =sigg(1)
371 wwa(k + 2) =sigg(2)
372 wwa(k + 3) =sigg(3)
373 wwa(k + 4) =sigg(4)
374 wwa(k + 5) =sigg(5)
375 ENDDO ! DO IL=1,NLAY
376 ELSE
377 DO il = 1,nlay
378 bufly => elbuf_tab(ng)%BUFLY(il)
379 imat = bufly%IMAT
380 nptt = bufly%NPTT
381 sig(1:5) = zero
382 k = 183 + (il-1)*5
383 DO ir=1,nptr
384 DO is=1,npts
385 DO it=1,nptt
386 lbuf => bufly%LBUF(ir,is,it)
387 DO j = 1,5
388 sig(j) = sig(j) + lbuf%SIG(kk(j) + i) / (nptt*npg)
389 ENDDO
390 ENDDO
391 ENDDO
392 ENDDO
393 mat_orth = matparam_tab(imat)%ORTHOTROPY
394 IF (mat_orth == 1) THEN
395 DO j = 1,5
396 wwa(k + j) = sig(j)
397 ENDDO
398 ELSE IF (mat_orth == 2) THEN ! rotate sig to global coord
399 dir_a => elbuf_tab(ng)%BUFLY(il)%DIRA
400 d1 = dir_a(i)
401 d2 = dir_a(i+nel)
402 d11 = d1*d1
403 d22 = d2*d2
404 d12 = d1*d2
405 wwa(k + 1) = d11*sig(1) + d22*sig(2) + two*d12 *sig(3)
406 wwa(k + 2) = d22*sig(1) + d11*sig(2) - two*d12 *sig(3)
407 wwa(k + 3) =-d12*sig(1) + d12*sig(2) +(d11-d22)*sig(3)
408 wwa(k + 4) =-d2 *sig(5) + d1 *sig(4)
409 wwa(k + 5) = d1 *sig(5) + d2 *sig(4)
410 END IF
411 ENDDO ! DO IL=1,NLAY
412 ENDIF ! idrape
413c------------ Viscous stress
414c
415 DO il = 1,nlay
416 bufly => elbuf_tab(ng)%BUFLY(il)
417 imat = bufly%IMAT
418 ivisc = matparam_tab(imat)%IVISC
419 nptt = bufly%NPTT
420 IF (ivisc > 0) THEN
421 k = 30382+(il-1)*5
422 func(1:5) = zero
423 DO ir=1,nptr
424 DO is=1,npts
425 DO it=1,nptt
426 lbuf => bufly%LBUF(ir,is,it)
427 DO j = 1,5
428 func(j) = func(j) + lbuf%VISC(kk(j) + i) / nptt
429 ENDDO
430 ENDDO
431 DO j = 1,5
432 wwa(k+j) = func(j) / npg
433 ENDDO
434 ENDDO ! DO IS=1,NPTS
435 ENDDO ! DO IR=1,NPTR
436c
437 ENDIF ! IVISC > 0
438 ENDDO ! DO IL=1,NLAY
439c
440c------------ Viscous stress
441c
442c
443c------------ Max/Min Plastic strain
444 IF (gbuf%G_PLA > 0) THEN
445 wwa(13) = ep30
446 wwa(14) = zero
447c
448c IF (NPG == 1 .and. NLAY == 1) THEN
449c DO IPT = 1, NPT
450c LBUF => ELBUF_TAB(NG)%BUFLY(1)%LBUF(1,1,IPT)
451c WWA(13) = MIN(WWA(13),LBUF%PLA(I))
452c WWA(14) = MAX(WWA(14),LBUF%PLA(I))
453c ENDDO
454c ELSEIF (NLAY > 1) THEN
455c DO IPT=1,NLAY
456c BUFLY => ELBUF_TAB(NG)%BUFLY(IPT)
457c IF (BUFLY%L_PLA > 0) THEN
458c WWA(13) = MIN(WWA(13),ABS(BUFLY%PLAPT(I)))
459c WWA(14) = MAX(WWA(14),ABS(BUFLY%PLAPT(I)))
460c ENDIF
461c ENDDO
462c ELSEIF (ELBUF_TAB(NG)%BUFLY(1)%LY_PLAPT > 0) THEN
463c DO IPT = 1, NPT
464c BUFLY => ELBUF_TAB(NG)%BUFLY(1)
465c JJ = (IPT-1)*NEL
466c WWA(13) = MIN(WWA(13),ABS(BUFLY%PLAPT(JJ+I)))
467c WWA(14) = MAX(WWA(14),ABS(BUFLY%PLAPT(JJ+I)))
468c ENDDO
469c ENDIF
470 IF (nlay > 1) THEN
471 IF (npg > 1) THEN
472 DO ipt = 1,nlay
473 bufly => elbuf_tab(ng)%BUFLY(ipt)
474 IF (bufly%L_PLA > 0) THEN
475 wwa(13) = min(wwa(13),abs(bufly%PLAPT(i)))
476 wwa(14) = max(wwa(14),abs(bufly%PLAPT(i)))
477 ENDIF
478 ENDDO
479 ELSE ! (NPG = 1)
480 DO ipt = 1,nlay
481 bufly => elbuf_tab(ng)%BUFLY(ipt)
482 nptt = bufly%NPTT
483 IF (bufly%L_PLA > 0) THEN
484 func(6) = zero
485 DO it=1,nptt
486 lbuf => bufly%LBUF(1,1,it)
487 func(6) = func(6) + abs(lbuf%PLA(i))/nptt
488 ENDDO
489 wwa(13) = min(wwa(13),func(6))
490 wwa(14) = max(wwa(14),func(6))
491 ENDIF
492 ENDDO
493 ENDIF ! IF (NPG > 1) THEN
494 ELSE ! (NLAY = 1)
495 IF (npg > 1) THEN
496 bufly => elbuf_tab(ng)%BUFLY(1)
497 nptt = bufly%NPTT
498 DO it=1,nptt
499 i1 = (it-1)*nel
500 IF (bufly%L_PLA > 0) THEN
501 wwa(13) = min(wwa(13),abs(bufly%PLAPT(i1+i)))
502 wwa(14) = max(wwa(14),abs(bufly%PLAPT(i1+i)))
503 ENDIF
504 ENDDO
505 ELSE ! (NPG = 1)
506 bufly => elbuf_tab(ng)%BUFLY(1)
507 nptt = bufly%NPTT
508 DO it=1,nptt
509 lbuf => bufly%LBUF(1,1,it)
510 IF (bufly%L_PLA > 0) THEN
511 wwa(13) = min(wwa(13),abs(lbuf%PLA(i)))
512 wwa(14) = max(wwa(14),abs(lbuf%PLA(i)))
513 ENDIF
514 ENDDO
515 ENDIF
516 ENDIF ! IF (NLAY > 1) THEN
517c----
518c------------ Plastic strain per layer
519c----
520 IF (mte == 25) THEN
521 IF (ifailure == 0)THEN
522 wwa(30279) = fail(i)
523 wwa(30280) = 100*fail(i)/npt
524 wwa(30281) = fail1(i)
525 wwa(30282) = fail2(i)
526 wwa(30283) = fail3(i)
527 ENDIF ! IFAILURE == 0
528C
529 DO ipt=1,nlay
530 IF(ipt > 99) EXIT
531 bufly => elbuf_tab(ng)%BUFLY(ipt)
532 nptt = bufly%NPTT
533 val_ly_average = zero
534 DO ir=1,nptr
535 DO is=1,npts
536 val_ly_ip = zero
537 DO it=1,nptt
538 lbuf => bufly%LBUF(ir,is,it)
539 val_ly_ip = val_ly_ip + lbuf%PLA(i)/nptt
540 ENDDO
541 val_ly_average = val_ly_average + val_ly_ip/npg
542 ENDDO ! DO IS=1,NPTS
543 ENDDO ! DO IR=1,NPTR
544 wwa(30283 + ipt ) = val_ly_average
545 ENDDO ! DO IPT=1,NLAY
546 ENDIF ! MTE == 25
547 ENDIF ! GBUF%G_PLA > 0
548c----
549c------------ Non-local plastic strain and non-local plastic strain rate
550c----
551 IF (gbuf%G_PLANL > 0) THEN
552 bufly => elbuf_tab(ng)%BUFLY(1)
553 wwa(37855) = zero
554 nptt = bufly%NPTT
555 DO ir = 1,nptr
556 DO is = 1,npts
557 DO it = 1,nptt
558 wwa(37855) = wwa(37855) +
559 . bufly%LBUF(ir,is,it)%PLANL(i)/(nptr*npts*nptt)
560 ENDDO
561 ENDDO
562 ENDDO
563 ENDIF
564 IF (gbuf%G_EPSDNL > 0) THEN
565 bufly => elbuf_tab(ng)%BUFLY(1)
566 wwa(37856) = zero
567 nptt = bufly%NPTT
568 DO ir = 1,nptr
569 DO is = 1,npts
570 DO it = 1,nptt
571 wwa(37856) = wwa(37856) +
572 . bufly%LBUF(ir,is,it)%EPSDNL(i)/(nptr*npts*nptt)
573 ENDDO
574 ENDDO
575 ENDDO
576 ENDIF
577c----
578c------------ User variables
579c----
580 IF ((mte>=29.AND.mte<=31).OR.
581 . mte==35.OR.mte==36.OR.mte==43.OR.
582 . mte==44.OR.mte==45.OR.mte==48.OR.mte>=50) THEN
583c
584 nuvar = elbuf_tab(ng)%BUFLY(1)%NVAR_MAT
585 ALLOCATE (var(nuvar,max(1,mpt)))
586 var = zero
587C---
588 IF (mte == 58 .or. mte == 158) THEN
589 IF (nlay > 1) THEN
590 DO il=1,nlay
591 nuvar = elbuf_tab(ng)%BUFLY(il)%NVAR_MAT
592 nptt = elbuf_tab(ng)%BUFLY(il)%NPTT
593 DO ir=1,nptr
594 DO is=1,npts
595 k = nptr*(is-1) + ir
596 DO it=1,nptt
597 uvar=>elbuf_tab(ng)%BUFLY(il)%MAT(ir,is,it)%VAR
598 DO j = 1, nuvar
599 i1 = (j-1)*nel
600 IF (j==4 .OR. j==5) THEN ! convert to eng strain
601 var(j,il) = var(j,il) + (exp(uvar(i1+i))-one)/npg
602 ELSE
603 var(j,il) = var(j,il) + uvar(i1+i)/npg
604 ENDIF
605 wwa(6518 + (il-1)*60*4 + (k-1)*60 + j) =
606 . uvar(i1 + i)
607 ENDDO
608 ENDDO
609 ENDDO
610 ENDDO
611 ENDDO
612 ELSE ! NLAY = 1
613 nuvar = elbuf_tab(ng)%BUFLY(1)%NVAR_MAT
614 nptt = elbuf_tab(ng)%BUFLY(1)%NPTT
615 DO ipt=1,nptt
616 DO ir=1,nptr
617 DO is=1,npts
618 uvar=>elbuf_tab(ng)%BUFLY(1)%MAT(ir,is,ipt)%VAR
619 k = nptr*(is-1) + ir
620 DO j = 1, nuvar
621 i1 = (j-1)*nel
622 IF (j==4 .OR. j==5) THEN ! convert to eng strain
623 var(j,ipt) = var(j, ipt) + (exp(uvar(i1+i))-one)/npg
624 ELSE
625 var(j,ipt) = var(j, ipt) + uvar(i1 + i)/npg
626 ENDIF
627 wwa(6518 + (ipt-1)*60*4 + (k-1)*60 + j) =
628 . uvar(i1 + i)
629 ENDDO
630 ENDDO
631 ENDDO
632 ENDDO
633 ENDIF ! NLAY
634 ELSE ! (MTE /= 58)
635 IF (nlay > 1) THEN
636 DO il=1,nlay
637 nuvar = elbuf_tab(ng)%BUFLY(il)%NVAR_MAT
638 nptt = elbuf_tab(ng)%BUFLY(il)%NPTT
639 DO ir=1,nptr
640 DO is=1,npts
641 k = nptr*(is-1) + ir
642 DO it=1,nptt
643 uvar=>elbuf_tab(ng)%BUFLY(il)%MAT(ir,is,it)%VAR
644 DO j = 1, nuvar
645 i1 = (j-1)*nel
646 var(j,il) = var(j,il) + uvar(i1+i)/npg
647 wwa(6518 + (il-1)*60*4 + (k-1)*60 + j) =
648 . uvar(i1 + i)
649 ENDDO
650 ENDDO
651 ENDDO
652 ENDDO
653 ENDDO
654 ELSE ! NLAY = 1
655 nuvar = elbuf_tab(ng)%BUFLY(1)%NVAR_MAT
656 nptt = elbuf_tab(ng)%BUFLY(1)%NPTT
657 DO ipt=1,nptt
658 DO ir=1,nptr
659 DO is=1,npts
660 uvar=>elbuf_tab(ng)%BUFLY(1)%MAT(ir,is,ipt)%VAR
661 k = nptr*(is-1) + ir
662 DO j = 1, nuvar
663 i1 = (j-1)*nel
664 var(j,ipt) = var(j, ipt) + uvar(i1 + i)/npg
665 wwa(6518 + (ipt-1)*60*4 + (k-1)*60 + j) =
666 . uvar(i1 + i)
667 ENDDO
668 ENDDO
669 ENDDO
670 ENDDO
671 ENDIF ! NLAY
672 ENDIF ! MTE=58
673C---
674 nu = min(60,nuvar)
675 DO j = 1, nu
676 wwa(23+j)=var(j,iabs(mpt)/2 + 1) ! UVAR(NU) - membrane
677 DO ipt = 1, mpt
678 IF (j <= 20) THEN
679 IF (ipt <= 5) THEN
680 iuv = 83
681 iaa = 5
682 ELSEIF(ipt > 5)THEN
683 iuv = 678
684 iaa = 94
685 ENDIF
686 ELSE
687 iuv = 2558
688 iaa = 99
689 ENDIF
690 wwa(iuv + (j - 1)*iaa + ipt) = var(j, ipt) ! UVAR(NU, IPT)
691 ENDDO ! ipt =1,mpt
692 ENDDO ! J = 1, NU
693c
694 DEALLOCATE (var)
695 ENDIF ! MTE user
696c----
697c------------ Strain
698c----
699 IF (istrain /= 0) THEN
700 wwa(15)=gbuf%STRA(kk(1)+i)
701 wwa(16)=gbuf%STRA(kk(2)+i)
702 wwa(17)=gbuf%STRA(kk(3)+i)
703 wwa(18)=gbuf%STRA(kk(4)+i)
704 wwa(19)=gbuf%STRA(kk(5)+i)
705 wwa(20)=gbuf%STRA(kk(6)+i)
706 wwa(21)=gbuf%STRA(kk(7)+i)
707 wwa(22)=gbuf%STRA(kk(8)+i)
708 ENDIF
709C pinching data
710 IF(ihbe ==11.AND.npinch > 0) THEN
711 wwa(37848:37853) = zero
712 DO ipg=1,4
713 wwa(37847+1) = wwa(37847+1) + fourth*gbuf%EPGPINCHXZ(4*(i-1)+ipg)
714 wwa(37847+2) = wwa(37847+2) + fourth*gbuf%EPGPINCHYZ(4*(i-1)+ipg)
715 wwa(37847+3) = wwa(37847+3) + fourth*gbuf%EPGPINCHZZ(4*(i-1)+ipg)
716 wwa(37847+4) = wwa(37847+4) + fourth*gbuf%FORPGPINCH(4*(i-1)+ipg)
717 wwa(37847+5) = wwa(37847+5) + fourth*gbuf%MOMPGPINCH(8*(i-1)+2*(ipg-1)+1)
718 wwa(37847+6) = wwa(37847+6) + fourth*gbuf%MOMPGPINCH(8*(i-1)+2*ipg)
719 ENDDO
720 wwa(37847+7) = gbuf%THK(i)
721 ENDIF
722C end of pinching data
723c---------------------------------
724 DO l=iadv,iadv+nvar-1
725 k = ithbuf(l)
726 ijk = ijk+1
727 wa(ijk)=wwa(k)
728 ENDDO
729 ijk = ijk+1
730 wa(ijk) = ii
731c--------
732 ENDIF ! K==N
733 ENDDO ! I=1,NEL
734c--------
735 ENDIF ! MTE /= 13
736 ENDIF ! ITY == ITYP
737 ENDDO ! NG=1,NGROUP
738! -------------------------------
739 ENDIF
740 666 continue
741 ENDDO
742
743C-----------
744 RETURN
745 END
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21
pinchtype_mod
Definition pinchtype_mod.F:69
pinchtype_mod::npinch
integer npinch
Definition pinchtype_mod.F:73
stack_mod
Definition stack_mod.F:72
thcoq
subroutine thcoq(elbuf_tab, matparam_tab, nthgrp2, ithgrp, iparg, ithbuf, wa, ipm, igeo, ixc, ixtg, pm, rthbuf, thke, stack)
Definition thcoq.F:37