OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
h3d_shell_tensor.F File Reference
#include "implicit_f.inc"
#include "mvsiz_p.inc"
#include "com01_c.inc"
#include "com04_c.inc"
#include "param_c.inc"
#include "tabsiz_c.inc"

Go to the source code of this file.

Functions/Subroutines

subroutine h3d_shell_tensor (elbuf_tab, shell_tensor, iparg, itens, invert, nelcut, el2fa, nbf, tens, epsdot, iadp, nbf_l, nbpart, iadg, x, ixc, igeo, ixtg, ipm, stack, id_elem, ity_elem, info1, info2, is_written_shell, ipartc, iparttg, layer_input, ipt_input, ply_input, gauss_input, iuvar_input, h3d_part, keyword, d, id, bufmat, mat_param, geo, drape_sh4n, drape_sh3n, drapeg)
subroutine sh4_tstrain (xn, yn, zn, dx, dy, dz, strain, nel)
subroutine c4sysg2l (xn, yn, zn, xl2, yl2, xl3, yl3, xl4, yl4, zl1, area, nel)
subroutine clsys3 (rx, ry, rz, sx, sy, sz, vq, det, nel)
subroutine sh3_tstrain (xn, yn, zn, dx, dy, dz, strain, nel, ihbe)
subroutine c3sysg2l (xn, yn, zn, xl2, yl2, xl3, yl3, area, nel)
subroutine u_from_f2 (f, strain, nel)

Function/Subroutine Documentation

◆ c3sysg2l()

subroutine c3sysg2l ( xn,
yn,
zn,
xl2,
yl2,
xl3,
yl3,
area,
integer nel )

Definition at line 2529 of file h3d_shell_tensor.F.

2530C-----------------------------------------------
2531C I m p l i c i t T y p e s
2532C-----------------------------------------------
2533#include "implicit_f.inc"
2534C---------+---------+---+---+--------------------------------------------
2535C VAR | SIZE |TYP| RW| DEFINITION
2536C---------+---------+---+---+--------------------------------------------
2537C NEL | 1 | I | R | SIZE OF THE ELEMENT GROUP NEL
2538C XN | 3*NEL | R | R | X-coordinate ARRAY (3n tria)
2539C YN | 3*NEL | R | R | Y-coordinate ARRAY (3n tria)
2540C ZN | 3*NEL | R | R | Z-coordinate ARRAY (3n tria)
2541C XL2 | NEL | R | W | Local X-coordinate of N2 (relative to N1)
2542C YL2 | NEL | R | W | Local Y-coordinate of N2 (relative to N1)
2543C XL3 | NEL | R | W | Local X-coordinate of N3 (relative to N1)
2544C YL3 | NEL | R | W | Local Y-coordinate of N3 (relative to N1)
2545C AREA | NEL | R | W | AREA of tria
2546C---------+---------+---+---+--------------------------------------------
2547C-----------------------------------------------
2548C D u m m y A r g u m e n t s
2549C-----------------------------------------------
2550 INTEGER NEL
2551 my_real
2552 . xn(3,*) , yn(3,*) , zn(3,*),
2553 . xl2(nel),xl3(nel),yl2(nel),yl3(nel),area(nel)
2554C-----------------------------------------------
2555C L o c a l V a r i a b l e s
2556C-----------------------------------------------
2557 INTEGER I
2558 my_real
2559 . rx(nel),ry(nel),rz(nel),sx(nel),sy(nel),sz(nel),
2560 . vq(3,3,nel), deta1(nel),xx,yy,zz
2561 DO i=1,nel
2562 rx(i)=xn(2,i)-xn(1,i)
2563 ry(i)=yn(2,i)-yn(1,i)
2564 rz(i)=zn(2,i)-zn(1,i)
2565 sx(i)=xn(3,i)-xn(1,i)
2566 sy(i)=yn(3,i)-yn(1,i)
2567 sz(i)=zn(3,i)-zn(1,i)
2568 ENDDO
2569 CALL clsys3(rx, ry, rz, sx, sy, sz,
2570 . vq, deta1,nel)
2571C------ Global -> Local Coordinate
2572 DO i=1,nel
2573 xx=rx(i)
2574 yy=ry(i)
2575 zz=rz(i)
2576 xl2(i)=vq(1,1,i)*xx+vq(2,1,i)*yy+vq(3,1,i)*zz
2577 yl2(i)=vq(1,2,i)*xx+vq(2,2,i)*yy+vq(3,2,i)*zz
2578C
2579 xx=sx(i)
2580 yy=sy(i)
2581 zz=sz(i)
2582 xl3(i)=vq(1,1,i)*xx+vq(2,1,i)*yy+vq(3,1,i)*zz
2583 yl3(i)=vq(1,2,i)*xx+vq(2,2,i)*yy+vq(3,2,i)*zz
2584 area(i)=0.5*deta1(i)
2585 ENDDO
2586c-----------
2587 RETURN
#define my_real
Definition cppsort.cpp:32
subroutine clsys3(rx, ry, rz, sx, sy, sz, vq, det, nel)
subroutine area(d1, x, x2, y, y2, eint, stif0)

◆ c4sysg2l()

subroutine c4sysg2l ( xn,
yn,
zn,
xl2,
yl2,
xl3,
yl3,
xl4,
yl4,
zl1,
area,
integer nel )

Definition at line 2254 of file h3d_shell_tensor.F.

2255C-----------------------------------------------
2256C I m p l i c i t T y p e s
2257C-----------------------------------------------
2258#include "implicit_f.inc"
2259C---------+---------+---+---+--------------------------------------------
2260C VAR | SIZE |TYP| RW| DEFINITION
2261C---------+---------+---+---+--------------------------------------------
2262C NEL | 1 | I | R | SIZE OF THE ELEMENT GROUP NEL
2263C XN | 4*NEL | R | R | X-coordinate ARRAY (4n quad)
2264C YN | 4*NEL | R | R | Y-coordinate ARRAY (4n quad)
2265C ZN | 4*NEL | R | R | Z-coordinate ARRAY (4n quad)
2266C XL2 | NEL | R | W | Local X-coordinate of N2 (relative to N1)
2267C YL2 | NEL | R | W | Local Y-coordinate of N2 (relative to N1)
2268C ZL1 | NEL | R | W | Local Z-coordinate of N1 (Z1 -Z1 Z1 -Z1)
2269C XL3 | NEL | R | W | Local X-coordinate of N3 (relative to N1)
2270C YL3 | NEL | R | W | Local Y-coordinate of N3 (relative to N1)
2271C XL4 | NEL | R | W | Local X-coordinate of N4 (relative to N1)
2272C YL4 | NEL | R | W | Local Y-coordinate of N4 (relative to N1)
2273C AREA | NEL | R | W | AREA of quad
2274C---------+---------+---+---+--------------------------------------------
2275C-----------------------------------------------
2276C D u m m y A r g u m e n t s
2277C-----------------------------------------------
2278 INTEGER NEL
2279 my_real
2280 . xn(4,*) , yn(4,*) , zn(4,*),
2281 . xl2(nel),xl3(nel),xl4(nel),yl2(nel),
2282 . yl3(nel),yl4(nel),zl1(nel),area(nel)
2283C-----------------------------------------------
2284C L o c a l V a r i a b l e s
2285C-----------------------------------------------
2286 INTEGER I
2287 my_real
2288 . rx(nel),ry(nel),rz(nel),sx(nel),sy(nel),sz(nel),
2289 . vq(3,3,nel), lxyz0(3),deta1(nel),xx,yy,zz
2290C-----------------------------------------------
2291 DO i=1,nel
2292 rx(i)=xn(2,i)+xn(3,i)-xn(1,i)-xn(4,i)
2293 ry(i)=yn(2,i)+yn(3,i)-yn(1,i)-yn(4,i)
2294 rz(i)=zn(2,i)+zn(3,i)-zn(1,i)-zn(4,i)
2295 sx(i)=xn(3,i)+xn(4,i)-xn(1,i)-xn(2,i)
2296 sy(i)=yn(3,i)+yn(4,i)-yn(1,i)-yn(2,i)
2297 sz(i)=zn(3,i)+zn(4,i)-zn(1,i)-zn(2,i)
2298 ENDDO
2299C------Local elem. base:
2300 CALL clsys3(rx, ry, rz, sx, sy, sz,
2301 . vq, deta1,nel)
2302C------ Global -> Local Coordinate FOURTH=0.25 ;
2303 DO i=1,nel
2304 lxyz0(1)=fourth*(xn(1,i)+xn(2,i)+xn(3,i)+xn(4,i))
2305 lxyz0(2)=fourth*(yn(1,i)+yn(2,i)+yn(3,i)+yn(4,i))
2306 lxyz0(3)=fourth*(zn(1,i)+zn(2,i)+zn(3,i)+zn(4,i))
2307 xx=xn(2,i)-xn(1,i)
2308 yy=yn(2,i)-yn(1,i)
2309 zz=zn(2,i)-zn(1,i)
2310 xl2(i)=vq(1,1,i)*xx+vq(2,1,i)*yy+vq(3,1,i)*zz
2311 yl2(i)=vq(1,2,i)*xx+vq(2,2,i)*yy+vq(3,2,i)*zz
2312 xx=xn(2,i)-lxyz0(1)
2313 yy=yn(2,i)-lxyz0(2)
2314 zz=zn(2,i)-lxyz0(3)
2315 zl1(i)=vq(1,3,i)*xx+vq(2,3,i)*yy+vq(3,3,i)*zz
2316C
2317 xx=xn(3,i)-xn(1,i)
2318 yy=yn(3,i)-yn(1,i)
2319 zz=zn(3,i)-zn(1,i)
2320 xl3(i)=vq(1,1,i)*xx+vq(2,1,i)*yy+vq(3,1,i)*zz
2321 yl3(i)=vq(1,2,i)*xx+vq(2,2,i)*yy+vq(3,2,i)*zz
2322C
2323 xx=xn(4,i)-xn(1,i)
2324 yy=yn(4,i)-yn(1,i)
2325 zz=zn(4,i)-zn(1,i)
2326 xl4(i)=vq(1,1,i)*xx+vq(2,1,i)*yy+vq(3,1,i)*zz
2327 yl4(i)=vq(1,2,i)*xx+vq(2,2,i)*yy+vq(3,2,i)*zz
2328 area(i)=fourth*deta1(i)
2329 ENDDO
2330c-----------
2331 RETURN

◆ clsys3()

subroutine clsys3 ( rx,
ry,
rz,
sx,
sy,
sz,
vq,
det,
integer nel )

Definition at line 2341 of file h3d_shell_tensor.F.

2342C-----------------------------------------------
2343C I m p l i c i t T y p e s
2344C-----------------------------------------------
2345#include "implicit_f.inc"
2346C-----------------------------------------------
2347C D u m m y A r g u m e n t s
2348C-----------------------------------------------
2349 INTEGER NEL
2350 my_real
2351 . rx(*) , ry(*) , rz(*),
2352 . sx(*) , sy(*) , sz(*),
2353 . det(*),vq(3,3,*)
2354C---------+---------+---+---+--------------------------------------------
2355C VAR | SIZE |TYP| RW| DEFINITION
2356C---------+---------+---+---+--------------------------------------------
2357C NEL | 1 | I | R | SIZE OF THE ELEMENT GROUP NEL
2358C RX | NEL | R | R | X-of covariant vecter g1
2359C RY | NEL | R | R | Y-of covariant vecter g1
2360C RZ | NEL | R | R | Z-of covariant vecter g1
2361C SX | NEL | R | R | X-of covariant vecter g2
2362C SY | NEL | R | R | Y-of covariant vecter g2
2363C SZ | NEL | R | R | Z-of covariant vecter g2
2364C VQ |3*3*NEL | R | W | Local elem sys bases
2365C DET | NEL | R | W | det of g1 ^ g2
2366C---------+---------+---+---+--------------------------------------------
2367C-----------------------------------------------
2368C L o c a l V a r i a b l e s
2369C-----------------------------------------------
2370 INTEGER I
2371C REAL
2372 my_real
2373 . e1x(nel), e1y(nel), e1z(nel),
2374 . e2x(nel), e2y(nel), e2z(nel),
2375 . e3x(nel), e3y(nel), e3z(nel),
2376 . c1,c2,cc,c1c1,c2c2,c1_1,c2_1
2377C~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
2378 DO i=1,nel
2379C---------E3------------
2380 e3x(i) = ry(i) * sz(i) - rz(i) * sy(i)
2381 e3y(i) = rz(i) * sx(i) - rx(i) * sz(i)
2382 e3z(i) = rx(i) * sy(i) - ry(i) * sx(i)
2383 det(i) = sqrt(e3x(i)*e3x(i) + e3y(i)*e3y(i) + e3z(i)*e3z(i))
2384C ----- EM20=1.0E-20
2385 det(i) = max(em20,det(i))
2386 e3x(i) = e3x(i) / det(i)
2387 e3y(i) = e3y(i) / det(i)
2388 e3z(i) = e3z(i) / det(i)
2389 ENDDO
2390C~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
2391 DO i=1,nel
2392 c1c1 = rx(i)*rx(i) + ry(i)*ry(i) + rz(i)*rz(i)
2393 c2c2 = sx(i)*sx(i) + sy(i)*sy(i) + sz(i)*sz(i)
2394C ----- ZERO=0., ONE=1.0
2395 IF(c1c1 /= zero) THEN
2396 c2_1 = sqrt(c2c2/max(em20,c1c1))
2397 c1_1 = one
2398 ELSEIF(c2c2 /= zero)THEN
2399 c2_1 = one
2400 c1_1 = sqrt(c1c1/max(em20,c2c2))
2401 END IF
2402 e1x(i) = rx(i)*c2_1+(sy(i)*e3z(i)-sz(i)*e3y(i))*c1_1
2403 e1y(i) = ry(i)*c2_1+(sz(i)*e3x(i)-sx(i)*e3z(i))*c1_1
2404 e1z(i) = rz(i)*c2_1+(sx(i)*e3y(i)-sy(i)*e3x(i))*c1_1
2405 ENDDO
2406C~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
2407 DO i=1,nel
2408 c1 = sqrt(e1x(i)*e1x(i) + e1y(i)*e1y(i) + e1z(i)*e1z(i))
2409 IF ( c1 /= zero) c1 = one / max(em20,c1)
2410 e1x(i) = e1x(i)*c1
2411 e1y(i) = e1y(i)*c1
2412 e1z(i) = e1z(i)*c1
2413 e2x(i) = e3y(i) * e1z(i) - e3z(i) * e1y(i)
2414 e2y(i) = e3z(i) * e1x(i) - e3x(i) * e1z(i)
2415 e2z(i) = e3x(i) * e1y(i) - e3y(i) * e1x(i)
2416 ENDDO
2417 DO i=1,nel
2418 vq(1,1,i)=e1x(i)
2419 vq(2,1,i)=e1y(i)
2420 vq(3,1,i)=e1z(i)
2421 vq(1,2,i)=e2x(i)
2422 vq(2,2,i)=e2y(i)
2423 vq(3,2,i)=e2z(i)
2424 vq(1,3,i)=e3x(i)
2425 vq(2,3,i)=e3y(i)
2426 vq(3,3,i)=e3z(i)
2427 ENDDO
2428c-----------
2429 RETURN
#define max(a, b)
Definition macros.h:21

◆ h3d_shell_tensor()

subroutine h3d_shell_tensor ( type (elbuf_struct_), dimension(ngroup), target elbuf_tab,
shell_tensor,
integer, dimension(nparg,*) iparg,
integer itens,
integer, dimension(*) invert,
integer nelcut,
integer, dimension(*) el2fa,
integer nbf,
tens,
epsdot,
integer, dimension(*) iadp,
integer nbf_l,
integer nbpart,
integer, dimension(nspmd,*) iadg,
x,
integer, dimension(nixc,*) ixc,
integer, dimension(npropgi,*) igeo,
integer, dimension(nixtg,*) ixtg,
integer, dimension(npropmi,*) ipm,
type (stack_ply) stack,
integer, dimension(*) id_elem,
integer, dimension(*) ity_elem,
integer info1,
integer info2,
integer, dimension(*) is_written_shell,
integer, dimension(*) ipartc,
integer, dimension(*) iparttg,
integer layer_input,
integer ipt_input,
integer ply_input,
integer gauss_input,
integer iuvar_input,
integer, dimension(*) h3d_part,
character(len=ncharline100) keyword,
d,
integer id,
dimension(sbufmat), target bufmat,
type (matparam_struct_), dimension(nummat), intent(in) mat_param,
dimension(npropg,numgeo), intent(in) geo,
type (drape_), dimension(numelc_drape), intent(in) drape_sh4n,
type (drape_), dimension(numeltg_drape), intent(in) drape_sh3n,
type (drapeg_), intent(in) drapeg )

Definition at line 44 of file h3d_shell_tensor.F.

51C-----------------------------------------------
52C M o d u l e s
53C-----------------------------------------------
54 USE elbufdef_mod
55 USE stack_mod
56 USE matparam_def_mod
58 USE drape_mod
59C-----------------------------------------------
60C I m p l i c i t T y p e s
61C-----------------------------------------------
62#include "implicit_f.inc"
63C-----------------------------------------------
64C C o m m o n B l o c k s
65C-----------------------------------------------
66#include "mvsiz_p.inc"
67C-----------------------------------------------
68#include "com01_c.inc"
69#include "com04_c.inc"
70#include "param_c.inc"
71#include "tabsiz_c.inc"
72C-----------------------------------------------
73C D u m m y A r g u m e n t s
74C-----------------------------------------------
75 INTEGER IPARG(NPARG,*),ITENS,INVERT(*),IUVAR_INPUT,
76 . EL2FA(*),IXC(NIXC,*), IGEO(NPROPGI,*),
77 . NELCUT,NBF,IADP(*),NBF_L,NBPART,IADG(NSPMD,*),
78 . IXTG(NIXTG,*),IPM(NPROPMI,*),ID_ELEM(*),ITY_ELEM(*),
79 . INFO1,INFO2,IS_WRITTEN_SHELL(*),IPARTC(*),IPARTTG(*),H3D_PART(*),
80 . LAYER_INPUT ,IPT_INPUT,GAUSS_INPUT,PLY_INPUT,ID
81 my_real, INTENT(IN),TARGET :: bufmat(sbufmat)
83 . tens(3,*),epsdot(6,*),x(3,*),shell_tensor(3,*),d(3,*)
84 my_real, INTENT(IN) :: geo(npropg,numgeo)
85 TYPE (ELBUF_STRUCT_), DIMENSION(NGROUP), TARGET :: ELBUF_TAB
86 TYPE (STACK_PLY) :: STACK
87 TYPE (MATPARAM_STRUCT_) ,DIMENSION(NUMMAT) ,INTENT(IN) :: MAT_PARAM
88 CHARACTER(LEN=NCHARLINE100):: KEYWORD
89 TYPE (DRAPE_) , INTENT(IN) :: DRAPE_SH4N(NUMELC_DRAPE), DRAPE_SH3N(NUMELTG_DRAPE)
90 TYPE (DRAPEG_), INTENT(IN) :: DRAPEG
91C-----------------------------------------------
92C L o c a l V a r i a b l e s
93C-----------------------------------------------
94 my_real :: a1,a2,a3,thk,chard,factor
95 my_real :: value(5)
96 INTEGER I,J,K,N,NG,NEL,NFT,ITY,NPT,MPT,IPT,NBFUNCT,NCHARD,MLW,
97 . ILAY,IR,IS,IT,NPTR,NPTS,NPTT,NLAY,NPG,IPLY,IDRAPE,
98 . IPID,NS1,NS2,ISTRE,IADBUF,NUPARAM,IMAT,NNI,N0,
99 . IHBE,IREP,BUF,ISROT,IVISC,IGTYP,ISUBSTACK,
100 . ID_PLY,IPANG,IPPOS,IPTHK,OFFSET,ISELECT,MAT_ORTH,
101 . IXLAY,IXFEM,LAYNPT_MAX,NUMEL_DRAPE,SEDRAPE,NLAY_MAX,
102 . IPT_ALL,ISLICE,PTS,IPG,LENS,MPT0
103 INTEGER NN1,NN2,NN3,NN4,NN5,NN6,NN7,NN8,NN9,NN10,I1
104 INTEGER PID(MVSIZ),MAT(MVSIZ),IOK_PART(MVSIZ),JJ(15)
105 my_real ,DIMENSION(3,MVSIZ) :: strain
106 my_real ,DIMENSION(4*MVSIZ) :: xn,yn,zn,dxn,dyn,dzn
107 my_real ,DIMENSION(:,:) , ALLOCATABLE :: sige,sigm,epsm
108
109 TYPE(BUF_LAY_) ,POINTER :: BUFLY
110 TYPE(G_BUFEL_) ,POINTER :: GBUF
111 TYPE(L_BUFEL_) ,POINTER :: LBUF
112 my_real, DIMENSION(:) ,POINTER :: uparam,dir_a,dir_b
113 !
114 INTEGER, DIMENSION(:) , ALLOCATABLE :: MATLY !! MATLY(MVSIZ*LAY_MAX)
115 my_real, DIMENSION(:) , ALLOCATABLE :: thkly !! THKLY(MVSIZ*LAY_MAX*LAYNPT_MAX)
116 my_real, DIMENSION(:,:), ALLOCATABLE :: posly,thk_ly
117C-----------------------------------------------
118 ! material orthotropy defined in MATPARAM data structure
119 ! MATPARAM%ORTHOTROPY
120 ! -> 1 : ISOTROPIC
121 ! -> 2 : IRTHOTROPIC
122 ! -> 3 : ANISOTROPIC
123
124 offset = 0
125 value(1:5) = zero
126 iselect = 1
127 id_ply = 0
128 npg = 1
129
130 nn1 = 1
131 nn2 = nn1
132 nn3 = nn2
133 nn4 = nn3 + numelq
134 nn5 = nn4 + numelc
135 nn6 = nn5 + numeltg
136 nn7 = nn6
137 nn8 = nn7
138 nn9 = nn8
139 nn10= nn9
140C
141 DO i=1,numelc+numeltg
142 is_written_shell(i) = 0
143 ENDDO
144C
145 DO 490 ng=1,ngroup
146C IF(ANIM_K == 0.AND.IPARG(8,NG) == 1)GOTO 490
147 mlw = iparg(1,ng)
148 nel = iparg(2,ng)
149 nft = iparg(3,ng)
150 ity = iparg(5,ng)
151 igtyp = iparg(38,ng)
152 isrot = iparg(41,ng)
153 ixfem = iparg(54,ng)
154 isubstack = iparg(71,ng)
155 idrape = elbuf_tab(ng)%IDRAPE
156 npt = iabs(iparg(6,ng))
157 iok_part(1:nel) = 0
158!
159 DO i=1,15
160 jj(i) = nel*(i-1)
161 ENDDO
162!
163 IF (mlw /= 13) THEN
164C-----------------------------------------------
165C QUAD
166C-----------------------------------------------
167 IF(ity == 2)THEN
168 DO i=1,nel
169 n = i + nft
170 shell_tensor(1,offset+nft+i) = zero
171 shell_tensor(2,offset+nft+i) = zero
172 shell_tensor(3,offset+nft+i) = zero
173 ENDDO
174C-----------------------------------------------
175C COQUES
176C-----------------------------------------------
177 ELSEIF (ity == 3 .OR. ity == 7) THEN
178 gbuf => elbuf_tab(ng)%GBUF
179 nptr = elbuf_tab(ng)%NPTR
180 npts = elbuf_tab(ng)%NPTS
181 nlay = elbuf_tab(ng)%NLAY
182 npg = nptr*npts
183C
184 ihbe = iparg(23,ng)
185 IF (ity == 3) THEN
186 n0 = 0
187 nni = nn4
188 IF (ihbe == 11) npg = 4
189 ELSE
190 n0 = numelc
191 nni = nn5
192 IF (ihbe == 30) npg = 3 !dkt18
193 ENDIF
194c
195 IF (ity == 3) offset = 0
196 IF (ity == 7) offset = numelc
197c
198 DO i=1,nel
199 IF (ity == 3) THEN
200 id_elem(offset+nft+i) = ixc(nixc,nft+i)
201 ity_elem(offset+nft+i) = 3
202 IF( h3d_part(ipartc(nft+i)) == 1) iok_part(i) = 1
203 ELSEIF (ity == 7) THEN
204 id_elem(offset+nft+i) = ixtg(nixtg,nft+i)
205 ity_elem(offset+nft+i) = 7
206 IF( h3d_part(iparttg(nft+i)) == 1) iok_part(i) = 1
207 ENDIF
208 ENDDO
209c
210 IF (mlw == 0) GOTO 490
211C
212 a1 = zero
213 a2 = zero
214 a3 = zero
215 istre = 1
216 npt = iabs(iparg(6,ng))
217 mpt = max(1,npt)
218 mpt0 = mpt
219 IF (npt == 0) mpt = 0
220C
221 IF (igtyp == 10 .OR. igtyp == 11 .OR. igtyp == 16 .OR. igtyp == 17) THEN
222 npt = 1
223 mpt = npt
224 ELSEIF (igtyp == 51 .OR. igtyp == 52) THEN
225 IF(layer_input == -2) THEN
226 npt= elbuf_tab(ng)%BUFLY(1)%NPTT
227 ELSEIF(layer_input == -3) THEN
228 npt= elbuf_tab(ng)%BUFLY(nlay)%NPTT
229 ELSEIF(layer_input > 0 .AND. layer_input <= nlay) THEN
230 npt= elbuf_tab(ng)%BUFLY(layer_input)%NPTT
231 ENDIF
232 IF (ply_input > 0) THEN
233 DO j=1,nlay
234 id_ply = 0
235 IF (igtyp == 51) THEN
236 id_ply = igeo(1,stack%IGEO(2+j,isubstack))
237 ELSEIF (igtyp == 52) THEN
238 id_ply = ply_info(1,stack%IGEO(2+j,isubstack)-numstack)
239 ENDIF
240 IF (id_ply == ply_input ) THEN
241 npt= elbuf_tab(ng)%BUFLY(j)%NPTT
242 EXIT
243 ENDIF
244 ENDDO
245 ENDIF
246 mpt = max(1,npt)
247 ENDIF
248c
249 ilay = layer_input
250 ipt = ipt_input
251 iply = ply_input
252 IF (ilay == -2) ilay = 1 ! Lower
253 IF (ilay == -3) ilay = nlay ! Upper
254 IF (ipt == -2) ipt = 1 ! Lower
255 IF (igtyp == 51 .OR. igtyp == 52) THEN
256 IF (ipt == -3 .AND. ilay > 0) ipt = max(1,elbuf_tab(ng)%BUFLY(ilay)%NPTT) ! Upper
257 ELSE
258 IF (ipt == -3) ipt = max(1,npt) ! Upper
259 ENDIF
260C------------------------
261C STRESS
262C------------------------
263 IF (keyword == 'TENS/STRESS/MEMB' .OR.
264 . keyword == 'TENS/STRESS/BEND' .OR.
265 . keyword == 'TENS/STRESS' .OR.
266 . keyword == 'TENS/STRAIN' .OR.
267 . keyword == 'TENS/MSTRAIN' ) THEN
268 IF (ity == 3) THEN
269 ipid = ixc(6,nft+1)
270 DO i=1,nel
271 mat(i)=ixc(1,nft+i)
272 pid(i)=ixc(6,nft+i)
273 ENDDO
274 ELSE ! ITY == 7
275 ipid = ixtg(5,nft+1)
276 DO i=1,nel
277 mat(i)=ixtg(1,nft+i)
278 pid(i)=ixtg(5,nft+i)
279 ENDDO
280 ENDIF
281c
282 irep = igeo(6,ipid)
283 ENDIF
284 IF( keyword == 'TENS/STRAIN' .OR. keyword == 'TENS/MSTRAIN') THEN
285 !Npt_max
286 laynpt_max = 1
287 ixlay = 0
288 IF(igtyp == 51 .OR. igtyp == 52) THEN
289 DO ilay=1,nlay
290 laynpt_max = max(laynpt_max ,elbuf_tab(ng)%BUFLY(ilay)%NPTT)
291 ENDDO
292 ENDIF
293 nlay_max = max(nlay,npt)
294 ALLOCATE(matly(mvsiz*nlay_max), thkly(mvsiz*nlay_max*laynpt_max),
295 . posly(mvsiz,nlay_max*laynpt_max),thk_ly(nel,nlay_max*laynpt_max))
296 matly = 0
297 thkly = zero
298 posly = zero
299 thk_ly = zero
300c
301 ! computing position of slice or Ply
302 IF(ity == 7) THEN
303 numel_drape = numeltg_drape
304 sedrape = stdrape
305 CALL layini(
306 . elbuf_tab(ng),1 ,nel ,geo ,igeo ,
307 . mat ,pid ,thkly ,matly ,posly ,
308 . igtyp ,ixfem ,ixlay ,nlay ,mpt0 ,
309 . isubstack ,stack ,drape_sh3n ,nft ,gbuf%THK,
310 . nel ,thk_ly ,drapeg%INDX_SH3N ,sedrape,numel_drape )
311 ELSE ! ITY = 3
312 numel_drape = numelc_drape
313 sedrape = scdrape
314 CALL layini(
315 . elbuf_tab(ng),1 ,nel ,geo ,igeo ,
316 . mat ,pid ,thkly ,matly ,posly ,
317 . igtyp ,ixfem ,ixlay ,nlay ,mpt0 ,
318 . isubstack ,stack ,drape_sh4n ,nft ,gbuf%THK ,
319 . nel ,thk_ly ,drapeg%INDX_SH4N,sedrape ,numel_drape )
320 ENDIF
321 ENDIF
322C--------------------------------------------------
323 IF (keyword == 'TENS/STRESS/MEMB') THEN
324C--------------------------------------------------
325 DO i=1,nel
326 value(1) = gbuf%FOR(jj(1)+i)
327 value(2) = gbuf%FOR(jj(2)+i)
328 value(3) = gbuf%FOR(jj(3)+i)
329 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,
330 . VALUE)
331 ENDDO
332C--------------------------------------------------
333 ELSEIF (keyword == 'TENS/STRESS/BEND') THEN
334C--------------------------------------------------
335 DO i=1,nel
336 value(1) = gbuf%MOM(jj(1)+i)
337 value(2) = gbuf%MOM(jj(2)+i)
338 value(3) = gbuf%MOM(jj(3)+i)
339 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,
340 . VALUE)
341 ENDDO
342
343C--------------------------------------------------
344 ELSEIF (keyword == 'TENS/STRESS') THEN ! stress output in element coordinate system
345C--------------------------------------------------
346c
347 iselect = 0
348 ALLOCATE (sige(nel,3))
349 sige(1:nel,1:3) = zero
350c
351 IF (mpt == 0) THEN ! ILAYER=NULL, NPT=NULL
352 iselect = 1
353 IF (ipt_input == -2 ) THEN
354 factor = -six ! lower
355 ELSEIF (ipt_input == -3) THEN
356 factor = six ! upper
357 ELSE
358 factor = zero ! mem
359 END IF
360 DO i=1,nel
361 sige(i,1) = gbuf%FOR(jj(1)+i) + gbuf%MOM(jj(1)+i) * factor
362 sige(i,2) = gbuf%FOR(jj(2)+i) + gbuf%MOM(jj(2)+i) * factor
363 sige(i,3) = gbuf%FOR(jj(3)+i) + gbuf%MOM(jj(3)+i) * factor
364 ENDDO
365c
366 ELSE IF (ilay == -1 .AND. iply == -1 .AND. ipt == -1) THEN
367 iselect = 1
368 DO i=1,nel
369 sige(i,1) = gbuf%FOR(jj(1)+i)
370 sige(i,2) = gbuf%FOR(jj(2)+i)
371 sige(i,3) = gbuf%FOR(jj(3)+i)
372 ENDDO
373c
374 ELSEIF (ilay == -1 .AND. iply > 0 .AND. ipt > 0) THEN
375c /TENS/STRESS/PLY=.../NPT=...
376 DO j=1,nlay ! shells type 17,19,51 and 52 only
377 IF (igtyp == 17 .OR. igtyp == 19 .OR. igtyp == 51) THEN
378 id_ply = igeo(1,stack%IGEO(2+j,isubstack))
379 ELSE IF (igtyp == 52) THEN
380 id_ply = ply_info(1,stack%IGEO(2+j,isubstack) - numstack)
381 END IF
382 ilay = j
383 IF (id_ply == iply .AND. ipt <= elbuf_tab(ng)%BUFLY(ilay)%NPTT) THEN
384 imat = elbuf_tab(ng)%BUFLY(ilay)%IMAT
385 ivisc = mat_param(imat)%IVISC
386 iselect = 1
387 DO i=1,nel
388 DO ir=1,nptr
389 DO is=1,npts
390 lbuf => elbuf_tab(ng)%BUFLY(ilay)%LBUF(ir,is,ipt)
391 sige(i,1) = sige(i,1) + lbuf%SIG(jj(1) + i) / npg
392 sige(i,2) = sige(i,2) + lbuf%SIG(jj(2) + i) / npg
393 sige(i,3) = sige(i,3) + lbuf%SIG(jj(3) + i) / npg
394 ENDDO
395 ENDDO
396 ENDDO
397 IF (ivisc > 0) THEN
398 DO i=1,nel
399 DO ir=1,nptr
400 DO is=1,npts
401 lbuf => elbuf_tab(ng)%BUFLY(ilay)%LBUF(ir,is,ipt)
402 sige(i,1) = sige(i,1) + lbuf%VISC(jj(1) + i) / npg
403 sige(i,2) = sige(i,2) + lbuf%VISC(jj(2) + i) / npg
404 sige(i,3) = sige(i,3) + lbuf%VISC(jj(3) + i) / npg
405 ENDDO
406 ENDDO
407 ENDDO
408 END IF
409 mat_orth = mat_param(imat)%ORTHOTROPY
410 IF (mat_orth > 0) THEN
411 IF (idrape > 0 .AND. (igtyp == 51 .OR. igtyp == 52) ) THEN
412 dir_a => elbuf_tab(ng)%BUFLY(ilay)%LBUF_DIR(ipt)%DIRA
413 dir_b => elbuf_tab(ng)%BUFLY(ilay)%LBUF_DIR(ipt)%DIRB
414 ELSE
415 dir_a => elbuf_tab(ng)%BUFLY(ilay)%DIRA
416 dir_b => elbuf_tab(ng)%BUFLY(ilay)%DIRB
417 ENDIF
418 END IF
419 IF (mat_orth == 2) THEN
420 CALL uroto_tens2d(nel,sige,dir_a)
421 ELSE IF (mat_orth == 3) THEN ! anisotropic (law 58,158 only)
422 CALL uroto_tens2d_aniso(nel,sige,dir_a,dir_b)
423 END IF
424 EXIT
425 ENDIF ! ID_PLY & ipt
426 ENDDO ! NLAY
427c
428 ELSEIF (ilay > 0 .AND. ilay <= nlay .AND. iply == -1) THEN
429 ! /TENS/STRESS/LAYER=
430 IF (igtyp == 10 .OR. igtyp == 11 .OR. igtyp == 16) THEN ! shells type 10,11,16 only
431 iselect = 1
432 imat = elbuf_tab(ng)%BUFLY(ilay)%IMAT
433 ivisc = mat_param(imat)%IVISC
434 DO i=1,nel
435 DO ir=1,nptr
436 DO is=1,npts
437 lbuf => elbuf_tab(ng)%BUFLY(ilay)%LBUF(ir,is,1)
438 sige(i,1) = sige(i,1) + lbuf%SIG(jj(1) + i) / npg
439 sige(i,2) = sige(i,2) + lbuf%SIG(jj(2) + i) / npg
440 sige(i,3) = sige(i,3) + lbuf%SIG(jj(3) + i) / npg
441 ENDDO
442 ENDDO
443 ENDDO
444 IF (ivisc > 0) THEN
445 DO i=1,nel
446 DO ir=1,nptr
447 DO is=1,npts
448 lbuf => elbuf_tab(ng)%BUFLY(ilay)%LBUF(ir,is,1)
449 sige(i,1) = sige(i,1) + lbuf%VISC(jj(1) + i) / npg
450 sige(i,2) = sige(i,2) + lbuf%VISC(jj(2) + i) / npg
451 sige(i,3) = sige(i,3) + lbuf%VISC(jj(3) + i) / npg
452 ENDDO
453 ENDDO
454 ENDDO
455 END IF
456 mat_orth = mat_param(imat)%ORTHOTROPY
457 IF (mat_orth > 0) THEN
458 dir_a => elbuf_tab(ng)%BUFLY(ilay)%DIRA
459 dir_b => elbuf_tab(ng)%BUFLY(ilay)%DIRB
460 END IF
461 IF (mat_orth == 2) THEN ! standard orthotropic rotation
462 CALL uroto_tens2d(nel,sige,dir_a)
463 ELSE IF (mat_orth == 3) THEN ! anisotropic (law 58,158 only)
464 CALL uroto_tens2d_aniso(nel,sige,dir_a,dir_b)
465 END IF
466 ENDIF
467
468 ELSEIF (ipt > 0 .AND. ilay ==-1 .AND. iply == -1) THEN ! shells type 1,9 only
469c /TENS/STRESS/NPT=
470 IF (igtyp == 1 .OR. igtyp == 9) THEN
471 IF (ipt <= elbuf_tab(ng)%BUFLY(1)%NPTT) THEN
472 iselect = 1
473 imat = elbuf_tab(ng)%BUFLY(1)%IMAT
474 ivisc = mat_param(imat)%IVISC
475 DO i=1,nel
476 DO ir=1,nptr
477 DO is=1,npts
478 lbuf => elbuf_tab(ng)%BUFLY(1)%LBUF(ir,is,ipt)
479 sige(i,1) = sige(i,1) + lbuf%SIG(jj(1) + i) / npg
480 sige(i,2) = sige(i,2) + lbuf%SIG(jj(2) + i) / npg
481 sige(i,3) = sige(i,3) + lbuf%SIG(jj(3) + i) / npg
482 ENDDO
483 ENDDO
484 ENDDO
485 IF (ivisc > 0) THEN
486 DO i=1,nel
487 DO ir=1,nptr
488 DO is=1,npts
489 lbuf => elbuf_tab(ng)%BUFLY(1)%LBUF(ir,is,ipt)
490 sige(i,1) = sige(i,1) + lbuf%VISC(jj(1) + i) / npg
491 sige(i,2) = sige(i,2) + lbuf%VISC(jj(2) + i) / npg
492 sige(i,3) = sige(i,3) + lbuf%VISC(jj(3) + i) / npg
493 ENDDO
494 ENDDO
495 ENDDO
496 END IF
497 mat_orth = mat_param(imat)%ORTHOTROPY
498 IF (mat_orth == 2) THEN
499 dir_a => elbuf_tab(ng)%BUFLY(1)%DIRA
500 CALL uroto_tens2d(nel,sige,dir_a)
501 END IF
502 ENDIF
503 ENDIF
504 ENDIF
505c---
506 IF (iselect == 1) THEN
508 . iok_part ,iselect ,nel ,offset ,nft ,
509 . is_written_shell,shell_tensor,sige )
510 END IF
511c
512 DEALLOCATE (sige)
513c--------------------------------------------------
514 ELSEIF (keyword == 'TENS/MSTRESS') THEN ! stress output in material (orthotropic) coordinates
515c--------------------------------------------------
516 iselect = 0
517 ALLOCATE (sigm(nel,3))
518 sigm(1:nel,1:3) = zero
519c
520 IF (ilay == -1 .AND. iply > 0 .AND. ipt > 0) THEN
521c /TENS/MSTRESS/PLY=.../NPT=...
522 DO j=1,nlay ! shells type 17,19,51 and 52 only
523 IF (igtyp == 17 .OR. igtyp == 19 .OR. igtyp == 51) THEN
524 id_ply = igeo(1,stack%IGEO(2+j,isubstack))
525 ELSE IF (igtyp == 52) THEN
526 id_ply = ply_info(1,stack%IGEO(2+j,isubstack) - numstack)
527 END IF
528
529 IF (id_ply == iply) THEN
530 ilay = j
531 imat = elbuf_tab(ng)%BUFLY(ilay)%IMAT
532 ivisc = mat_param(imat)%IVISC
533 IF (ipt <= elbuf_tab(ng)%BUFLY(ilay)%NPTT) THEN
534 iselect = 1
535 DO i=1,nel
536 DO ir=1,nptr
537 DO is=1,npts
538 lbuf => elbuf_tab(ng)%BUFLY(ilay)%LBUF(ir,is,ipt)
539 sigm(i,1) = sigm(i,1) + lbuf%SIG(jj(1) + i) / npg
540 sigm(i,2) = sigm(i,2) + lbuf%SIG(jj(2) + i) / npg
541 sigm(i,3) = sigm(i,3) + lbuf%SIG(jj(3) + i) / npg
542 ENDDO
543 ENDDO
544 ENDDO
545 IF (ivisc > 0) THEN
546 DO i=1,nel
547 DO ir=1,nptr
548 DO is=1,npts
549 lbuf => elbuf_tab(ng)%BUFLY(ilay)%LBUF(ir,is,ipt)
550 sigm(i,1) = sigm(i,1) + lbuf%VISC(jj(1) + i) / npg
551 sigm(i,2) = sigm(i,2) + lbuf%VISC(jj(2) + i) / npg
552 sigm(i,3) = sigm(i,3) + lbuf%VISC(jj(3) + i) / npg
553 ENDDO
554 ENDDO
555 ENDDO
556 END IF
557 ENDIF
558 ENDIF
559 ENDDO ! NLAY
560c
561 ELSEIF (ilay > 0 .AND. ilay <= nlay .AND. iply == -1 .AND. ipt == -1) THEN
562c /TENS/MSTRESS/LAYER=
563 IF (igtyp == 10 .OR. igtyp == 11 .OR. igtyp == 16) THEN ! shells type 10,11,16 only
564 iselect = 1
565 imat = elbuf_tab(ng)%BUFLY(ilay)%IMAT
566 ivisc = mat_param(imat)%IVISC
567 DO i=1,nel
568 DO ir=1,nptr
569 DO is=1,npts
570 lbuf => elbuf_tab(ng)%BUFLY(ilay)%LBUF(ir,is,1)
571 sigm(i,1) = sigm(i,1) + lbuf%SIG(jj(1) + i) / npg
572 sigm(i,2) = sigm(i,2) + lbuf%SIG(jj(2) + i) / npg
573 sigm(i,3) = sigm(i,3) + lbuf%SIG(jj(3) + i) / npg
574 ENDDO
575 ENDDO
576 ENDDO
577 IF (ivisc > 0) THEN
578 DO i=1,nel
579 DO ir=1,nptr
580 DO is=1,npts
581 lbuf => elbuf_tab(ng)%BUFLY(ilay)%LBUF(ir,is,ipt)
582 sigm(i,1) = sigm(i,1) + lbuf%VISC(jj(1) + i) / npg
583 sigm(i,2) = sigm(i,2) + lbuf%VISC(jj(2) + i) / npg
584 sigm(i,3) = sigm(i,3) + lbuf%VISC(jj(3) + i) / npg
585 ENDDO
586 ENDDO
587 ENDDO
588 END IF
589 ENDIF
590
591 ELSEIF (ipt > 0 .AND. ilay ==-1 .AND. iply == -1) THEN ! shells type 1,9 only
592c /TENS/MSTRESS/NPT=
593 IF (igtyp == 1 .OR. igtyp == 9) THEN
594 IF (ipt <= elbuf_tab(ng)%BUFLY(1)%NPTT) THEN
595 iselect = 1
596 imat = elbuf_tab(ng)%BUFLY(ilay)%IMAT
597 ivisc = mat_param(imat)%IVISC
598 DO i=1,nel
599 DO ir=1,nptr
600 DO is=1,npts
601 lbuf => elbuf_tab(ng)%BUFLY(1)%LBUF(ir,is,ipt)
602 sigm(i,1) = sigm(i,1) + lbuf%SIG(jj(1) + i) / npg
603 sigm(i,2) = sigm(i,2) + lbuf%SIG(jj(2) + i) / npg
604 sigm(i,3) = sigm(i,3) + lbuf%SIG(jj(3) + i) / npg
605 ENDDO
606 ENDDO
607 ENDDO
608 IF (ivisc > 0) THEN
609 DO i=1,nel
610 DO ir=1,nptr
611 DO is=1,npts
612 lbuf => elbuf_tab(ng)%BUFLY(1)%LBUF(ir,is,ipt)
613 sigm(i,1) = sigm(i,1) + lbuf%VISC(jj(1) + i) / npg
614 sigm(i,2) = sigm(i,2) + lbuf%VISC(jj(2) + i) / npg
615 sigm(i,3) = sigm(i,3) + lbuf%VISC(jj(3) + i) / npg
616 ENDDO
617 ENDDO
618 ENDDO
619 END IF
620
621 ENDIF
622 ENDIF
623 ENDIF
624c---
626 . iok_part ,iselect ,nel ,offset ,nft ,
627 . is_written_shell,shell_tensor,sigm )
628
629 DEALLOCATE (sigm)
630C--------------------------------------------------
631 ELSE IF (keyword == 'TENS/STRAIN/MEMB') THEN
632C--------------------------------------------------
633 DO i=1,nel
634 n = i + nft
635 thk = gbuf%THK(i)
636 j = el2fa(nni+n)
637 value(1) = gbuf%STRA(jj(1)+i)
638 value(2) = gbuf%STRA(jj(2)+i)
639 value(3) = gbuf%STRA(jj(3)+i)
640 value(3) = value(3) * half
641 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,
642 . shell_tensor,i,offset,nft,VALUE)
643 ENDDO
644C--------------------------------------------------
645 ELSEIF (keyword == 'TENS/STRAIN/BEND') THEN ! bend
646C--------------------------------------------------
647 DO i=1,nel
648 n = i + nft
649 thk = gbuf%THK(i)
650 j = el2fa(nni+n)
651 value(1) = gbuf%STRA(jj(6)+i) * thk
652 value(2) = gbuf%STRA(jj(7)+i) * thk
653 value(3) = gbuf%STRA(jj(8)+i) * thk
654 value(3) = value(3) * half
655 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,
656 . shell_tensor,i,offset,nft,VALUE)
657 ENDDO
658C--------------------------------------------------
659 ELSEIF (keyword == 'TENS/STRAIN') THEN ! strain tensor output in element coordinate system
660C--------------------------------------------------
661 IF (mpt == 0) THEN ! ILAYER=NULL, NPT=NULL
662 iselect = 1
663 DO i=1,nel
664 IF (ipt == 1) THEN
665 factor = -half*gbuf%THK(i)
666 ELSE
667 factor = half*gbuf%THK(i)
668 ENDIF
669 value(1) = gbuf%STRA(jj(1)+i) + factor*gbuf%STRA(jj(6)+i)
670 value(2) = gbuf%STRA(jj(2)+i) + factor*gbuf%STRA(jj(7)+i)
671 value(3) = gbuf%STRA(jj(3)+i) + factor*gbuf%STRA(jj(8)+i)
672 value(3) = value(3) * half
673 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,
674 . shell_tensor,i,offset,nft,VALUE)
675 ENDDO
676c
677 ELSE IF (ilay == -1 .AND. iply == -1 .AND. ipt == -1) THEN
678 DO i=1,nel
679 value(1) = gbuf%STRA(jj(1)+i)
680 value(2) = gbuf%STRA(jj(2)+i)
681 value(3) = gbuf%STRA(jj(3)+i)
682 value(3) = value(3) * half
683 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,
684 . shell_tensor,i,offset,nft,VALUE)
685 ENDDO
686c
687c PLY=IPLY NPT=IPT
688
689 ELSE IF (iply > 0 .AND. ipt > 0) THEN
690 IF (igtyp == 17 .OR. igtyp == 19 .OR. igtyp == 51 .OR. igtyp == 52) THEN
691 ipang = 1
692 ipthk = ipang + nlay
693 ippos = ipthk + nlay
694 ipt_all = 0
695 DO j=1,nlay
696 bufly => elbuf_tab(ng)%BUFLY(j)
697 nptt = bufly%NPTT
698 IF (igtyp == 17 .OR. igtyp == 51) THEN
699 id_ply = igeo(1,stack%IGEO(2+j,isubstack))
700 ELSEIF (igtyp == 52) THEN
701 id_ply = ply_info(1,stack%IGEO(2+j,isubstack)-numstack)
702 ENDIF
703 IF (id_ply == iply .AND. ipt <= nptt) THEN
704 islice = ipt_all + ipt
705 IF(npg > 1) THEN
706 lens = nel*gbuf%G_STRPG/npg
707 DO i=1,nel
708 n = i + nft
709 thk = gbuf%THK(i)
710 factor = posly(i,islice)
711 value(1:3) = zero
712 DO ipg = 1, npg
713 pts = (ipg -1)*lens
714 value(1) = value(1)+ gbuf%STRPG(pts + jj(1)+i) + factor*gbuf%STRPG(pts + jj(6)+i) * thk
715 value(2) = value(2)+ gbuf%STRPG(pts + jj(2)+i) + factor*gbuf%STRPG(pts + jj(7)+i) * thk
716 value(3) = value(3)+ gbuf%STRPG(pts + jj(3)+i) + factor*gbuf%STRPG(pts + jj(8)+i) * thk
717 ENDDO
718 value(3) = value(3) * half
719 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,
720 . shell_tensor,i,offset,nft,VALUE)
721 ENDDO
722 ELSE
723 DO i=1,nel
724 n = i + nft
725 thk = gbuf%THK(i)
726 factor = posly(i,islice)
727 value(1) = gbuf%STRA(jj(1)+i) + factor*gbuf%STRA(jj(6)+i) * thk
728 value(2) = gbuf%STRA(jj(2)+i) + factor*gbuf%STRA(jj(7)+i) * thk
729 value(3) = gbuf%STRA(jj(3)+i) + factor*gbuf%STRA(jj(8)+i) * thk
730 value(3) = value(3) * half
731 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,
732 . shell_tensor,i,offset,nft,VALUE)
733 ENDDO
734 ENDIF ! NPG
735 ENDIF
736 ipt_all = ipt_all + nptt
737 ENDDO
738 ENDIF
739c
740 ELSEIF (ilay > 0 .AND. ilay <= nlay .AND. iply == -1) THEN
741c ILAYER=IL PLY=NULL NPT=NULL
742 IF (igtyp == 10 .OR. igtyp == 11 .OR. igtyp == 16) THEN
743 IF (npg > 1) THEN
744 lens = nel*gbuf%G_STRPG/npg
745 DO i=1,nel
746 n = i + nft
747 thk = gbuf%THK(i)
748 factor = posly(i,ilay)
749 j = el2fa(nni+n)
750 value(1:3) = zero
751 DO ipg = 1, npg
752 pts = (ipg -1)*lens
753 value(1) = value(1)+ gbuf%STRPG(pts + jj(1)+i) + factor*gbuf%STRPG(pts + jj(6)+i) * thk
754 value(2) = value(2)+ gbuf%STRPG(pts + jj(2)+i) + factor*gbuf%STRPG(pts + jj(7)+i) * thk
755 value(3) = value(3)+ gbuf%STRPG(pts + jj(3)+i) + factor*gbuf%STRPG(pts + jj(8)+i) * thk
756 ENDDO
757 value(1) = value(1)/npg
758 value(2) = value(2)/npg
759 value(3) = value(3)/npg
760 value(3) = value(3) * half
761 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,
762 . shell_tensor,i,offset,nft,VALUE)
763 ENDDO
764 ELSE
765 DO i=1,nel
766 n = i + nft
767 thk = gbuf%THK(i)
768 factor = posly(i,ilay)
769 j = el2fa(nni+n)
770 value(1) = gbuf%STRA(jj(1)+i) + factor*gbuf%STRA(jj(6)+i) * thk
771 value(2) = gbuf%STRA(jj(2)+i) + factor*gbuf%STRA(jj(7)+i) * thk
772 value(3) = gbuf%STRA(jj(3)+i) + factor*gbuf%STRA(jj(8)+i) * thk
773 value(3) = value(3) * half
774 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,
775 . shell_tensor,i,offset,nft,VALUE)
776 ENDDO
777 ENDIF ! NPG
778 ENDIF
779c
780 ELSEIF (ipt <= mpt .AND. ipt > 0) THEN
781c NPT=IPT PLY=NULL ILAYER=NULL
782 IF (igtyp == 1 .OR. igtyp == 9) THEN
783 IF(npg > 1) THEN
784 lens = nel*gbuf%G_STRPG/npg
785 DO i=1,nel
786 n = i + nft
787 thk = gbuf%THK(i)
788 factor = posly(i,ipt)
789 j = el2fa(nni+n)
790 value(1:3) = zero
791 DO ipg =1,npg
792 pts = (ipg -1)*lens
793 value(1) = VALUE(1)+ gbuf%STRPG(pts + jj(1)+i) + factor*gbuf%STRPG(pts + jj(6)+i) * thk
794 value(2) = value(2)+ gbuf%STRPG(pts + jj(2)+i) + factor*gbuf%STRPG(pts + jj(7)+i) * thk
795 value(3) = value(3)+ gbuf%STRPG(pts + jj(3)+i) + factor*gbuf%STRPG(pts + jj(8)+i) * thk
796 ENDDO
797 value(3) = value(3) * half
798 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,
799 . shell_tensor,i,offset,nft,VALUE)
800 ENDDO
801 ELSE
802 DO i=1,nel
803 n = i + nft
804 thk = gbuf%THK(i)
805 factor = posly(i,ipt)
806 j = el2fa(nni+n)
807 value(1) = gbuf%STRA(jj(1)+i) + factor*gbuf%STRA(jj(6)+i) * thk
808 value(2) = gbuf%STRA(jj(2)+i) + factor*gbuf%STRA(jj(7)+i) * thk
809 value(3) = gbuf%STRA(jj(3)+i) + factor*gbuf%STRA(jj(8)+i) * thk
810 value(3) = value(3) * half
811 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,
812 . shell_tensor,i,offset,nft,VALUE)
813 ENDDO
814 ENDIF ! NPG
815 ENDIF
816 ENDIF
817 DEALLOCATE(matly, thkly,posly,thk_ly)
818C--------------------------------------------------
819 ELSEIF (keyword == 'TENS/MSTRAIN') THEN ! strain tensor output in material coordinate system
820C--------------------------------------------------
821 ALLOCATE (epsm(nel,3))
822 epsm(:,:) = zero
823c
824 IF (iply > 0 .AND. ipt > 0) THEN
825
826 IF (igtyp == 17 .OR. igtyp == 19 .OR. igtyp == 51 .OR. igtyp == 52) THEN
827 ipang = 1
828 ipthk = ipang + nlay
829 ippos = ipthk + nlay
830 ipt_all = 0
831 DO j=1,nlay
832 bufly => elbuf_tab(ng)%BUFLY(j)
833 nptt = bufly%NPTT
834 IF (igtyp == 17 .OR. igtyp == 19 .OR. igtyp == 51) THEN
835 id_ply = igeo(1,stack%IGEO(2+j,isubstack))
836 ELSEIF (igtyp == 52) THEN
837 id_ply = ply_info(1,stack%IGEO(2+j,isubstack)-numstack)
838 ENDIF
839 IF (id_ply == iply .AND. ipt <= nptt) THEN
840 ilay = j
841 islice = ipt_all + ipt
842 IF(npg > 1) THEN
843 lens = nel*gbuf%G_STRPG/npg
844 DO i=1,nel
845 thk = gbuf%THK(i)
846 factor = posly(i,islice)
847 epsm(i,1:3) = zero
848 DO ipg = 1, npg
849 pts = (ipg -1)*lens
850 epsm(i,1) = epsm(i,1)+ gbuf%STRPG(pts + jj(1)+i) + factor*gbuf%STRPG(pts + jj(6)+i) * thk
851 epsm(i,2) = epsm(i,2)+ gbuf%STRPG(pts + jj(2)+i) + factor*gbuf%STRPG(pts + jj(7)+i) * thk
852 epsm(i,3) = epsm(i,3)+ gbuf%STRPG(pts + jj(3)+i) + factor*gbuf%STRPG(pts + jj(8)+i) * thk
853 ENDDO
854 epsm(i,1) = epsm(i,1)/npg
855 epsm(i,2) = epsm(i,2)/npg
856 epsm(i,3) = half*epsm(i,3)/npg
857 ENDDO
858 ELSE
859 DO i=1,nel
860 thk = gbuf%THK(i)
861 factor = posly(i,islice)
862 epsm(i,1) = gbuf%STRA(jj(1)+i) + factor*gbuf%STRA(jj(6)+i) * thk
863 epsm(i,2) = gbuf%STRA(jj(2)+i) + factor*gbuf%STRA(jj(7)+i) * thk
864 epsm(i,3) = gbuf%STRA(jj(3)+i) + factor*gbuf%STRA(jj(8)+i) * thk
865 epsm(i,3) = epsm(i,3) * half
866 ENDDO
867 ENDIF ! NPG
868 imat = elbuf_tab(ng)%BUFLY(ilay)%IMAT
869 mat_orth = mat_param(imat)%ORTHOTROPY
870 IF (mat_orth > 0) THEN
871 IF (idrape > 0 .AND. (igtyp == 51 .OR. igtyp == 52) ) THEN
872 dir_a => elbuf_tab(ng)%BUFLY(ilay)%LBUF_DIR(ipt)%DIRA
873 dir_b => elbuf_tab(ng)%BUFLY(ilay)%LBUF_DIR(ipt)%DIRB
874 ELSE
875 dir_a => elbuf_tab(ng)%BUFLY(ilay)%DIRA
876 dir_b => elbuf_tab(ng)%BUFLY(ilay)%DIRB
877 ENDIF
878 END IF
879 IF (mat_orth == 2) THEN
880 CALL roto_tens2d(nel,epsm,dir_a)
881 ELSE IF (mat_orth == 3) THEN ! anisotropic (law 58,158 only)
882 CALL roto_tens2d_aniso(nel,epsm,dir_a,dir_b)
883 END IF
884 ENDIF
885 ipt_all = ipt_all + nptt
886 ENDDO
887 ENDIF
888c
889 ELSEIF (ilay > 0 .AND. ilay <= nlay .AND. iply == -1 .AND. ipt == -1) THEN
890c PLY=NULL ILAYER=IL NPT=NULL
891 IF (igtyp == 10 .OR. igtyp == 11 .OR. igtyp == 16) THEN
892 IF(npg > 1) THEN
893 lens = nel*gbuf%G_STRPG/npg
894 DO i=1,nel
895 thk = gbuf%THK(i)
896 factor = posly(i,ilay)
897 epsm(i,1:3) = zero
898 DO ipg = 1, npg
899 pts = (ipg -1)*lens
900 epsm(i,1) = epsm(i,1)+ gbuf%STRPG(pts + jj(1)+i) + factor*gbuf%STRPG(pts + jj(6)+i) * thk
901 epsm(i,2) = epsm(i,2)+ gbuf%STRPG(pts + jj(2)+i) + factor*gbuf%STRPG(pts + jj(7)+i) * thk
902 epsm(i,3) = epsm(i,3)+ gbuf%STRPG(pts + jj(3)+i) + factor*gbuf%STRPG(pts + jj(8)+i) * thk
903 ENDDO
904 epsm(i,1) = epsm(i,1)/npg
905 epsm(i,2) = epsm(i,2)/npg
906 epsm(i,3) = half*epsm(i,3)/npg
907 ENDDO
908 ELSE
909 DO i=1,nel
910 thk = gbuf%THK(i)
911 factor = posly(i,ilay)
912 epsm(i,1) = gbuf%STRA(jj(1)+i) + factor*gbuf%STRA(jj(6)+i) * thk
913 epsm(i,2) = gbuf%STRA(jj(2)+i) + factor*gbuf%STRA(jj(7)+i) * thk
914 epsm(i,3) = gbuf%STRA(jj(3)+i) + factor*gbuf%STRA(jj(8)+i) * thk
915 epsm(i,3) = epsm(i,3) * half
916 ENDDO
917 ENDIF
918 imat = elbuf_tab(ng)%BUFLY(ilay)%IMAT
919 mat_orth = mat_param(imat)%ORTHOTROPY
920 IF (mat_orth > 0) THEN
921 dir_a => elbuf_tab(ng)%BUFLY(ilay)%DIRA
922 dir_b => elbuf_tab(ng)%BUFLY(ilay)%DIRB
923 END IF
924 IF (mat_orth == 2) THEN
925 CALL roto_tens2d(nel,epsm,dir_a)
926 ELSE IF (mat_orth == 3) THEN ! anisotropic (law 58,158 only)
927 CALL roto_tens2d_aniso(nel,epsm,dir_a,dir_b)
928 END IF
929 ENDIF
930c
931 ELSEIF (ipt > 0 .AND. ipt <= mpt .AND. iply == -1 .AND. ilay == -1) THEN
932c PLY=NULL ILAYER=NULL NPT=IPT
933 IF (igtyp == 1 .OR. igtyp == 9) THEN
934 IF(npg > 1) THEN
935 lens = nel*gbuf%G_STRPG/npg
936 DO i=1,nel
937 thk = gbuf%THK(i)
938 factor = posly(i,ipt)
939 epsm(i,1:3) = zero
940 DO ipg = 1, npg
941 pts = (ipg -1)*lens
942 epsm(i,1) = epsm(i,1)+ gbuf%STRPG(pts + jj(1)+i) + factor*gbuf%STRPG(pts + jj(6)+i) * thk
943 epsm(i,2) = epsm(i,2)+ gbuf%STRPG(pts + jj(2)+i) + factor*gbuf%STRPG(pts + jj(7)+i) * thk
944 epsm(i,3) = epsm(i,3)+ gbuf%STRPG(pts + jj(3)+i) + factor*gbuf%STRPG(pts + jj(8)+i) * thk
945 ENDDO
946 epsm(i,1) = epsm(i,1)/npg
947 epsm(i,2) = epsm(i,2)/npg
948 epsm(i,3) = half*epsm(i,3)/npg
949 ENDDO
950 ELSE
951 DO i=1,nel
952 thk = gbuf%THK(i)
953 factor = posly(i,ipt)
954 epsm(i,1) = gbuf%STRA(jj(1)+i) + factor*gbuf%STRA(jj(6)+i) * thk
955 epsm(i,2) = gbuf%STRA(jj(2)+i) + factor*gbuf%STRA(jj(7)+i) * thk
956 epsm(i,3) = gbuf%STRA(jj(3)+i) + factor*gbuf%STRA(jj(8)+i) * thk
957 epsm(i,3) = epsm(i,3) * half
958 ENDDO
959 ENDIF
960 imat = elbuf_tab(ng)%BUFLY(1)%IMAT
961 mat_orth = mat_param(imat)%ORTHOTROPY
962 IF (mat_orth == 2) THEN
963 dir_a => elbuf_tab(ng)%BUFLY(1)%DIRA
964 CALL roto_tens2d(nel,epsm,dir_a)
965 END IF
966 ENDIF
967 ENDIF
968c---
970 . iok_part ,iselect ,nel ,offset ,nft ,
971 . is_written_shell,shell_tensor,epsm )
972
973 DEALLOCATE (epsm)
974 DEALLOCATE(matly, thkly,posly,thk_ly)
975C--------------------------------------------------
976 ELSEIF (keyword == 'TENS/EPSDOT/MEMB') THEN
977C--------------------------------------------------
978 a1 = one
979 a2 = zero
980 DO i=1,nel
981 thk = gbuf%THK(i)
982 value(1) = a1*epsdot(1,i+nft+offset) + a2*epsdot(4,i+nft+offset)*thk
983 value(2) = a1*epsdot(2,i+nft+offset) + a2*epsdot(5,i+nft+offset)*thk
984 value(3) = (a1*epsdot(3,i+nft+offset) + a2*epsdot(6,i+nft+offset)*thk)* half
985 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,
986 . VALUE)
987 ENDDO
988C--------------------------------------------------
989 ELSEIF (keyword == 'TENS/EPSDOT/BEND') THEN
990C--------------------------------------------------
991 DO i=1,nel
992 thk = gbuf%THK(i)
993 value(1) = epsdot(4,i+nft+offset)
994 value(2) = epsdot(5,i+nft+offset)
995 value(3) = epsdot(6,i+nft+offset) * half
996 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,
997 . VALUE)
998 ENDDO
999C--------------------------------------------------
1000 ELSEIF (keyword == 'TENS/EPSDOT') THEN
1001C--------------------------------------------------
1002c ILAYER=NULL NPT=NULL
1003 IF ( ilay == -1 .AND. ipt == -1 .AND. iply == -1) THEN
1004 a1 = one
1005 a2 = zero
1006 DO i=1,nel
1007 thk = gbuf%THK(i)
1008 value(1) = a1*epsdot(1,i+nft+offset) + a2*epsdot(4,i+nft+offset)*thk
1009 value(2) = a1*epsdot(2,i+nft+offset) + a2*epsdot(5,i+nft+offset)*thk
1010 value(3) = (a1*epsdot(3,i+nft+offset) + a2*epsdot(6,i+nft+offset)*thk)* half
1011 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,
1012 . shell_tensor,i,offset,nft,VALUE)
1013 ENDDO
1014c PLY=IPLY NPT=IPT
1015 ELSEIF ( iply > 0 .AND. ipt <= mpt .AND. ipt > 0 ) THEN
1016 IF (igtyp == 17 .OR. igtyp == 51 .OR. igtyp == 52) THEN
1017 ipang = 1
1018 ipthk = ipang + nlay
1019 ippos = ipthk + nlay
1020 DO j=1,nlay
1021 IF (igtyp == 17 .OR. igtyp == 51) THEN
1022 id_ply = igeo(1,stack%IGEO(2+j,isubstack))
1023 ELSEIF (igtyp == 52) THEN
1024 id_ply = ply_info(1,stack%IGEO(2+j,isubstack)-numstack)
1025 ENDIF
1026 bufly => elbuf_tab(ng)%BUFLY(j)
1027 nptt = bufly%NPTT
1028 IF (id_ply == iply .AND. ipt <= nptt) THEN
1029 a2 = stack%GEO(ippos+j,isubstack)+
1030 . half*(((2*ipt-one)/nptt)-one) *
1031 . stack%GEO(ipthk+j,isubstack)
1032 DO i=1,nel
1033 thk = gbuf%THK(i)
1034 value(1) = epsdot(1,i+nft+offset) + a2*epsdot(4,i+nft+offset)*thk
1035 value(2) = epsdot(2,i+nft+offset) + a2*epsdot(5,i+nft+offset)*thk
1036 value(3) =(epsdot(3,i+nft+offset) + a2*epsdot(6,i+nft+offset)*thk)* half
1037 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,
1038 . shell_tensor,i,offset,nft,VALUE)
1039 ENDDO
1040 ENDIF
1041 ENDDO
1042 ENDIF
1043
1044c PLY=NULL ILAYER=ILAY NPT=IPT
1045 ELSEIF (iply == -1 .AND. ilay <= nlay .AND. ilay > 0 .AND. ipt <= mpt .AND. ipt > 0 ) THEN
1046 IF (igtyp == 51 .OR. igtyp == 52) THEN
1047 a1 = zero
1048 a2 = zero
1049 ns1 = 8
1050 ns2 = 8
1051 ipang = 1
1052 ipthk = ipang + nlay
1053 ippos = ipthk + nlay
1054 IF (igtyp == 17 .OR. igtyp == 51) THEN
1055 id_ply = igeo(1,stack%IGEO(2+ilay,isubstack))
1056 ELSEIF (igtyp == 52) THEN
1057 id_ply = ply_info(1,stack%IGEO(2+ilay,isubstack)-numstack)
1058 ENDIF
1059 bufly => elbuf_tab(ng)%BUFLY(ilay)
1060 nptt = bufly%NPTT
1061 IF (ipt <= nptt) THEN
1062 a1 = one
1063 a2 = stack%GEO(ippos+ilay,isubstack)+
1064 . half*(((2*ipt-one)/nptt)-one) *
1065 . stack%GEO(ipthk+ilay,isubstack)
1066 DO i=1,nel
1067 n = i + nft
1068 thk = gbuf%THK(i)
1069 value(1) = a1*epsdot(1,i+nft+offset) + a2*epsdot(4,i+nft+offset)*thk
1070 value(2) = a1*epsdot(2,i+nft+offset) + a2*epsdot(5,i+nft+offset)*thk
1071 value(3) = (a1*epsdot(3,i+nft+offset) + a2*epsdot(6,i+nft+offset)*thk)* half
1072 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,
1073 . VALUE)
1074 ENDDO
1075 ENDIF
1076 ENDIF
1077c PLY=NULL ILAYER=IL NPT=NULL
1078 ELSEIF (iply == -1 .AND. ilay <= nlay .AND. ilay > 0 .AND. ipt == -1 ) THEN
1079 IF (igtyp == 10 .OR. igtyp == 11 .OR. igtyp == 16 .OR. igtyp == 17) THEN
1080 a1 = one
1081 a2 = half*(((2*ilay-one)/nlay)-one)
1082 DO i=1,nel
1083 n = i + nft
1084 thk = gbuf%THK(i)
1085 value(1) = a1*epsdot(1,i+nft+offset) + a2*epsdot(4,i+nft+offset)*thk
1086 value(2) = a1*epsdot(2,i+nft+offset) + a2*epsdot(5,i+nft+offset)*thk
1087 value(3) = (a1*epsdot(3,i+nft+offset) + a2*epsdot(6,i+nft+offset)*thk)* half
1088 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,
1089 . VALUE)
1090 ENDDO
1091 ELSEIF (igtyp == 51 .OR. igtyp == 52) THEN
1092 a1 = one
1093 a2 = stack%GEO(ippos+ilay,isubstack)+
1094 . half*(((2*ipt-one)/nptt)-one) *
1095 . stack%GEO(ipthk+ilay,isubstack)
1096 DO i=1,nel
1097 n = i + nft
1098 thk = gbuf%THK(i)
1099 value(1) = a1*epsdot(1,i+nft+offset) + a2*epsdot(4,i+nft+offset)*thk
1100 value(2) = a1*epsdot(2,i+nft+offset) + a2*epsdot(5,i+nft+offset)*thk
1101 value(3) = (a1*epsdot(3,i+nft+offset) + a2*epsdot(6,i+nft+offset)*thk)* half
1102 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,
1103 . VALUE)
1104 ENDDO
1105 ENDIF
1106c PLY=NULL ILAYER=NULL NPT=IPT
1107 ELSEIF ( ipt <= mpt .AND. ipt > 0) THEN
1108 a1 = one
1109 a2 = half*(((2*ipt-one)/mpt)-one)
1110 IF (igtyp == 1 .OR. igtyp == 9) THEN
1111 DO i=1,nel
1112 thk = gbuf%THK(i)
1113 value(1) = a1*epsdot(1,i+nft+offset) + a2*epsdot(4,i+nft+offset)*thk
1114 value(2) = a1*epsdot(2,i+nft+offset) + a2*epsdot(5,i+nft+offset)*thk
1115 value(3) = (a1*epsdot(3,i+nft+offset) + a2*epsdot(6,i+nft+offset)*thk)* half
1116 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,
1117 . VALUE)
1118 ENDDO
1119 ENDIF
1120 ENDIF
1121C--------------------------------------------------
1122 ELSE IF (keyword == 'TENS/STRAIN_ENG') THEN
1123C--------------------------------------------------
1124 IF (ity == 3 ) THEN !4n
1125 DO i=1,nel
1126 n = i + nft
1127 nni = ixc(2,n)
1128 j = 4*(i-1) +1
1129 xn(j)=x(1,nni)
1130 yn(j)=x(2,nni)
1131 zn(j)=x(3,nni)
1132 dxn(j)=d(1,nni)
1133 dyn(j)=d(2,nni)
1134 dzn(j)=d(3,nni)
1135 nni = ixc(3,n)
1136 xn(j+1)=x(1,nni)
1137 yn(j+1)=x(2,nni)
1138 zn(j+1)=x(3,nni)
1139 dxn(j+1)=d(1,nni)
1140 dyn(j+1)=d(2,nni)
1141 dzn(j+1)=d(3,nni)
1142 nni = ixc(4,n)
1143 xn(j+2)=x(1,nni)
1144 yn(j+2)=x(2,nni)
1145 zn(j+2)=x(3,nni)
1146 dxn(j+2)=d(1,nni)
1147 dyn(j+2)=d(2,nni)
1148 dzn(j+2)=d(3,nni)
1149 nni = ixc(5,n)
1150 xn(j+3)=x(1,nni)
1151 yn(j+3)=x(2,nni)
1152 zn(j+3)=x(3,nni)
1153 dxn(j+3)=d(1,nni)
1154 dyn(j+3)=d(2,nni)
1155 dzn(j+3)=d(3,nni)
1156 strain(1:3,i)=zero
1157 ENDDO
1158 CALL sh4_tstrain(xn,yn,zn,dxn,dyn,dzn,strain,nel)
1159 DO i=1,nel
1160 value(1:3)= strain(1:3,i)
1161 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,
1162 . VALUE)
1163 ENDDO
1164 ELSEIF (ity == 7) THEN
1165 DO i=1,nel
1166 n = i + nft
1167 nni = ixtg(2,n)
1168 j = 3*(i-1) +1
1169 xn(j)=x(1,nni)
1170 yn(j)=x(2,nni)
1171 zn(j)=x(3,nni)
1172 dxn(j)=d(1,nni)
1173 dyn(j)=d(2,nni)
1174 dzn(j)=d(3,nni)
1175 nni = ixtg(3,n)
1176 xn(j+1)=x(1,nni)
1177 yn(j+1)=x(2,nni)
1178 zn(j+1)=x(3,nni)
1179 dxn(j+1)=d(1,nni)
1180 dyn(j+1)=d(2,nni)
1181 dzn(j+1)=d(3,nni)
1182 nni = ixtg(4,n)
1183 xn(j+2)=x(1,nni)
1184 yn(j+2)=x(2,nni)
1185 zn(j+2)=x(3,nni)
1186 dxn(j+2)=d(1,nni)
1187 dyn(j+2)=d(2,nni)
1188 dzn(j+2)=d(3,nni)
1189 strain(1:3,i)=zero
1190 ENDDO
1191 CALL sh3_tstrain(xn,yn,zn,dxn,dyn,dzn,strain,nel,ihbe)
1192 DO i=1,nel
1193 value(1:3)= strain(1:3,i)
1194 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,
1195 . VALUE)
1196 ENDDO
1197 END IF
1198C--------------------------------------------------
1199C-----------------------------------------------
1200 ELSEIF (keyword == 'TENS/STRESS/TMAX') THEN
1201C---------- -------------------------------------
1202 DO i=1,nel
1203 value(1:3) = gbuf%TM_SIG1(jj(1:3) + i)
1204 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,
1205 . VALUE)
1206 ENDDO
1207C---------- -------------------------------------
1208 ELSEIF (keyword == 'TENS/STRESS/TMIN') THEN
1209C---------- -------------------------------------
1210 DO i=1,nel
1211 value(1:3) = gbuf%TM_SIG3(jj(1:3) + i)
1212 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,
1213 . VALUE)
1214 ENDDO
1215C---------- -------------------------------------
1216 ELSEIF (keyword == 'TENS/STRAIN/TMAX') THEN
1217C---------- -------------------------------------
1218 DO i=1,nel
1219 value(1:3) = gbuf%TM_STRA1(jj(1:3) + i)
1220 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,
1221 . VALUE)
1222 ENDDO
1223C---------- -------------------------------------
1224 ELSEIF (keyword == 'TENS/STRAIN/TMIN') THEN
1225C---------- -------------------------------------
1226 DO i=1,nel
1227 value(1:3) = gbuf%TM_STRA3(jj(1:3) + i)
1228 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,
1229 . VALUE)
1230 ENDDO
1231C--------------------------------------------------
1232 ELSEIF (keyword == 'TENS/BSTRESS') THEN
1233C--------------------------------------------------
1234 IF (mlw == 87) THEN
1235 imat = ixc(1,nft+1)
1236 iadbuf = ipm(7,imat)
1237 nuparam= ipm(9,imat)
1238 uparam => bufmat(iadbuf:iadbuf+nuparam)
1239 nbfunct = uparam(25)
1240 nchard = 34 + 2*nbfunct + 22
1241 chard = uparam(nchard)
1242 ELSEIF (mlw == 36) THEN
1243 imat = ixc(1,nft+1)
1244 iadbuf = ipm(7,imat)
1245 nuparam= ipm(9,imat)
1246 uparam => bufmat(iadbuf:iadbuf+nuparam)
1247 nbfunct = uparam(1)
1248 nchard = 2*nbfunct + 14
1249 chard = uparam(nchard)
1250 ENDIF
1251 IF ( ilay == -1 .AND. ipt == -1 .AND. iply == -1) THEN !global value = mean on gauss points IPs and layers
1252 IF(id == -1) THEN ! sum of all backstresses
1253 IF(mlw == 36 .AND. chard > zero) THEN
1254 DO i=1,nel
1255 DO j=1,nlay
1256 bufly => elbuf_tab(ng)%BUFLY(j)
1257 nptt = bufly%NPTT
1258 DO ir=1,nptr
1259 DO is=1,npts
1260 DO it=1,nptt
1261 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,it)
1262 DO k=1,3
1263 value(k) = value(k) + lbuf%SIGB(jj(k) + i)/npg/nptt/nlay
1264 ENDDO !K
1265 ENDDO !IT
1266 ENDDO !IS
1267 ENDDO !IR
1268 ENDDO !Jlay
1269 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE )
1270 ENDDO !I
1271 ELSEIF(mlw == 78) THEN
1272 DO i=1,nel
1273 DO j=1,nlay
1274 bufly => elbuf_tab(ng)%BUFLY(j)
1275 nptt = bufly%NPTT
1276 DO ir=1,nptr
1277 DO is=1,npts
1278 DO it=1,nptt
1279 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,it)
1280 DO k=1,3
1281 value(k) = value(k) + (lbuf%SIGA(jj(k) + i)+lbuf%SIGB(jj(k) + i))/npg/nptt/nlay
1282 ENDDO !K
1283 ENDDO !IT
1284 ENDDO !IS
1285 ENDDO !IR
1286 ENDDO !Jlay
1287 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE )
1288 ENDDO !I
1289 ELSEIF(mlw == 87 .AND. chard > zero) THEN
1290 !SIGBXX(I) = SIGB(I,1) + SIGB(I,4) + SIGB(I,7) + SIGB(I,10)
1291 !SIGBYY(I) = SIGB(I,2) + SIGB(I,5) + SIGB(I,8) + SIGB(I,11)
1292 !SIGBXY(I) = SIGB(I,3) + SIGB(I,6) + SIGB(I,9) + SIGB(I,12)
1293 DO i=1,nel
1294 DO j=1,nlay
1295 bufly => elbuf_tab(ng)%BUFLY(j)
1296 nptt = bufly%NPTT
1297 DO ir=1,nptr
1298 DO is=1,npts
1299 DO it=1,nptt
1300 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,it)
1301 DO k=1,3
1302 value(k) = value(k) + (lbuf%SIGB(jj(k) + i )
1303 . +lbuf%SIGB(jj(k+3) + i )
1304 . +lbuf%SIGB(jj(k+6) + i )
1305 . +lbuf%SIGB(jj(k+9) + i ))/npg/nptt/nlay
1306 ENDDO !K
1307 ENDDO !IT
1308 ENDDO !IS
1309 ENDDO !IR
1310 ENDDO !Jlay
1311 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE )
1312 ENDDO !I
1313 ENDIF ! MLW==
1314 ELSEIF(id > 0) THEN !!
1315 IF(mlw == 36.AND. chard > zero) THEN ! forcement ID=1 y a qu une BS
1316 DO i=1,nel
1317 DO j=1,nlay
1318 bufly => elbuf_tab(ng)%BUFLY(j)
1319 nptt = bufly%NPTT
1320 DO ir=1,nptr
1321 DO is=1,npts
1322 DO it=1,nptt
1323 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,it)
1324 DO k=1,3
1325 value(k) = value(k) + lbuf%SIGB(jj(k) + i)/npg/nptt/nlay
1326 ENDDO !K
1327 ENDDO !IT
1328 ENDDO !IS
1329 ENDDO !IR
1330 ENDDO !Jlay
1331 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE )
1332 ENDDO !I
1333 ELSEIF(mlw == 78) THEN
1334 IF(id == 1) THEN
1335 DO i=1,nel
1336 DO j=1,nlay
1337 bufly => elbuf_tab(ng)%BUFLY(j)
1338 nptt = bufly%NPTT
1339 DO ir=1,nptr
1340 DO is=1,npts
1341 DO it=1,nptt
1342 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,it)
1343 DO k=1,3
1344 value(k) = value(k) + lbuf%SIGA(jj(k) + i) /npg/nptt/nlay
1345 ENDDO !K
1346 ENDDO !IT
1347 ENDDO !IS
1348 ENDDO !IR
1349 ENDDO !Jlay
1350 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE )
1351 ENDDO !I
1352 ELSEIF(id == 2) THEN
1353 DO i=1,nel
1354 DO j=1,nlay
1355 bufly => elbuf_tab(ng)%BUFLY(j)
1356 nptt = bufly%NPTT
1357 DO ir=1,nptr
1358 DO is=1,npts
1359 DO it=1,nptt
1360 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,it)
1361 DO k=1,3
1362 value(k) = value(k) + lbuf%SIGB(jj(k) + i) /npg/nptt/nlay
1363 ENDDO !K
1364 ENDDO !IT
1365 ENDDO !IS
1366 ENDDO !IR
1367 ENDDO !Jlay
1368 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE )
1369 ENDDO !I
1370 ELSEIF(id == 3) THEN
1371 DO i=1,nel
1372 DO j=1,nlay
1373 bufly => elbuf_tab(ng)%BUFLY(j)
1374 nptt = bufly%NPTT
1375 DO ir=1,nptr
1376 DO is=1,npts
1377 DO it=1,nptt
1378 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,it)
1379 DO k=1,3
1380 value(k) = value(k) + lbuf%SIGC(jj(k) + i) /npg/nptt/nlay
1381 ENDDO !K
1382 ENDDO !IT
1383 ENDDO !IS
1384 ENDDO !IR
1385 ENDDO !Jlay
1386 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE )
1387 ENDDO !I
1388 ENDIF !ID ==
1389 ELSEIF(mlw == 87.AND. chard > zero) THEN
1390 IF(id == 1) THEN
1391 DO i=1,nel
1392 DO j=1,nlay
1393 bufly => elbuf_tab(ng)%BUFLY(j)
1394 nptt = bufly%NPTT
1395 DO ir=1,nptr
1396 DO is=1,npts
1397 DO it=1,nptt
1398 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,it)
1399 DO k=1,3
1400 value(k) = value(k) + lbuf%SIGB(jj(k) + i ) /npg/nptt/nlay
1401 ENDDO !K
1402 ENDDO !IT
1403 ENDDO !IS
1404 ENDDO !IR
1405 ENDDO !Jlay
1406 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE )
1407 ENDDO !I
1408 ELSEIF(id == 2) THEN
1409 DO i=1,nel
1410 DO j=1,nlay
1411 bufly => elbuf_tab(ng)%BUFLY(j)
1412 nptt = bufly%NPTT
1413 DO ir=1,nptr
1414 DO is=1,npts
1415 DO it=1,nptt
1416 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,it)
1417 DO k=1,3
1418 value(k) = value(k) + lbuf%SIGB(jj(k+3) + i) /npg/nptt/nlay
1419 ENDDO !K
1420 ENDDO !IT
1421 ENDDO !IS
1422 ENDDO !IR
1423 ENDDO !Jlay
1424 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE )
1425 ENDDO !I
1426 ELSEIF(id == 3) THEN
1427 DO i=1,nel
1428 DO j=1,nlay
1429 bufly => elbuf_tab(ng)%BUFLY(j)
1430 nptt = bufly%NPTT
1431 DO ir=1,nptr
1432 DO is=1,npts
1433 DO it=1,nptt
1434 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,it)
1435 DO k=1,3
1436 value(k) = value(k) + lbuf%SIGB(jj(k+6) + i) /npg/nptt/nlay
1437 ENDDO !K
1438 ENDDO !it
1439 ENDDO !IS
1440 ENDDO !IR
1441 ENDDO !Jlay
1442 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE )
1443 ENDDO !I
1444 ELSEIF( id == 4) THEN
1445 DO i=1,nel
1446 DO j=1,nlay
1447 bufly => elbuf_tab(ng)%BUFLY(j)
1448 nptt = bufly%NPTT
1449 DO ir=1,nptr
1450 DO is=1,npts
1451 DO it=1,nptt
1452 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,it)
1453 DO k=1,3
1454 value(k) = value(k) + lbuf%SIGB(jj(k+9) + i)/npg/nptt/nlay
1455 ENDDO !K
1456 ENDDO !IT
1457 ENDDO !IS
1458 ENDDO !IR
1459 ENDDO !Jlay
1460 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE )
1461 ENDDO !I
1462 ENDIF !(ID == )
1463 endif!(MLW ==
1464 ENDIF !(ID >0 )
1465C--------------------------------------------------
1466c PLY=IPLY NPT=IPT pour prop avec ply 17-51-52
1467C--------------------------------------------------
1468 ELSEIF ( iply > 0 .AND. ipt <= mpt .AND. ipt > 0 ) THEN
1469 DO j=1,nlay
1470 id_ply = 0
1471 IF (igtyp == 17 .OR. igtyp == 51) THEN
1472 id_ply = igeo(1,stack%IGEO(2+j,isubstack))
1473 ELSEIF (igtyp == 52) THEN
1474 id_ply = ply_info(1,stack%IGEO(2+j,isubstack) - numstack)
1475 ENDIF
1476 IF (id_ply == iply) THEN
1477 bufly => elbuf_tab(ng)%BUFLY(j)
1478 !--------
1479 ! LAW36
1480 !--------
1481 IF (mlw == 36 .AND.( id == -1 .OR. id == 1).AND. chard > zero) THEN
1482 DO i=1,nel
1483 DO ir=1,nptr
1484 DO is=1,npts
1485 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1486 DO k=1,3
1487 value(k) = value(k) + lbuf%SIGB(jj(k) + i)/npg
1488 ENDDO !k
1489 ENDDO !IS
1490 ENDDO !IR
1491 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1492 ENDDO ! I
1493 !--------
1494 ! LAW78
1495 !--------
1496 ELSEIF (mlw == 78) THEN
1497 IF(id == -1) THEN ! somme of all backstresses
1498 DO i=1,nel
1499 DO ir=1,nptr
1500 DO is=1,npts
1501 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1502 DO k=1,3
1503 value(k) = value(k) + (lbuf%SIGA(jj(k) + i)+lbuf%SIGB(jj(k) + i))/npg
1504 ENDDO !k
1505 ENDDO !IS
1506 ENDDO !IR
1507 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1508 ENDDO ! I
1509 ELSEIF(id ==1 ) THEN
1510 DO i=1,nel
1511 DO ir=1,nptr
1512 DO is=1,npts
1513 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1514 DO k=1,3
1515 value(k) = value(k) + lbuf%SIGA(jj(k) + i)/npg
1516 ENDDO !k
1517 ENDDO !IS
1518 ENDDO !IR
1519 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1520 ENDDO ! I
1521 ELSEIF(id ==2 ) THEN
1522 DO i=1,nel
1523 DO ir=1,nptr
1524 DO is=1,npts
1525 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1526 DO k=1,3
1527 value(k) = value(k) + lbuf%SIGB(jj(k) + i)/npg
1528 ENDDO !k
1529 ENDDO !IS
1530 ENDDO !IR
1531 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1532 ENDDO ! I
1533
1534 ELSEIF(id ==3 ) THEN
1535 DO i=1,nel
1536 DO ir=1,nptr
1537 DO is=1,npts
1538 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1539 DO k=1,3
1540 value(k) = value(k) + lbuf%SIGC(jj(k) + i)/npg
1541 ENDDO !k
1542 ENDDO !IS
1543 ENDDO !IR
1544 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1545 ENDDO ! I
1546 ENDIF !ID == -1
1547 !--------
1548 ! LAW87
1549 !--------
1550 ELSEIF( mlw == 87 .AND. chard > zero) THEN
1551 IF(id == -1) THEN ! somme of all backstresses
1552 DO i=1,nel
1553 DO ir=1,nptr
1554 DO is=1,npts
1555 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1556 DO k=1,3
1557 value(k) = value(k) + (lbuf%SIGB(jj(k) + i )
1558 . +lbuf%SIGB(jj(k+3) + i )
1559 . +lbuf%SIGB(jj(k+6) + i )
1560 . +lbuf%SIGB(jj(k+9) + i ))/npg
1561
1562 ENDDO !k
1563 ENDDO !IS
1564 ENDDO !IR
1565 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1566 ENDDO ! I
1567 ELSEIF(id ==1 ) THEN
1568 DO i=1,nel
1569 DO ir=1,nptr
1570 DO is=1,npts
1571 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1572 DO k=1,3
1573 value(k) = value(k) + lbuf%SIGB(jj(k) + i)/npg
1574 ENDDO !k
1575 ENDDO !IS
1576 ENDDO !IR
1577 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1578 ENDDO ! I
1579 ELSEIF(id ==2 ) THEN
1580 DO i=1,nel
1581 DO ir=1,nptr
1582 DO is=1,npts
1583 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1584 DO k=1,3
1585 value(k) = value(k) + lbuf%SIGB(jj(k+3) + i)/npg
1586 ENDDO !k
1587 ENDDO !IS
1588 ENDDO !IR
1589 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1590 ENDDO ! I
1591 ELSEIF(id ==3 ) THEN
1592 DO i=1,nel
1593 DO ir=1,nptr
1594 DO is=1,npts
1595 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1596 DO k=1,3
1597 value(k) = value(k) + lbuf%SIGB(jj(k+6) + i)/npg
1598 ENDDO !k
1599 ENDDO !IS
1600 ENDDO !IR
1601 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1602 ENDDO ! I
1603 ELSEIF(id ==4 ) THEN
1604 DO i=1,nel
1605 DO ir=1,nptr
1606 DO is=1,npts
1607 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1608 DO k=1,3
1609 value(k) = value(k) + lbuf%SIGB(jj(k+9) + i)/npg
1610 ENDDO !k
1611 ENDDO !IS
1612 ENDDO !IR
1613 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1614 ENDDO ! I
1615 ENDIF !ID == -1
1616 ENDIF !(MLW ==
1617 END IF !(ID_PLY == IPLY
1618 ENDDO !JLAY
1619C--------------------------------------------------
1620c PLY=NULL ILAYER=IL NPT=IPT
1621C--------------------------------------------------
1622 ELSEIF (ilay > 0 .AND. ilay <= nlay .AND. ipt <= mpt .AND. ipt > 0 ) THEN
1623 j = ilay
1624 IF(igtyp == 9) j = 1
1625 bufly => elbuf_tab(ng)%BUFLY(j)
1626 !--------
1627 ! LAW36
1628 !--------
1629 IF (mlw == 36.AND. (id==-1 . or .id==1).AND. chard > zero) THEN
1630 DO i=1,nel
1631 DO ir=1,nptr
1632 DO is=1,npts
1633 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1634 DO k=1,3
1635 value(k) = value(k) + lbuf%SIGB(jj(k) + i)/npg
1636 ENDDO !k
1637 ENDDO !IS
1638 ENDDO !IR
1639 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1640 ENDDO ! I
1641 !--------
1642 ! LAW78
1643 !--------
1644 ELSEIF (mlw == 78) THEN
1645 IF(id == -1) THEN ! somme of all backstresses
1646 DO i=1,nel
1647 DO ir=1,nptr
1648 DO is=1,npts
1649 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1650 DO k=1,3
1651 value(k) = value(k) + (lbuf%SIGA(jj(k) + i)+lbuf%SIGB(jj(k) + i))/npg
1652 ENDDO !k
1653 ENDDO !IS
1654 ENDDO !IR
1655 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1656 ENDDO ! I
1657 ELSEIF(id ==1 ) THEN
1658 DO i=1,nel
1659 DO ir=1,nptr
1660 DO is=1,npts
1661 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1662 DO k=1,3
1663 value(k) = value(k) + lbuf%SIGA(jj(k) + i)/npg
1664 ENDDO !k
1665 ENDDO !IS
1666 ENDDO !IR
1667 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1668 ENDDO ! I
1669 ELSEIF(id ==2 ) THEN
1670 DO i=1,nel
1671 DO ir=1,nptr
1672 DO is=1,npts
1673 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1674 DO k=1,3
1675 value(k) = value(k) + lbuf%SIGB(jj(k) + i)/npg
1676 ENDDO !k
1677 ENDDO !IS
1678 ENDDO !IR
1679 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1680 ENDDO ! i
1681 ELSEIF(id ==3 ) THEN
1682 DO i=1,nel
1683 DO ir=1,nptr
1684 DO is=1,npts
1685 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1686 DO k=1,3
1687 value(k) = value(k) + lbuf%SIGC(jj(k) + i)/npg
1688 ENDDO !k
1689 ENDDO !IS
1690 ENDDO !IR
1691 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1692 ENDDO ! I
1693 ENDIF !ID == -1
1694 !--------
1695 ! LAW87
1696 !--------
1697 ELSEIF( mlw == 87 .AND. chard > zero) THEN
1698 IF(id == -1) THEN ! somme of all backstresses
1699 DO i=1,nel
1700 DO ir=1,nptr
1701 DO is=1,npts
1702 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1703 DO k=1,3
1704 value(k) = value(k) + (lbuf%SIGB(jj(k) + i )
1705 . +lbuf%SIGB(jj(k+3) + i )
1706 . +lbuf%SIGB(jj(k+6) + i )
1707 . +lbuf%SIGB(jj(k+9) + i ))/npg
1708
1709 ENDDO !k
1710 ENDDO !IS
1711 ENDDO !IR
1712 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1713 ENDDO ! I
1714 ELSEIF(id ==1 ) THEN
1715 DO i=1,nel
1716 DO ir=1,nptr
1717 DO is=1,npts
1718 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1719 DO k=1,3
1720 value(k) = value(k) + lbuf%SIGB(jj(k) + i)/npg
1721 ENDDO !k
1722 ENDDO !IS
1723 ENDDO !IR
1724 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1725 ENDDO ! I
1726 ELSEIF(id ==2 ) THEN
1727 DO i=1,nel
1728 DO ir=1,nptr
1729 DO is=1,npts
1730 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1731 DO k=1,3
1732 value(k) = value(k) + lbuf%SIGB(jj(k+3) + i)/npg
1733 ENDDO !k
1734 ENDDO !IS
1735 ENDDO !IR
1736 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1737 ENDDO ! I
1738 ELSEIF(id ==3 ) THEN
1739 DO i=1,nel
1740 DO ir=1,nptr
1741 DO is=1,npts
1742 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1743 DO k=1,3
1744 VALUE(k) = value(k) + lbuf%SIGB(jj(k+6) + i)/npg
1745 ENDDO !k
1746 ENDDO !IS
1747 ENDDO !IR
1748 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1749 ENDDO ! I
1750 ELSEIF(id ==4 ) THEN
1751 DO i=1,nel
1752 DO ir=1,nptr
1753 DO is=1,npts
1754 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1755 DO k=1,3
1756 value(k) = value(k) + lbuf%SIGB(jj(k+9) + i)/npg
1757 ENDDO !k
1758 ENDDO !is
1759 ENDDO !IR
1760 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1761 ENDDO ! I
1762 ENDIF !ID == -1
1763 ENDIF !(MLW ==
1764C--------------------------------------------------
1765c PLY=NULL ILAYER=IL NPT=NULL ! prop 9-10-11 have layers and not PLY
1766C--------------------------------------------------
1767 ELSEIF (iply == -1 .AND. ilay <= nlay .AND. ilay > 0 .AND. ipt == -1 ) THEN
1768 IF (igtyp == 9 .OR.igtyp == 10 .OR. igtyp == 11 ) THEN
1769 ! output in orthotopic frame / elementary
1770 j = ilay
1771 IF(igtyp == 9) j = 1
1772 bufly => elbuf_tab(ng)%BUFLY(j)
1773 nptt = bufly%NPTT
1774 !--------
1775 ! LAW36
1776 !--------
1777 IF (mlw == 36.AND. (id==-1 .OR. id==1) .AND. chard > zero) THEN ! only one Bstress
1778 DO i=1,nel
1779 DO ir=1,nptr
1780 DO is=1,npts
1781 DO it=1,nptt
1782 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,it)
1783 DO k=1,3
1784 value(k) = value(k) + lbuf%SIGB(jj(k) + i)/npg/nptt
1785 ENDDO !K
1786 enddo!IT
1787 ENDDO !IS
1788 enddo!IR
1789 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1790 value(1:3) = zero
1791 ENDDO !I
1792 !--------
1793 ! LAW78
1794 !--------
1795 ELSEIF (mlw == 78) THEN
1796 IF(id == -1) THEN ! somme of all backstresses
1797 DO i=1,nel
1798 DO ir=1,nptr
1799 DO is=1,npts
1800 DO it=1,nptt
1801 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,it)
1802 DO k=1,3
1803 value(k) = value(k) + (lbuf%SIGA(jj(k) + i)+lbuf%SIGB(jj(k) + i))/npg/nptt
1804 ENDDO !k
1805 enddo!IT
1806 ENDDO !IS
1807 ENDDO !IR
1808 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1809 ENDDO ! I
1810 ELSEIF(id ==1 ) THEN
1811 DO i=1,nel
1812 DO ir=1,nptr
1813 DO is=1,npts
1814 DO it=1,nptt
1815 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,it)
1816 DO k=1,3
1817 value(k) = value(k) + lbuf%SIGA(jj(k) + i)/npg/nptt
1818 ENDDO !k
1819 enddo!IT
1820 ENDDO !is
1821 ENDDO !IR
1822 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1823 ENDDO ! I
1824 ELSEIF(id ==2 ) THEN
1825 DO i=1,nel
1826 DO ir=1,nptr
1827 DO is=1,npts
1828 DO it=1,nptt
1829 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,it)
1830 DO k=1,3
1831 value(k) = value(k) + lbuf%SIGB(jj(k) + i)/npg/nptt
1832 ENDDO !k
1833 enddo!IT
1834 ENDDO !IS
1835 ENDDO !IR
1836 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1837 ENDDO ! I
1838 ELSEIF(id ==3 ) THEN
1839 DO i=1,nel
1840 DO ir=1,nptr
1841 DO is=1,npts
1842 DO it=1,nptt
1843 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,it)
1844 DO k=1,3
1845 value(k) = value(k) + lbuf%SIGC(jj(k) + i)/npg/nptt
1846 ENDDO !k
1847 enddo!IT
1848 ENDDO !IS
1849 ENDDO !IR
1850 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1851 ENDDO ! I
1852 ENDIF !ID == -1
1853 !--------
1854 ! LAW87
1855 !--------
1856 ELSEIF( mlw == 87 .AND. chard > zero) THEN
1857 IF(id == -1) THEN ! somme of all backstresses
1858 DO i=1,nel
1859 DO ir=1,nptr
1860 DO is=1,npts
1861 DO it=1,nptt
1862 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,it)
1863 DO k=1,3
1864 value(k) = value(k) + (lbuf%SIGB(jj(k) + i )
1865 . +lbuf%SIGB(jj(k+3) + i )
1866 . +lbuf%SIGB(jj(k+6) + i )
1867 . +lbuf%SIGB(jj(k+9) + i ))/npg/nptt
1868
1869 ENDDO !k
1870 enddo!IT
1871 ENDDO !IS
1872 ENDDO !IR
1873 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1874 ENDDO ! I
1875 ELSEIF(id ==1 ) THEN
1876 DO i=1,nel
1877 DO ir=1,nptr
1878 DO is=1,npts
1879 DO it=1,nptt
1880 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,it)
1881 DO k=1,3
1882 value(k) = value(k) + lbuf%SIGB(jj(k) + i)/npg/nptt
1883 ENDDO !k
1884 enddo!IT
1885 ENDDO !IS
1886 ENDDO !IR
1887 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1888 ENDDO ! I
1889 ELSEIF(id ==2 ) THEN
1890 DO i=1,nel
1891 DO ir=1,nptr
1892 DO is=1,npts
1893 DO it=1,nptt
1894 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,it)
1895 DO k=1,3
1896 value(k) = value(k) + lbuf%SIGB(jj(k+3) + i )/npg/nptt
1897 ENDDO !k
1898 enddo!IT
1899 ENDDO !IS
1900 ENDDO !IR
1901 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1902 ENDDO ! I
1903 ELSEIF(id ==3 ) THEN
1904 DO i=1,nel
1905 DO ir=1,nptr
1906 DO is=1,npts
1907 DO it=1,nptt
1908 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,it)
1909 DO k=1,3
1910 value(k) = value(k) + lbuf%SIGB(jj(k+6) + i )/npg/nptt
1911 ENDDO !k
1912 enddo!IT
1913 ENDDO !IS
1914 ENDDO !IR
1915 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1916 ENDDO ! I
1917 ELSEIF(id ==4 ) THEN
1918 DO i=1,nel
1919 DO ir=1,nptr
1920 DO is=1,npts
1921 DO it=1,nptt
1922 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,it)
1923 DO k=1,3
1924 value(k) = value(k) + lbuf%SIGB(jj(k+9) + i )/npg
1925 ENDDO !k
1926 enddo!IT
1927 ENDDO !IS
1928 ENDDO !IR
1929 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1930 ENDDO ! I
1931 ENDIF !ID == -1
1932 ENDIF !(MLW ==
1933 ENDIF !IGTYP ==
1934C---------------------------------------------------------------------------
1935c ILAYER=-1 NPT=IPT IPLY=-1 ID=-1
1936C---------------------------------------------------------------------------
1937 ELSE IF(ilay == -1 .AND. ipt > 0 .AND. ipt<=mpt .AND. iply == -1 ) THEN ! output for each layer/PLY as if LAYER==ALL
1938 DO j=1,nlay
1939 bufly => elbuf_tab(ng)%BUFLY(j)
1940 nptt = bufly%NPTT
1941 IF (ipt <= nptt ) THEN
1942 !--------
1943 ! LAW36
1944 !--------
1945 IF (mlw == 36.AND. (id==-1 .OR. id==1) .AND. chard > zero) THEN
1946 DO i=1,nel
1947 DO ir=1,nptr
1948 DO is=1,npts
1949 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1950 DO k=1,3
1951 value(k) = value(k) + lbuf%SIGB(jj(k) + i)/npg
1952 ENDDO !k
1953 ENDDO !IS
1954 ENDDO !IR
1955 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1956 ENDDO ! I
1957 !--------
1958 ! LAW78
1959 !--------
1960 ELSEIF (mlw == 78) THEN
1961 IF(id == -1) THEN ! somme of all backstresses
1962 DO i=1,nel
1963 DO ir=1,nptr
1964 DO is=1,npts
1965 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1966 DO k=1,3
1967 value(k) = value(k) + (lbuf%SIGA(jj(k) + i)+lbuf%SIGB(jj(k) + i))/npg
1968 ENDDO !k
1969 ENDDO !IS
1970 ENDDO !IR
1971 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1972 ENDDO ! I
1973 ELSEIF(id ==1 ) THEN
1974 DO i=1,nel
1975 DO ir=1,nptr
1976 DO is=1,npts
1977 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1978 DO k=1,3
1979 value(k) = value(k) + lbuf%SIGA(jj(k) + i)/npg
1980 ENDDO !k
1981 ENDDO !IS
1982 ENDDO !IR
1983 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1984 ENDDO ! I
1985 ELSEIF(id ==2 ) THEN
1986 DO i=1,nel
1987 DO ir=1,nptr
1988 DO is=1,npts
1989 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
1990 DO k=1,3
1991 value(k) = VALUE(k) + lbuf%SIGB(jj(k) + i)/npg
1992 ENDDO !k
1993 ENDDO !IS
1994 ENDDO !IR
1995 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
1996 ENDDO ! I
1997 ELSEIF(id ==3 ) THEN
1998 DO i=1,nel
1999 DO ir=1,nptr
2000 DO is=1,npts
2001 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
2002 DO k=1,3
2003 value(k) = value(k) + lbuf%SIGC(jj(k) + i)/npg
2004 ENDDO !k
2005 ENDDO !IS
2006 ENDDO !IR
2007 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
2008 ENDDO ! I
2009 ENDIF !ID == -1
2010 !--------
2011 ! LAW87
2012 !--------
2013 ELSEIF( mlw == 87 .AND. chard > zero) THEN
2014 IF(id == -1) THEN ! somme of all backstresses
2015 DO i=1,nel
2016 DO ir=1,nptr
2017 DO is=1,npts
2018 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
2019 DO k=1,3
2020 value(k) = value(k) + (lbuf%SIGB(jj(k) + i )
2021 . +lbuf%SIGB(jj(k+3) + i )
2022 . +lbuf%SIGB(jj(k+6) + i )
2023 . +lbuf%SIGB(jj(k+9) + i ))/npg
2024
2025 ENDDO !k
2026 ENDDO !IS
2027 ENDDO !IR
2028 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
2029 ENDDO ! I
2030 ELSEIF(id ==1 ) THEN
2031 DO i=1,nel
2032 DO ir=1,nptr
2033 DO is=1,npts
2034 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
2035 DO k=1,3
2036 value(k) = value(k) + lbuf%SIGB(jj(k) + i )/npg
2037 ENDDO !k
2038 ENDDO !IS
2039 ENDDO !IR
2040 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
2041 ENDDO ! I
2042 ELSEIF(id ==2 ) THEN
2043 DO i=1,nel
2044 DO ir=1,nptr
2045 DO is=1,npts
2046 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
2047 DO k=1,3
2048 value(k) = value(k) + lbuf%SIGB(jj(k+3) + i )/npg
2049 ENDDO !k
2050 ENDDO !IS
2051 ENDDO !IR
2052 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
2053 ENDDO ! I
2054 ELSEIF(id ==3 ) THEN
2055 DO i=1,nel
2056 DO ir=1,nptr
2057 DO is=1,npts
2058 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
2059 DO k=1,3
2060 value(k) = value(k) + lbuf%SIGB(jj(k+6) + i )/npg
2061 ENDDO !k
2062 ENDDO !IS
2063 ENDDO !IR
2064 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
2065 ENDDO ! I
2066 ELSEIF(id ==4 ) THEN
2067 DO i=1,nel
2068 DO ir=1,nptr
2069 DO is=1,npts
2070 lbuf => elbuf_tab(ng)%BUFLY(j)%LBUF(ir,is,ipt)
2071 DO k=1,3
2072 value(k) = value(k) + lbuf%SIGB(jj(k+9) + i )/npg
2073 ENDDO !k
2074 ENDDO !IS
2075 ENDDO !IR
2076 CALL h3d_write_sh_tensor(iok_part,iselect,is_written_shell,shell_tensor,i,offset,nft,VALUE)
2077 ENDDO ! I
2078 ENDIF !ID == -1
2079 ENDIF !(MLW ==
2080 endif! (IPT <= NPTT )
2081 ENDDO !JLAY
2082
2083 END IF
2084c ........
2085C---------------------------------------------------------------------------
2086c ELSEIF (KEYWORD == 'NEWKEY') THEN ! New Output Example
2087C---------------------------------------------------------------------------
2088c ILAYER=NULL NPT=NULL
2089c IF ( ILAY == -1 .AND. IPT == -1 .AND. IPLY == -1) THEN
2090c DO I=1,NEL
2091c VALUE(I) =
2092c ENDDO
2093c PLY=IPLY NPT=IPT
2094c ELSEIF ( IPLY > 0 .AND. IPT <= MPT .AND. IPT > 0 ) THEN
2095c IF (IGTYP == 17 .OR. IGTYP == 51 .OR. IGTYP == 52) THEN
2096c
2097c ENDIF
2098c
2099c PLY=NULL ILAYER=ILAY NPT=IPT
2100c ELSEIF (IPLY == -1 .AND. ILAY <= NLAY .AND. ILAY > 0 .AND. IPT <= MPT .AND. IPT > 0 ) THEN
2101c IF (IGTYP == 51 .OR. IGTYP == 52) THEN
2102c
2103c ENDIF
2104c PLY=NULL ILAYER=IL NPT=NULL
2105c ELSEIF (IPLY == -1 .AND. ILAY <= NLAY .AND. ILAY > 0 .AND. IPT == -1 ) THEN
2106c IF (IGTYP == 10 .OR. IGTYP == 11 .OR. IGTYP == 16 .OR. IGTYP == 17) THEN
2107c
2108c ELSEIF (IGTYP == 51 .OR. IGTYP == 52) THEN
2109c
2110c ENDIF
2111c PLY=NULL ILAYER=NULL NPT=IPT
2112c ELSEIF ( IPT <= MPT .AND. IPT > 0) THEN
2113c IF (IGTYP == 1 .OR. IGTYP == 9) THEN
2114c
2115c ENDIF
2116c ENDIF
2117 ENDIF
2118
2119C-----------------------------------------------
2120C RNUR
2121C-----------------------------------------------
2122 ELSEIF (ity == 50) THEN
2123c DO I=1,NEL
2124c N = I + NFT
2125c TENS(1,EL2FA(NN9+N)) = ZERO
2126c TENS(2,EL2FA(NN9+N)) = ZERO
2127c TENS(3,EL2FA(NN9+N)) = ZERO
2128c ENDDO
2129C-----------------------------------------------
2130 ELSE
2131 ENDIF ! IF (MLW /= 13)
2132 ENDIF ! IF(ITY == 2)
2133 490 CONTINUE
2134 500 CONTINUE
2135C-----------------------------------------------
2136C
2137 RETURN
subroutine sh3_tstrain(xn, yn, zn, dx, dy, dz, strain, nel, ihbe)
subroutine sh4_tstrain(xn, yn, zn, dx, dy, dz, strain, nel)
subroutine h3d_write_sh_tensor(iok_part, iselect, is_written, tensor, i, offset, nft, value)
subroutine h3d_write_sh_tensor_array(iok_part, iselect, nel, offset, nft, is_written, tensor, value)
subroutine layini(elbuf_str, jft, jlt, geo, igeo, mat, pid, thkly, matly, posly, igtyp, ixfem, ixlay, nlay, npt, isubstack, stack, drape, nft, thk, nel, ratio_thkly, indx_drape, sedrape, numel_drape)
Definition layini.F:47
initmumps id
integer numeltg_drape
Definition drape_mod.F:92
integer scdrape
Definition drape_mod.F:92
integer stdrape
Definition drape_mod.F:92
integer numelc_drape
Definition drape_mod.F:92
integer, parameter ncharline100
integer, dimension(:,:), allocatable ply_info
Definition stack_mod.F:133
subroutine nbfunct(nfunct, ntable, npts, lsubmodel)
Definition nbfunc.F:37
subroutine roto_tens2d(nel, sig, dir)
Definition roto_tens2d.F:29
subroutine roto_tens2d_aniso(nel, tens, dir_a, dir_b)
subroutine uroto_tens2d(nel, sig, dir)
subroutine uroto_tens2d_aniso(nel, tens, dir_a, dir_b)

◆ sh3_tstrain()

subroutine sh3_tstrain ( xn,
yn,
zn,
dx,
dy,
dz,
strain,
integer nel,
integer ihbe )

Definition at line 2439 of file h3d_shell_tensor.F.

2440C---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8
2441#include "implicit_f.inc"
2442C---------+---------+---+---+--------------------------------------------
2443C VAR | SIZE |TYP| RW| DEFINITION
2444C---------+---------+---+---+--------------------------------------------
2445C NEL | 1 | I | R | SIZE OF THE ELEMENT GROUP NEL
2446C XN | 3*NEL | R | R | X-coordinate ARRAY (3n tria)
2447C YN | 3*NEL | R | R | Y-coordinate ARRAY (3n tria)
2448C ZN | 3*NEL | R | R | Z-coordinate ARRAY (3n tria)
2449C DX | 3*NEL | R | R | X-Displ ARRAY (3n tria)
2450C DY | 3*NEL | R | R | Y-Displ ARRAY (3n tria)
2451C DZ | 3*NEL | R | R | D-Displ ARRAY (3n tria)
2452C STRAIN | 3*NEL | R | W | STRAIN ARRAY
2453C---------+---------+---+---+--------------------------------------------
2454C-----------------------------------------------
2455C D U M M Y A R G U M E N T S
2456C-----------------------------------------------
2457 INTEGER NEL,IHBE
2458 my_real
2459 . xn(3,nel) , yn(3,nel) , zn(3,nel),
2460 . dx(3,nel) , dy(3,nel) , dz(3,nel),strain(3,*)
2461C-----------------------------------------------
2462C L O C A L V A R I A B L E S
2463C-----------------------------------------------
2464 INTEGER I,J,NNOD
2465 parameter(nnod = 3)
2466 my_real
2467 . x0n(nnod,nel) , y0n(nnod,nel) , z0n(nnod,nel),
2468 . deta1(nel),xx,yy,zz,
2469 . xl2(nel),xl3(nel),yl2(nel),yl3(nel),
2470 . ux1,uy1,ux2,ux3,uy2,uy3,
2471 . px2(nel),px3(nel),py2(nel),py3(nel),
2472 . ux21(nel),ux31(nel),uy21(nel),uy31(nel),
2473 . x0l2(nel),x0l3(nel),y0l2(nel),
2474 . y0l3(nel),area(nel),area_i(nel),fxx,fyy,fxy,fyx,f(2,2,nel)
2475C----------------------------------------------
2476C------Compute coordinates in elementary local sys: actual configuration first
2477 CALL c3sysg2l(xn,yn,zn,xl2,yl2,xl3,yl3,area,nel)
2478C------initial configuration :
2479 DO i=1,nel
2480 x0n(1:nnod,i) = xn(1:nnod,i)-dx(1:nnod,i)
2481 y0n(1:nnod,i) = yn(1:nnod,i)-dy(1:nnod,i)
2482 z0n(1:nnod,i) = zn(1:nnod,i)-dz(1:nnod,i)
2483 ENDDO
2484 CALL c3sysg2l(x0n,y0n,z0n,x0l2,y0l2,x0l3,y0l3,area,nel)
2485C-----------[B0]-----------------
2486 DO i=1,nel
2487 area_i(i)=0.5/area(i)
2488 px2(i)= y0l3(i)*area_i(i)
2489 py2(i)=-x0l3(i)*area_i(i)
2490 px3(i)=-y0l2(i)*area_i(i)
2491 py3(i)= x0l2(i)*area_i(i)
2492 ENDDO
2493C------ objective disp (free rigide rotation)
2494 ux1= zero
2495 uy1= zero
2496 DO i=1,nel
2497 ux2 = xl2(i)-x0l2(i)
2498 uy2 = yl2(i)-y0l2(i)
2499 ux3 = xl3(i)-x0l3(i)
2500 uy3 = yl3(i)-y0l3(i)
2501 ux21(i)=ux2-ux1
2502 ux31(i)=ux3-ux1
2503 uy21(i)=uy2-uy1
2504 uy31(i)=uy3-uy1
2505 ENDDO
2506C---------------
2507C MEMBRANE [F]-1 = [B0]{d}, [e] ({d} has been projected to free rigide rotation)
2508 DO i=1,nel
2509 f(1,1,i) = px2(i)*ux21(i)+px3(i)*ux31(i)
2510 f(2,2,i) = py2(i)*uy21(i)+py3(i)*uy31(i)
2511 f(2,1,i) = px2(i)*uy21(i)+px3(i)*uy31(i)
2512 f(1,2,i) = py2(i)*ux21(i)+py3(i)*ux31(i)
2513 strain(1,i) = f(1,1,i)
2514 strain(2,i) = f(2,2,i)
2515 strain(3,i) = 0.5*(f(1,2,i) + f(2,1,i))
2516 ENDDO
2517C--- local sys for 3N isn't really material ----
2518 IF (ihbe==3) CALL u_from_f2(f,strain,nel )
2519C
2520 RETURN
subroutine u_from_f2(f, strain, nel)
subroutine c3sysg2l(xn, yn, zn, xl2, yl2, xl3, yl3, area, nel)

◆ sh4_tstrain()

subroutine sh4_tstrain ( xn,
yn,
zn,
dx,
dy,
dz,
strain,
integer nel )

Definition at line 2146 of file h3d_shell_tensor.F.

2147C---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8
2148#include "implicit_f.inc"
2149C---------+---------+---+---+--------------------------------------------
2150C VAR | SIZE |TYP| RW| DEFINITION
2151C---------+---------+---+---+--------------------------------------------
2152C NEL | 1 | I | R | SIZE OF THE ELEMENT GROUP NEL
2153C XN | 4*NEL | R | R | X-coordinate ARRAY (4n quad)
2154C YN | 4*NEL | R | R | Y-coordinate ARRAY (4n quad)
2155C ZN | 4*NEL | R | R | Z-coordinate ARRAY (4n quad)
2156C DX | 4*NEL | R | R | X-Displ ARRAY (4n quad)
2157C DY | 4*NEL | R | R | Y-Displ ARRAY (4n quad)
2158C DZ | 4*NEL | R | R | D-Displ ARRAY (4n quad)
2159C STRAIN | 3*NEL | R | W | STRAIN ARRAY
2160C---------+---------+---+---+--------------------------------------------
2161C-----------------------------------------------
2162C D U M M Y A R G U M E N T S
2163C-----------------------------------------------
2164 INTEGER NEL
2165 my_real
2166 . xn(4,nel) , yn(4,nel) , zn(4,nel),
2167 . dx(4,nel) , dy(4,nel) , dz(4,nel),strain(3,*)
2168C-----------------------------------------------
2169C L O C A L V A R I A B L E S
2170C-----------------------------------------------
2171 INTEGER I,J,NNOD
2172 parameter(nnod = 4)
2173 my_real
2174 . x0n(nnod,nel) , y0n(nnod,nel) , z0n(nnod,nel),
2175 . lxyz0(3),deta1(nel),corel(2,nnod),xx,yy,zz,
2176 . xl2(nel),xl3(nel),xl4(nel),yl2(nel),
2177 . yl3(nel),yl4(nel),zl1(nel),zl(nel),
2178 . x13,x24,y13,y24,ux1,uy1,ux2,ux3,ux4,uy2,uy3,uy4,
2179 . px1(nel),px2(nel),py1(nel),py2(nel),
2180 . ux13(nel),ux24(nel),uy13(nel),uy24(nel),
2181 . x0l2(nel),x0l3(nel),x0l4(nel),y0l2(nel),
2182 . y0l3(nel),y0l4(nel),area(nel),area_i(nel),fxx,fyy,fxy,fyx
2183C----------------------------------------------
2184C------Compute coordinates in elementary local sys: actual configuration first
2185 CALL c4sysg2l(xn,yn,zn,xl2,yl2,xl3,yl3,xl4,yl4,zl1,area,nel)
2186C------initial configuration :
2187 DO i=1,nel
2188 x0n(1:nnod,i) = xn(1:nnod,i)-dx(1:nnod,i)
2189 y0n(1:nnod,i) = yn(1:nnod,i)-dy(1:nnod,i)
2190 z0n(1:nnod,i) = zn(1:nnod,i)-dz(1:nnod,i)
2191 ENDDO
2192 CALL c4sysg2l(x0n,y0n,z0n,x0l2,y0l2,x0l3,y0l3,x0l4,y0l4,zl,area,nel)
2193C-------[B0]---remove the origine to the center--------------
2194 DO i=1,nel
2195C-----------EM20=1.0E-20,FOURTH=0.25,HALF=0.5,ZERO=0.--------------
2196 area_i(i) = one/max(em20,area(i))
2197 lxyz0(1)=fourth*(x0l2(i)+x0l3(i)+x0l4(i))
2198 lxyz0(2)=fourth*(y0l2(i)+y0l3(i)+y0l4(i))
2199 corel(1,1)=-lxyz0(1)
2200 corel(1,2)=x0l2(i)-lxyz0(1)
2201 corel(1,3)=x0l3(i)-lxyz0(1)
2202 corel(1,4)=x0l4(i)-lxyz0(1)
2203 corel(2,1)=-lxyz0(2)
2204 corel(2,2)=y0l2(i)-lxyz0(2)
2205 corel(2,3)=y0l3(i)-lxyz0(2)
2206 corel(2,4)=y0l4(i)-lxyz0(2)
2207C----
2208 x13 =(corel(1,1)-corel(1,3))*half
2209 x24 =(corel(1,2)-corel(1,4))*half
2210 y13 =(corel(2,1)-corel(2,3))*half
2211 y24 =(corel(2,2)-corel(2,4))*half
2212 py2(i) =x13*area_i(i)
2213 py1(i) =-x24*area_i(i)
2214 px2(i) =-y13*area_i(i)
2215 px1(i) =y24*area_i(i)
2216 ENDDO
2217C------ objective disp (or projected disp to free rigide rotation)
2218 ux1= zero
2219 uy1= zero
2220 DO i=1,nel
2221 ux2 = xl2(i)-x0l2(i)
2222 uy2 = yl2(i)-y0l2(i)
2223 ux3 = xl3(i)-x0l3(i)
2224 uy3 = yl3(i)-y0l3(i)
2225 ux4 = xl4(i)-x0l4(i)
2226 uy4 = yl4(i)-y0l4(i)
2227 ux13(i)=ux1-ux3
2228 ux24(i)=ux2-ux4
2229 uy13(i)=uy1-uy3
2230 uy24(i)=uy2-uy4
2231 ENDDO
2232C---------------
2233C MEMBRANE [F]-1 = [B0]{d}, [e] ({d} has been projected to free rigide rotation)
2234C---------------
2235 DO i=1,nel
2236 fxx = px1(i)*ux13(i)+px2(i)*ux24(i)
2237 fyy = py1(i)*uy13(i)+py2(i)*uy24(i)
2238 fyx = px1(i)*uy13(i)+px2(i)*uy24(i)
2239 fxy = py1(i)*ux13(i)+py2(i)*ux24(i)
2240 strain(1,i) = fxx
2241 strain(2,i) = fyy
2242 strain(3,i) = 0.5*(fxy+fyx)
2243 ENDDO
2244C
2245 RETURN
subroutine c4sysg2l(xn, yn, zn, xl2, yl2, xl3, yl3, xl4, yl4, zl1, area, nel)

◆ u_from_f2()

subroutine u_from_f2 ( f,
strain,
integer nel )

Definition at line 2594 of file h3d_shell_tensor.F.

2595C-----------------------------------------------
2596C I m p l i c i t T y p e s
2597C-----------------------------------------------
2598#include "implicit_f.inc"
2599C-----------------------------------------------
2600C D u m m y A r g u m e n t s
2601C-----------------------------------------------
2602 INTEGER NEL
2603 my_real
2604 . f(2,2,nel), strain(3,*)
2605C-----------------------------------------------
2606C L o c a l V a r i a b l e s
2607C-----------------------------------------------
2608 INTEGER I
2609 DOUBLE PRECISION
2610 . FMAT(2,2),UM(2,2),IC,I2C,I3C,IU,I2U,I3U,
2611 . C11,C21,C31,C12,C22,C32,C13,C23,C33,DETJ0,DETM1,
2612 . CC11,CC21,CC31,CC12,CC22,CC32,CC13,CC23,CC33,
2613 . A,B,ZZ,A1,A2,A3,A4,ALPHA,C_A,S_A,C_A2,S_A2
2614C-----------------------------------------------
2615 DO i=1,nel
2616 fmat(1,1)= f(1,1,i)+1.0
2617 fmat(2,2)= f(2,2,i)+1.0
2618 fmat(1,2)= f(1,2,i)
2619 fmat(2,1)= f(2,1,i)
2620 c11 = fmat(1,1)*fmat(1,1)+fmat(2,1)*fmat(2,1)
2621 c12 = fmat(1,1)*fmat(1,2)+fmat(2,1)*fmat(2,2)
2622 c22 = fmat(1,2)*fmat(1,2)+fmat(2,2)*fmat(2,2)
2623 cc12 = 0.5*(c11-c22)
2624 ic = 0.5*(c11+c22)
2625 b = sqrt(cc12*cc12+c12*c12)
2626 cc11 = sqrt(ic + b)
2627 cc22 = sqrt(ic - b)
2628 IF (abs(cc12)<em20) THEN
2629 c_a = 0.0
2630 s_a = 1.0
2631 ELSE
2632 alpha =0.5*(atan(c12/cc12))
2633 c_a = cos(alpha)
2634 s_a = sin(alpha)
2635 END IF
2636 c_a2 = c_a*c_a
2637 s_a2 = s_a*s_a
2638 um(1,1) = cc11*c_a2+cc22*s_a2
2639 um(2,2) = cc11*s_a2+cc22*c_a2
2640 um(1,2) = (cc11-cc22)*s_a*c_a
2641 strain(1,i) = um(1,1)-one
2642 strain(2,i) = um(2,2)-one
2643 strain(3,i) = um(1,2)
2644 END DO
2645C
2646 RETURN
#define alpha
Definition eval.h:35