43
44
45
46 USE elbufdef_mod
47 use element_mod , only : nixc
48
49#include "implicit_f.inc"
50
51
52
53#include "mvsiz_p.inc"
54#include "param_c.inc"
55#include "impl1_c.inc"
56#include "comlock.inc"
57#include "units_c.inc"
58#include "scr17_c.inc"
59
60
61
62 INTEGER IXC(NIXC,*),JFT,JLT,IREP,NPT,ISMSTR,PID(*),MAT(*),NGL(*),NEL
63 INTEGER NPLAT,NLAY,IPLAT(*)
65 . pm(npropm,*),geo(npropg,*), x(3,*),
66 . mx23(*),my13(*),my23(*),my34(*),
67 . x13(*),x24(*),y13(*),y24(*),mx13(*),
68 . vq(mvsiz,3,3),
area(*),z1(*),mx34(*),vqn(mvsiz,3,4),area_i(*)
69
71 . dir_a(nel,*),dir_b(nel,*),offg(*),off(*),
72 . corelv(mvsiz,2,4),thk(*)
73 double precision
74 . smstr(*)
75 TYPE(ELBUF_STRUCT_) :: ELBUF_STR
76
77
78
79 INTEGER NNOD,I,J,K,L,M,II(6),SPLAT,JPLAT(MVSIZ),MAT_1
80 parameter(nnod = 4)
82 . lxyz0(3),deta1(mvsiz),rx(mvsiz), ry(mvsiz), rz(mvsiz),
83 . sx(mvsiz),sy(mvsiz),r11(mvsiz),r12(mvsiz),r13(mvsiz),
84 . r21(mvsiz),r22(mvsiz),r23(mvsiz),r31(mvsiz),r32(mvsiz),
85 . r33(mvsiz),xl2(mvsiz),xl3(mvsiz),xl4(mvsiz),yl2(mvsiz),
86 . yl3(mvsiz),yl4(mvsiz),ssz(mvsiz),
87 . l13(mvsiz),l24(mvsiz),ll(mvsiz),
88 . tol,x13_2(mvsiz),y13_2(mvsiz),x24_2(mvsiz),y24_2(mvsiz),
89 . z2(mvsiz),a_4,sz,sz1,sz2,sl,s1,
90 .
91 . j0,j1,j2
93 . xx1,xx2,xx3,xx4,yy1,yy2,yy3,yy4,zz1,zz2,zz3,zz4
94
95 DO i=1,6
96 ii(i) = nel*(i-1)
97 ENDDO
98
99 tol=em8
100 mat_1 = ixc(1,jft)
101 DO i=jft,jlt
102 mat(i) = mat_1
103 pid(i) = ixc(6,i)
104 ngl(i) = ixc(7,i)
105 ENDDO
106
107
108
109 DO i=jft,jlt
110 rx(i)=x(1,ixc(3,i))+x(1,ixc(4,i))-x(1,ixc(2,i))-x(1,ixc(5,i))
111 sx(i)=x(1,ixc(4,i))+x(1,ixc(5,i))-x(1,ixc(2,i))-x(1,ixc(3,i))
112 ry(i)=x(2,ixc(3,i))+x(2,ixc(4,i))-x(2,ixc(2,i))-x(2,ixc(5,i))
113 sy(i)=x(2,ixc(4,i))+x(2,ixc(5,i))-x(2,ixc(2,i))-x(2,ixc(3,i))
114 rz(i)=x(3,ixc(3,i))+x(3,ixc(4,i))-x(3,ixc(2,i))-x(3,ixc(5,i))
115 ssz(i)=x(3,ixc(4,i))+x(3,ixc(5,i))-x(3,ixc(2,i))-x(3,ixc(3,i))
116 ENDDO
117 k = 0
119 . rx, ry, rz,
120 . sx, sy, ssz,
121 . r11,r12,r13,r21,r22,r23,r31,r32,r33,deta1,offg )
122 DO i=jft,jlt
123 area(i)=fourth*deta1(i)
124 area_i(i)=one/
area(i)
125 vq(i,1,1)=r11(i)
126 vq(i,2,1)=r21(i)
127 vq(i,3,1)=r31(i)
128 vq(i,1,2)=r12(i)
129 vq(i,2,2)=r22(i)
130 vq(i,3,2)=r32(i)
131 vq(i,1,3)=r13(i)
132 vq(i,2,3)=r23(i)
133 vq(i,3,3)=r33(i)
134 ENDDO
135
136
137
138 DO i=jft,jlt
139 j=ixc(2,i)
140 k=ixc(3,i)
141 l=ixc(4,i)
142 m=ixc(5,i)
143 lxyz0(1)=fourth*(x(1,l)+x(1,m)+x(1,j)+x(1,k))
144 lxyz0(2)=fourth*(x(2,l)+x(2,m)+x(2,j)+x(2,k))
145 lxyz0(3)=fourth*(x(3,l)+x(3,m)+x(3,j)+x(3,k))
146
147 xx1=x(1,k)-x(1,j)
148 yy1=x(2,k)-x(2,j)
149 zz1=x(3,k)-x(3,j)
150
151 xl2(i)=r11(i)*xx1+r21(i)*yy1+r31(i)*zz1
152 yl2(i)=r12(i)*xx1+r22(i)*yy1+r32(i)*zz1
153
154 xx2=x(1,j)-lxyz0(1)
155 yy2=x(2,j)-lxyz0(2)
156 zz2=x(3,j)-lxyz0(3)
157 z1(i)=r13(i)*xx2+r23(i)*yy2+r33(i)*zz2
158
159 xx3=x(1,l)-x(1,j)
160 yy3=x(2,l)-x(2,j)
161 zz3=x(3,l)-x(3,j)
162 xl3(i)=r11(i)*xx3+r21(i)*yy3+r31(i)*zz3
163 yl3(i)=r12(i)*xx3+r22(i)*yy3+r32(i)*zz3
164
165 xx4=x(1,m)-x(1,j)
166 yy4=x(2,m)-x(2,j)
167 zz4=x(3,m)-x(3,j)
168 xl4(i)=r11(i)*xx4+r21(i)*yy4+r31(i)*zz4
169 yl4(i)=r12(i)*xx4+r22(i)*yy4+r32(i)*zz4
170 ENDDO
171
172
173
174 IF(ismstr==1.OR.ismstr==2)THEN
175 DO i=jft,jlt
176 IF(abs(offg(i))==two)THEN
177 xl2(i)=smstr(ii(1)+i)
178 yl2(i)=smstr(ii(2)+i)
179 xl3(i)=smstr(ii(3)+i)
180 yl3(i)=smstr(ii(4)+i)
181 xl4(i)=smstr(ii(5)+i)
182 yl4(i)=smstr(ii(6)+i)
183 z1(i)=zero
185 . ((xl2(i)-xl4(i))*yl3(i)-xl3(i)*(yl2(i)-yl4(i)))
186 area_i(i)=one/
max(em20,
area(i))
187 ELSE
188 smstr(ii(1)+i)=xl2(i)
189 smstr(ii(2)+i)=yl2(i)
190 smstr(ii(3)+i)=xl3(i)
191 smstr(ii(4)+i)=yl3(i)
192 smstr(ii(5)+i)=xl4(i)
193 smstr(ii(6)+i)=yl4(i)
194 ENDIF
195 ENDDO
196 ENDIF
197 IF(ismstr==1)THEN
198 DO i=jft,jlt
199 IF(offg(i)==one)offg(i)=two
200 ENDDO
201 ENDIF
202
203
204
205 IF (irep > 0) THEN
206 CALL cortdir3(elbuf_str,dir_a,dir_b ,jft ,jlt ,
207 . nlay ,irep ,rx ,ry ,rz ,
208 . sx ,sy ,ssz ,r11 ,r21 ,
209 . r31 ,r12 ,r22 ,r32 ,nel )
210 ENDIF
211
212 DO i=jft,jlt
213 lxyz0(1)=fourth*(xl2(i)+xl3(i)+xl4(i))
214 lxyz0(2)=fourth*(yl2(i)+yl3(i)+yl4(i))
215 corelv(i,1,1)=-lxyz0(1)
216 corelv(i,1,2)=xl2(i)-lxyz0(1)
217 corelv(i,1,3)=xl3(i)-lxyz0(1)
218 corelv(i,1,4)=xl4(i)-lxyz0(1)
219 corelv(i,2,1)=-lxyz0(2)
220 corelv(i,2,2)=yl2(i)-lxyz0(2)
221 corelv(i,2,3)=yl3(i)-lxyz0(2)
222 corelv(i,2,4)=yl4(i)-lxyz0(2)
223 x13(i)=(corelv(i,1,1)-corelv(i,1,3))*half
224 x24(i)=(corelv(i,1,2)-corelv(i,1,4))*half
225 y13(i)=(corelv(i,2,1)-corelv(i,2,3))*half
226 y24(i)=(corelv(i,2,2)-corelv(i,2,4))*half
227
228 mx13(i)=(corelv(i,1,1)+corelv(i,1,3))*half
229 mx23(i)=(corelv(i,1,2)+corelv(i,1,3))*half
230 mx34(i)=(corelv(i,1,3)+corelv(i,1,4))*half
231 my13(i)=(corelv(i,2,1)+corelv(i,2,3))*half
232 my23(i)=(corelv(i,2,2)+corelv(i,2,3))*half
233 my34(i)=(corelv(i,2,3)+corelv(i,2,4))*half
234
235 x13_2(i) =x13(i)*x13(i)
236 y13_2(i) =y13(i)*y13(i)
237 x24_2(i) =x24(i)*x24(i)
238 y24_2(i) =y24(i)*y24(i)
239 l13(i)=x13_2(i)+y13_2(i)
240 l24(i)=x24_2(i)+y24_2(i)
241 ll(i)=half*(l13(i)+l24(i))
242 s1=em01*thk(i)*thk(i)
244 ENDDO
245 IF (imp_chk > 0) THEN
246 s1 =.577350269189626
247 DO i=jft,jlt
248 j1=(mx23(i)*my13(i)-mx13(i)*my23(i))*s1
249 j2=-(mx13(i)*my34(i)-mx34(i)*my13(i))*s1
251 xx1=j0+j2-j1
252 xx2=j0+j2+j1
253 xx3=j0-j2+j1
254 xx4=j0-j2-j1
255 IF(offg(i)/=zero)THEN
256 IF(xx1<=zero.OR.xx2<=zero.OR.
257 . xx3<=zero.OR.xx4<=zero)THEN
258#include "lockon.inc"
259 WRITE(iout ,2001) ngl(i)
260#include "lockoff.inc"
261 idel7nok = 1
262 imp_ir = imp_ir + 1
263 ENDIF
264 ENDIF
265 ENDDO
266 ENDIF
267
268 DO i=jft,jlt
269 z2(i)=z1(i)*z1(i)
270 IF (z2(i)<ll(i)*tol.OR.npt==1) THEN
271 z1(i)=zero
272 ELSE
273
274
275
276
278
279 sz1=mx13(i)*y24(i)-my13(i)*x24(i)
280 sz2=a_4+sz1
281 sz=z2(i)*l24(i)
282 sl=one/sqrt(sz+sz2*sz2)
283 vqn(i,1,1)=-z1(i)*y24(i)
284 vqn(i,2,1)= z1(i)*x24(i)
285 vqn(i,3,1)= sz2*sl
286 vqn(i,1,3)=-vqn(i,1,1)
287 vqn(i,2,3)=-vqn(i,2,1)
288 vqn(i,1,1)= vqn(i,1,1)*sl
289 vqn(i,2,1)= vqn(i,2,1)*sl
290
291 sz2=a_4-sz1
292 sl=one/sqrt(sz+sz2*sz2)
293 vqn(i,1,3)= vqn(i,1,3)*sl
294 vqn(i,2,3)= vqn(i,2,3)*sl
295 vqn(i,3,3)= sz2*sl
296
297 sz1=mx13(i)*y13(i)-my13(i)*x13(i)
298 sz2=a_4+sz1
299 sz=z2(i)*l13(i)
300 sl=one/sqrt(sz+sz2*sz2)
301 vqn(i,1,2)=-z1(i)*y13(i)
302 vqn(i,2,2)= z1(i)*x13(i)
303 vqn(i,3,2)= sz2*sl
304 vqn(i,1,4)=-vqn(i,1,2)
305 vqn(i,2,4)=-vqn(i,2,2)
306 vqn(i,1,2)= vqn(i,1,2)*sl
307 vqn(i,2,2)= vqn(i,2,2)*sl
308
309 sz2=a_4-sz1
310 sl=one/sqrt(sz+sz2*sz2)
311 vqn(i,1,4)= vqn(i,1,4)*sl
312 vqn(i,2,4)= vqn(i,2,4)*sl
313 vqn(i,3,4)= sz2*sl
314 ENDIF
315 ENDDO
316
317 nplat=jft-1
318 splat= 0
319 DO i=jft,jlt
320 IF (z1(i)==zero) THEN
321 nplat=nplat+1
322 iplat(nplat)=i
323 ELSE
324 splat=splat+1
325 jplat(splat)=i
326 ENDIF
327 ENDDO
328 DO i=nplat+1,jlt
329 iplat(i)=jplat(i-nplat)
330 ENDDO
331
332 DO i=jft,jlt
333 off(i)=offg(i)
334 ENDDO
335
336 RETURN
337 2001 FORMAT(/' ZERO OR NEGATIVE SHELL SUB-AREA : ELEMENT NB:',i8/)
subroutine clskew3(jft, jlt, irep, rx, ry, rz, sx, sy, sz, e1x, e2x, e3x, e1y, e2y, e3y, e1z, e2z, e3z, det)
subroutine cortdir3(elbuf_str, dir_a, dir_b, jft, jlt, nlay, irep, rx, ry, rz, sx, sy, sz, e1x, e1y, e1z, e2x, e2y, e2z, nel)
subroutine area(d1, x, x2, y, y2, eint, stif0)