46 1 X ,IPARG ,IPM ,IGEO ,IXC ,
47 2 IXTG ,WA,WAP0 ,IPARTC, IPARTTG,
48 3 IPART_STATE,STAT_INDXC,STAT_INDXTG,THKE,SIZP0,
49 4 GEO ,STACK,DRAPE_SH4N,DRAPE_SH3N,DRAPEG)
57 use element_mod ,
only : nixc,nixtg
61#include "implicit_f.inc"
77 INTEGER IXC(NIXC,*),IXTG(NIXTG,*),
78 . IPARG(NPARG,*),IPM(NPROPMI,*),IGEO(NPROPGI,*),
79 . IPARTC(*), IPARTTG(*), IPART_STATE(*),
80 . stat_indxc(*), stat_indxtg(*)
82 . thke(*),x(3,*),geo(*)
83 TYPE (ELBUF_STRUCT_),
DIMENSION(NGROUP),
TARGET :: ELBUF_TAB
84 TYPE (STACK_PLY) :: STACK
85 TYPE (DRAPE_) :: DRAPE_SH4N(NUMELC_DRAPE), DRAPE_SH3N(NUMELTG_DRAPE)
86 TYPE (DRAPEG_) :: DRAPEG
87 double precision WA(*),WAP0(*)
91 INTEGER I,J,K,N,II,JJ,LEN,IOFF,IE,NG,NEL,NFT,LFT,NPT,
92 . LLT,ITY,MLW,IH,IHBE, ID, IPRT0, IPRT,IR,IS,IT,J1,J2,
93 . NPG,IPG,MPT,IPT,NPTR,NPTS,NPTT,NLAY,,ITHK,NF3,
94 . IGTYP,NPT_ALL,IL,KK(12),NF1,IREL,IBID0,MAT_1,PID_1,ILAY,IDRAPE,
96 INTEGER,
DIMENSION(:),
ALLOCATABLE :: PTWA
97 INTEGER,
DIMENSION(:),
ALLOCATABLE :: PTWA_P0
99 . THK, EM, EB, H1, H2, H3
101 . pg,mpg,qpg(2,4),thkq,
102 . sk(2),st(2),mk(2),mt(2),shk(2),sht(2),zz
103 CHARACTER*100 DELIMIT,LINE
104 TYPE(g_bufel_) ,
POINTER :: GBUF
105 TYPE(L_BUFEL_) ,
POINTER :: LBUF
106 TYPE(BUF_LAY_) ,
POINTER :: BUFLY
107 INTEGER LAYNPT_MAX,NLAY_MAX,ISUBSTACK,IPT_ALL,JDIR,L_DIRA,L_DIRB,IREP,
110 .
DIMENSION(:),
POINTER :: dir_a,dir_b
112 . qt(9,mvsiz),tens(6),zh,thkp ,thk0(mvsiz)
114 INTEGER,
DIMENSION(:) ,
ALLOCATABLE :: MATLY
115 my_real,
DIMENSION(:) ,
ALLOCATABLE :: thkly
116 my_real,
DIMENSION(:,:) ,
ALLOCATABLE :: posly,thk_ly
117 my_real,
ALLOCATABLE,
DIMENSION(:) ,
TARGET :: dira,dirb
119 parameter(pg = .577350269189626)
120 parameter(mpg=-.577350269189626)
121 DATA qpg/mpg,mpg,pg,mpg,pg,pg,mpg,pg/
123 ./
'#---1----|----2----|----3----|----4----|----5----|----6----|'/
125 ./
'----7----|----8----|----9----|----10---|'/
127 CALL my_alloc(ptwa,
max(stat_numelc ,stat_numeltg))
128 ALLOCATE(ptwa_p0(0:
max(1,stat_numelc_g,stat_numeltg_g)))
133 IF (stat_numelc==0)
GOTO 200
139 gbuf => elbuf_tab(ng)%GBUF
148 isubstack=iparg(71,ng)
150 nptr = elbuf_tab(ng)%NPTR
151 npts = elbuf_tab(ng)%NPTS
152 nptt = elbuf_tab(ng)%NPTT
153 nlay = elbuf_tab(ng)%NLAY
156 IF (ihbe == 23) npg=4
160 IF (ihbe>10.OR.igtyp==16.OR.ishfram ==0)
THEN
162 ELSEIF (ishfram ==1)
THEN
176 thk0(lft:llt) = gbuf%THK(lft:llt)
178 thk0(lft:llt) = thke(lft+nft:llt+nft)
182 IF(igtyp == 51 .OR. igtyp == 52)
THEN
183 DO ilay=1, elbuf_tab(ng
184 laynpt_max =
max(laynpt_max , elbuf_tab(ng)%BUFLY(ilay)%NPTT)
187 nlay_max =
max(nlay,npt, elbuf_tab(ng)%NLAY)
188 ALLOCATE(matly(mvsiz*nlay_max), thkly(mvsiz*nlay_max*laynpt_max),
189 . posly(mvsiz,nlay_max*laynpt_max),thk_ly(nel,nlay_max*laynpt_max))
194 numel_drape = numelc_drape
196 CALL layini(elbuf_tab(ng),lft ,llt ,geo ,igeo ,
197 . mat_1 ,pid_1 ,thkly ,matly ,posly ,
198 . igtyp ,ibid0 ,ibid0 ,nlay ,mpt ,
199 . isubstack,stack ,drape_sh4n ,nft ,thke ,
200 . nel ,thk_ly ,drapeg%INDX_SH4N,sedrape,numel_drape)
201 l_dira = elbuf_tab(ng)%BUFLY(1)%LY_DIRA
202 l_dirb = elbuf_tab(ng)%BUFLY(1)%LY_DIRB
203 ALLOCATE(dira(nlay*nel*l_dira))
204 ALLOCATE(dirb(nlay*nel*l_dirb))
207 IF (l_dira == 0)
THEN
209 ELSEIF (irep == 0)
THEN
210 IF(idrape > 0 .AND. (igtyp == 51 .OR. igtyp == 52))
THEN
212 j1 = 1+(j-1)*l_dira*nel
214 dira(j1:j2) = elbuf_tab(ng)%BUFLY(j)%LBUF_DIR(1)%DIRA(1:nel*l_dira)
218 j1 = 1+(j-1)*l_dira*nel
220 dira(j1:j2) = elbuf_tab(ng)%BUFLY(j)%DIRA(1:nel*l_dira)
224 dir_a => dira(1:nlay*nel*l_dira)
225 dir_b => dirb(1:nlay*nel*l_dirb)
226 CALL get_q4lsys(lft ,llt ,ixc(1,nf1),x ,gbuf%OFF,
227 . irel ,qt ,nlay ,irep ,nel ,
228 . dir_a ,dir_b,elbuf_tab(ng))
234 npt_all = npt_all + elbuf_tab(ng)%BUFLY(il)%NPTT
237 IF (iparg(6,ng) == 0) mpt=0
244 IF (ipart_state(iprt)==0) cycle
246 IF (mlw /= 0 .AND. mlw /= 13)
THEN
260 IF (mlw /= 0 .AND. mlw /= 13)
THEN
268 IF (mlw /= 0 .AND. mlw /= 13)
THEN
269 wa(jj) = gbuf%EINT(i)
274 IF (mlw /= 0 .AND. mlw /= 13)
THEN
275 wa(jj) = gbuf%EINT(i+llt)
280 IF (ihbe==11 .or. ihbe==23 .or. mlw == 0 .or. mlw == 13)
THEN
289 wa(jj) = gbuf%HOURG(kk(1)+i)
291 wa(jj) = gbuf%HOURG(kk(2)+i)
293 wa(jj) = gbuf%HOURG(kk(3)+i)
298 IF (mlw == 0 .or. mlw == 13)
THEN
305 ELSEIF (npg == 1)
THEN
306 tens(1:5) = gbuf%FOR(kk(1:5)+i)
313 tens(1:3) = gbuf%MOM(kk(1:3)+i)
321 IF (gbuf%G_PLA > 0)
THEN
329 lbuf => elbuf_tab(ng)%BUFLY(1)%LBUF(ir,is,1)
330 ipg = nptr*(is-1) + ir
333 tens(1:5) = gbuf%FORPG(k+kk(1:5)+i)
341 IF (gbuf%G_PLA > 0)
THEN
348 tens(1:3) = gbuf%MOMPG(k+kk(1:3)+i)
358 ELSEIF (mlw == 0 .or. mlw == 13)
THEN
370 bufly => elbuf_tab(ng)%BUFLY(il)
373 jdir = 1 + (il-1)*nel*2
379 lbuf => bufly%LBUF(ir,is,it)
380 tens(1:5) = lbuf%SIG(kk(1:5)+i)
381 CALL orth2loc(tens,dir_a,dir_b,ii,ilaw,igtyp,nel)
388 IF (bufly%L_PLA > 0)
THEN
394 wa(jj) = posly(i,ipt)*two
398 ipt_all = ipt_all + nptt
404 IF (mlw==0 .or. mlw==13)
THEN
433 st(1) = gbuf%HOURG(kk(1)+i)
434 st(2) =-gbuf%HOURG(kk(2)+i)
435 mt(1) = gbuf%HOURG(kk(3)+i)
436 mt(2) =-gbuf%HOURG(kk(4)+i)
437 sk(1) =-gbuf%HOURG(kk(7)+i)
438 sk(2) = gbuf%HOURG(kk(8)+i)
439 mk(1) =-gbuf%HOURG(kk(9)+i)
440 mk(2) = gbuf%HOURG(kk(10)+i)
441 sht(1)= gbuf%HOURG(kk(5)+i)
442 sht(2)=-gbuf%HOURG(kk(6)+i)
444 shk(2)= gbuf%HOURG(kk(12)+i)
447 IF (mpt == 0 .and. mlw /= 0 .and. mlw /= 13)
THEN
449 tens(1:2) = gbuf%FOR(kk(1:2)+i)
450 . + st(1:2)*qpg(2,ipg)+sk(1:2)*qpg(1,ipg)
451 tens(3) = gbuf%FOR(kk(3)+i)
452 tens(4) = gbuf%FOR(kk(4)+i)
453 . + sht(2)*qpg(2,ipg)+shk(2)*qpg(1,ipg)
454 tens(5) = gbuf%FOR(kk(5)+i)
455 . + sht(1)*qpg(2,ipg)+shk(1)*qpg(1,ipg)
461 tens(1:2) = gbuf%MOM(kk(1:2)+i)
462 . + mt(1:2)*qpg(2,ipg)+mk(1:2)*qpg(1,ipg)
463 tens(3) = gbuf%MOM(kk(3)+i)
471 IF (gbuf%G_PLA > 0)
THEN
477 ELSEIF (mlw /= 0 .and. mlw /= 13)
THEN
480 bufly => elbuf_tab(ng)%BUFLY(il)
483 jdir = 1 + (il-1)*nel*2
487 lbuf => bufly%LBUF(1,1,it)
489 zz = posly(i,ipt)*thkq
491 tens(1:2) = lbuf%SIG(kk(1:2)+i)
492 . + (st(1:2)+zz*mt(1:2))*qpg(2,ipg)
493 . + (sk(1:2)+zz*mk(1:2))*qpg(1,ipg)
494 tens(3) = lbuf%SIG(kk(3)+i)
495 tens(4) = lbuf%SIG(kk(4)+i)
496 . + sht(2)*qpg(2,ipg)+shk(2)*qpg(1,ipg)
497 tens(5) = lbuf%SIG(kk(5)+i)
498 . + sht(1)*qpg(2,ipg)+shk(1)*qpg(1,ipg)
499 CALL orth2loc(tens,dir_a,dir_b,ii,ilaw,igtyp,nel
515 ipt_all = ipt_all + nptt
525 IF(
ALLOCATED(dirb))
DEALLOCATE(dirb)
526 IF(
ALLOCATED(dira))
DEALLOCATE(dira)
527 DEALLOCATE(matly, thkly, posly, thk_ly)
551 IF (ispmd == 0.AND.len
THEN
559 ioff = nint(wap0(j + 1))
561 iprt = nint(wap0(j + 2))
562 IF (iprt /= iprt0)
THEN
563 IF (izipstrs == 0)
THEN
564 WRITE(iugeo,
'(A)') delimit
565 WRITE(iugeo,
'(A)')
'/INISHE/STRS_F/GLOB'
567 .
'#------------------------ REPEAT --------------------------'
569 .
'# SHELLID NPT NPG THK'
570 WRITE(iugeo,
'(A)')
'# EM, EB, H1, H2, H3'
571 WRITE(iugeo,
'(A/A/A/A/A)')
572 .
'# IF(NPT == 0), REPEAT I=1,NPG :',
576 .
'# M12,M23,M31,EPSP '
577 WRITE(iugeo,
'(A/A/A)')
578 .
'# IF(NPT /= 0) REPEAT K=1,NPT : REPEAT I=1,NPG :',
580 .
'# S12, S23, S31, EPSP, T '
582 .
'#---------------------- END REPEAT ------------------------'
583 WRITE(iugeo,
'(A)') delimit
585 WRITE(line,
'(A)') delimit
587 WRITE(line,
'(A)')
'/INISHE/STRS_F/GLOB'
590 .
'#------------------------ REPEAT --------------------------'
593 .
'# SHELLID NPT NPG THK'
595 WRITE(line,
'(A)')
'# EM, EB, H1, H2, H3'
597 WRITE(line,
'(A)')
'# IF(NPT == 0), REPEAT I=1,NPG :'
599 WRITE(line,
'(A)')
'# N1, N2, N3 '
601 WRITE(line,
'(A)')
'# N12, N23, N31'
603 WRITE(line,
'(A)')
'# M1, M2, M3 '
605 WRITE(line,
'(A)')
'# M12, M23, M31, EPSP'
608 .
'# IF(NPT /= 0) REPEAT K=1,NPT : REPEAT I=1,NPG :'
610 WRITE(line,
'(A)')
'# S1, S2, S3'
612 WRITE(line,
'(A)')
'# S12,S23,S31, EPSP, T '
615 .
'#---------------------- END REPEAT ------------------------'
617 WRITE(line,
'(A)') delimit
624 npt = nint(wap0(j + 4))
625 npg = nint(wap0(j + 5))
633 IF (izipstrs == 0)
THEN
634 WRITE(iugeo,
'(3I10,1PE20.13)')id,npt,npg,thk
635 WRITE(iugeo,
'(1P5E20.13)')em,eb,h1,h2,h3
637 WRITE(line,
'(3I10,1PE20.13)')id,npt,npg,thk
639 WRITE(line,
'(1P5E20.13)')em,eb,h1,h2,h3
644 IF (izipstrs == 0)
THEN
645 WRITE(iugeo,
'(1P3E20.13)')(wap0(j + k),k=1,9)
646 WRITE(iugeo,
'(1P4E20.13)')(wap0(j + k),k=10,13
656 IF (izipstrs == 0)
THEN
657 WRITE(iugeo,
'(1P3E20.13)')(wap0(j + k),k=1,3)
658 WRITE(iugeo,
'(1P5E20.13)')(wap0(j + k),k
674 IF (stat_numeltg==0)
GOTO 300
680 gbuf => elbuf_tab(ng)%GBUF
689 isubstack=iparg(71,ng)
690 nptr = elbuf_tab(ng)%NPTR
691 npts = elbuf_tab(ng)%NPTS
692 nptt = elbuf_tab(ng)%NPTT
693 nlay = elbuf_tab(ng)%NLAY
710 pid_1 = ixtg(nixtg-1,nf1)
712 thk0(lft:llt) = gbuf%THK(lft:llt)
715 thk0(lft:llt) = thke(lft+nf3:llt+nf3)
719 IF(igtyp == 51 .OR. igtyp == 52)
THEN
720 DO ilay=1, elbuf_tab(ng)%NLAY
721 laynpt_max =
max(laynpt_max , elbuf_tab(ng)%BUFLY(ilay)%NPTT)
724 nlay_max =
max(nlay,npt, elbuf_tab(ng)%NLAY)
725 ALLOCATE(matly(mvsiz*nlay_max), thkly(mvsiz*nlay_max*laynpt_max),
726 . posly(mvsiz,nlay_max*laynpt_max),thk_ly(nel,nlay_max*laynpt_max))
731 numel_drape = numeltg_drape
733 CALL layini(elbuf_tab(ng),lft ,llt ,geo ,igeo
734 . mat_1 ,pid_1 ,thkly ,matly
735 . igtyp ,ibid0 ,ibid0 ,nlay ,mpt ,
736 . isubstack,stack ,drape_sh3n ,nft ,thke ,
737 . nel ,thk_ly ,drapeg%INDX_SH3N, sedrape,numel_drape
739 l_dira = elbuf_tab(ng)%BUFLY(1)%LY_DIRA
740 l_dirb = elbuf_tab(ng)%BUFLY(1)%LY_DIRB
741 ALLOCATE(dira(nlay*nel*l_dira))
742 ALLOCATE(dirb(nlay*nel*l_dirb))
745 IF (l_dira == 0)
THEN
747 ELSEIF (irep == 0)
THEN
749 j1 = 1+(j-1)*l_dira*nel
751 dira(j1:j2) = elbuf_tab(ng)%BUFLY(j)%DIRA(1:nel*l_dira)
754 dir_a => dira(1:nlay*nel*l_dira)
755 dir_b => dirb(1:nlay*nel*l_dirb)
756 CALL get_t3lsys(lft ,llt ,ixtg(1,nf1),x ,gbuf%OFF,
757 . irel ,qt ,nlay ,irep ,nel ,
758 . dir_a ,dir_b,elbuf_tab(ng))
764 npt_all = npt_all + elbuf_tab(ng)%BUFLY(il)%NPTT
767 IF (iparg(6,ng) == 0) mpt=0
774 IF (ipart_state(iprt) == 0) cycle
776 IF (mlw /= 0 .AND. mlw /= 13)
THEN
784 wa(jj) = ixtg(nixtg,n)
790 IF (mlw /= 0 .AND. mlw /= 13)
THEN
798 IF (mlw /= 0 .AND. mlw /= 13)
THEN
799 wa(jj) = gbuf%EINT(i)
804 IF (mlw /= 0 .AND. mlw /= 13)
THEN
805 wa(jj) = gbuf%EINT(i+llt)
817 IF (mlw == 0 .or. mlw == 13)
THEN
824 ELSEIF (npg == 1)
THEN
825 tens(1:5) = gbuf%FOR(kk(1:5)+i)
832 tens(1:3) = gbuf%MOM(kk(1:3)+i)
840 IF (gbuf%G_PLA > 0)
THEN
847 lbuf => elbuf_tab(ng)%BUFLY(1)%LBUF(ipg,1,1)
850 tens(1:5) = gbuf%FORPG(k+kk(1:5)+i)
858 IF (gbuf%G_PLA > 0)
THEN
865 tens(1:3) = gbuf%MOMPG(k+kk(1:3)+i)
874 IF (mlw == 0 .or. mlw == 13)
THEN
886 bufly => elbuf_tab(ng)%BUFLY(il)
889 jdir = 1 + (il-1)*nel*2
894 lbuf => bufly%LBUF(ipg,1,it)
895 tens(1:5) = lbuf%SIG(kk(1:5)+i)
896 CALL orth2loc(tens,dir_a,dir_b,ii,ilaw,igtyp,nel)
903 IF (bufly%L_PLA > 0)
THEN
909 wa(jj) = posly(i,ipt)*two
912 ipt_all = ipt_all + nptt
922 IF(
ALLOCATED(dirb))
DEALLOCATE(dirb)
923 IF(
ALLOCATED(dira))
DEALLOCATE(dira)
924 DEALLOCATE(matly, thkly, posly, thk_ly)
946 IF (ispmd == 0.AND.len > 0)
THEN
948 DO n=1,stat_numeltg_g
954 ioff = nint(wap0(j + 1))
956 iprt = nint(wap0(j + 2))
957 IF (iprt /= iprt0)
THEN
958 IF (izipstrs == 0)
THEN
959 WRITE(iugeo,
'(A)') delimit
960 WRITE(iugeo,
'(A)')
'/INISH3/STRS_F/GLOB'
962 .
'#------------------------ REPEAT --------------------------'
964 .
'# SH3NID NPT NPG THK'
966 .
'# EM, EB, H1, H2, H3'
967 WRITE(iugeo,
'(A/A/A/A/A)')
968 .
'# IF(NPT == 0), REPEAT I=1,NPG :',
972 .
'# M12,M23,M31,EPSP '
973 WRITE(iugeo,
'(A/A/A)')
974 .
'# IF(NPT /= 0) REPEAT K=1,NPT : REPEAT I=1,NPG :',
976 .
'# S12,S23,S31, EPSP, T '
978 .
'#---------------------- END REPEAT ------------------------'
979 WRITE(iugeo,
'(A)') delimit
981 WRITE(line,
'(A)') delimit
983 WRITE(line,
'(A)')
'/INISH3/STRS_F/GLOB'
986 .
'#------------------------ REPEAT --------------------------'
989 .
'# SH3NID NPT NPG THK'
992 .
'# EM, EB, H1, H2, H3'
995 .
'# IF(NPT == 0), REPEAT I=1,NPG :'
997 WRITE(line,
'(A)')
'# N1, N2, N3'
999 WRITE(line,
'(A)')
'# N12, N23, N31'
1001 WRITE(line,
'(A)')
'# M1, M2, M3 '
1003 WRITE(line,
'(A)')
'# M12, M23, M31,EPSP '
1006 .
'# IF(NPT /= 0) REPEAT K=1,NPT : REPEAT I=1,NPG :'
1008 WRITE(line,
'(A)')
'# S1, S2, S3 '
1010 WRITE(line,
'(A)')
'# S12, S23, S31, EPSP, T '
1013 .
'#---------------------- END REPEAT ------------------------'
1015 WRITE(line,
'(A)') delimit
1020 id = nint(wap0(j + 3))
1021 npt = nint(wap0(j + 4))
1022 npg = nint(wap0(j + 5))
1030 IF (izipstrs == 0)
THEN
1031 WRITE(iugeo,
'(3I10,1PE20.13)')id,npt,npg,thk
1032 WRITE(iugeo,
'(1P5E20.13)')em,eb,h1,h2,h3
1034 WRITE(line,
'(3I10,1PE20.13)')id,npt,npg,thk
1036 WRITE(line,
'(1P5E20.13)')em,eb,h1,h2,h3
1041 IF (izipstrs == 0)
THEN
1042 WRITE(iugeo,
'(1P3E20.13)')(wap0(j + k),k=1,9)
1043 WRITE(iugeo,
'(1P4E20.13)')(wap0(j + k),k=10,13)
1053 IF (izipstrs == 0)
THEN
1054 WRITE(iugeo,
'(1P3E20.13)')(wap0(j + k),k=1,3)
1055 WRITE(iugeo,
'(1P5E20.13)')(wap0(j + k),k=4,8)