37
38
39
40 USE elbufdef_mod
43 use element_mod , only : nixc
44
45
46
47#include "implicit_f.inc"
48
49
50
51#include "param_c.inc"
52#include "com01_c.inc"
53#include "com04_c.inc"
54#include "scr03_c.inc"
55
56
57
58 TYPE(ELBUF_STRUCT_), TARGET, DIMENSION(NGROUP) :: ELBUF_TAB
59 TYPE (NLOCAL_STR_) , TARGET :: NLOC_DMG
60 INTEGER IXC(NIXC,*),NFT,NEL,NG,IPM(NPROPMI,*)
61 my_real ,
DIMENSION(NUMELC+NUMELTG),
INTENT(IN) ::
64 . x(3,*),xrefc(4,3,*),dt_nl,bufmat(*),time
65 LOGICAL :: FAILURE
66
67
68
69 INTEGER :: IMAT,NDOF,L_NLOC,N1,N2,N3,N4,K,I,NDNOD
70 INTEGER, DIMENSION(NEL) :: POS1,POS2,POS3,POS4
72 . le_min,len,damp,dens,ntn_unl,ntn_vnl,
73 . ntvar,z01(11,11),wf1(11,11),zn1(12,11),b1,b2,
74 . b3,b4,nth1,nth2,bth1,bth2,k1,k12,k2,sspnl,le_max
75 my_real,
DIMENSION(:,:),
ALLOCATABLE :: var_reg,vpred
77 . DIMENSION(:), POINTER :: fnl,unl,vnl,dnl,mnl,thck
78 my_real,
DIMENSION(NEL) :: x1,x2,x3,x4,
79 . y1,y2,y3,y4,px1,px2,py1,py2,e1x,e2x,e3x,
80 . e1y,e2y,e3y,e1z,e2z,e3z,x2l,y2l,x3l,y3l,
81 . x4l,y4l,z1,z2,z3,z4,surf,offg,vols,btb11,
82 . btb12,btb22
83 TYPE(BUF_NLOC_),POINTER :: BUFNL
84 my_real,
DIMENSION(:,:),
POINTER ::
85 . massth,fnlth,vnlth,unlth
86
87 DATA z01/
88 1 0. ,0. ,0. ,0. ,0. ,
89 1 0. ,0. ,0. ,0. ,0. ,0. ,
90 2 -.5 ,0.5 ,0. ,0. ,0. ,
91 2 0. ,0. ,0. ,0. ,0. ,0. ,
92 3 -.5 ,0. ,0.5 ,0. ,0. ,
93 3 0. ,0. ,0. ,0. ,0. ,0. ,
94 4 -.5 ,-.1666667,0.1666667,0.5 ,0. ,
95 4 0. ,0. ,0. ,0. ,0. ,0. ,
96 5 -.5 ,-.25 ,0. ,0.25 ,0.5 ,
97 5 0. ,0. ,0. ,0. ,0. ,0. ,
98 6 -.5 ,-.3 ,-.1 ,0.1 ,0.3 ,
99 6 0.5 ,0. ,0. ,0. ,0. ,0. ,
100 7 -.5 ,-.3333333,-.1666667,0.0 ,0.1666667,
101 7 0.3333333,0.5 ,0. ,0. ,0. ,0. ,
102 8 -.5 ,-.3571429,-.2142857,-.0714286,0.0714286,
103 8 0.2142857,0.3571429,0.5 ,0. ,0. ,0. ,
104 9 -.5 ,-.375 ,-.25 ,-.125 ,0.0 ,
105 9 0.125 ,0.25 ,0.375 ,0.5 ,0. ,0. ,
106 a -.5 ,-.3888889,-.2777778,-.1666667,-.0555555,
107 a 0.0555555,0.1666667,0.2777778,0.3888889,0.5 ,0. ,
108 b -.5 ,-.4 ,-.3 ,-.2 ,-.1 ,
109 b 0. ,0.1 ,0.2 ,0.3 ,0.4 ,0.5 /
110
111 DATA wf1/
112 1 1. ,0. ,0. ,0. ,0. ,
113 1 0. ,0. ,0. ,0. ,0. ,0. ,
114 2 0.5 ,0.5 ,0. ,0. ,0. ,
115 2 0. ,0. ,0. ,0. ,0. ,0. ,
116 3 0.25 ,0.5 ,0.25 ,0. ,0. ,
117 3 0. ,0. ,0. ,0. ,0. ,0. ,
118 4 0.1666667,0.3333333,0.3333333,0.1666667,0. ,
119 4 0. ,0. ,0. ,0. ,0. ,0. ,
120 5 0.125 ,0.25 ,0.25 ,0.25 ,0.125 ,
121 5 0. ,0. ,0. ,0. ,0. ,0. ,
122 6 0.1 ,0.2 ,0.2 ,0.2 ,0.2 ,
123 6 0.1 ,0. ,0. ,0. ,0. ,0. ,
124 7 0.0833333,0.1666667,0.1666667,0.1666667,0.1666667,
125 7 0.1666667,0.0833333,0. ,0. ,0. ,0. ,
126 8 0.0714286,0.1428571,0.1428571,0.1428571,0.1428571,
127 8 0.1428571,0.1428571,0.0714286,0. ,0. ,0. ,
128 9 0.0625 ,0.125 ,0.125 ,0.125 ,0.125 ,
129 9 0.125 ,0.125 ,0.125 ,0.0625 ,0. ,0. ,
130 a 0.0555556,0.1111111,0.1111111,0.1111111,0.1111111,
131 a 0.1111111,0.1111111,0.1111111,0.1111111,0.0555556,0. ,
132 b 0.05 ,0.1 ,0.1 ,0.1 ,0.1 ,
133 b 0.1 ,0.1 ,0.1 ,0.1 ,0.1 ,0.05 /
134
135 DATA zn1/
136 1 0. ,0. ,0. ,0. ,0. ,0. ,
137 1 0. ,0. ,0. ,0. ,0. ,0. ,
138 2 -.5 ,0.5 ,0. ,0. ,0. ,0. ,
139 2 0. ,0. ,0. ,0. ,0. ,0. ,
140 3 -.5 ,-.25 ,0.25 ,0.5 ,0. ,0. ,
141 3 0. ,0. ,0. ,0. ,0. ,0. ,
142 4 -.5 ,-.3333333,0. ,0.3333333,0.5 ,0. ,
143 4 0. ,0. ,0. ,0. ,0. ,0. ,
144 5 -.5 ,-.375 ,-0.125 ,0.125 ,0.375 ,0.5 ,
145 5 0. ,0. ,0. ,0. ,0. ,0. ,
146 6 -.5 ,-.4 ,-.2 ,0.0 ,0.2 ,0.4 ,
147 6 0.5 ,0. ,0. ,0. ,0. ,0. ,
148 7 -.5 ,-.4166667,-.25 ,-.0833333,0.0833333,0.25 ,
149 7 0.4166667,0.5 ,0. ,0. ,0. ,0. ,
150 8 -.5 ,-.4285715,-.2857143,-.1428572,0.0 ,0.1428572,
151 8 0.2857143,0.4285715,0.5 ,0. ,0. ,0. ,
152 9 -.5 ,-.4375 ,-.3125 ,-.1875 ,-.0625 ,0.0625 ,
153 9 0.1875 ,0.3125 ,0.4375 ,0.5 ,0. ,0. ,
154 a -.5 ,-.4444444,-.3333333,-.2222222,-.1111111,0. ,
155 a 0.1111111,0.2222222,0.3333333,0.4444444,0.5 ,0. ,
156 b -.5 ,-.45 ,-.35 ,-.25 ,-.15 ,-.05 ,
157 b 0.05 ,0.15 ,0.25 ,0.35 ,0.45 ,0.5 /
158
159
160 l_nloc = nloc_dmg%L_NLOC
161
162 fnl => nloc_dmg%FNL(1:l_nloc,1)
163 vnl => nloc_dmg%VNL(1:l_nloc)
164 dnl => nloc_dmg%DNL(1:l_nloc)
165 unl => nloc_dmg%UNL(1:l_nloc)
166 mnl => nloc_dmg%MASS(1:l_nloc)
167
168 imat = ixc(1,1+nft)
169
170 le_min = sqrt(minval(
area(nft+1:nft+nel)))
171
172 ndof = elbuf_tab(ng)%BUFLY(1)%NPTT
173
174 thck => elbuf_tab(ng)%GBUF%THK(1:nel)
175
176 IF (ndof>1) THEN
177 le_min =
min(le_min,minval(thck(1:nel))/ndof)
178 ENDIF
179
180 len = nloc_dmg%LEN(imat)
181
182 le_max = nloc_dmg%LE_MAX(imat)
183
184 damp = nloc_dmg%DAMP(imat)
185
186 dens = nloc_dmg%DENS(imat)
187
188 sspnl = nloc_dmg%SSPNL(imat)
189
190 dt_nl =
min(dt_nl,0.5d0*((two*
min(le_min,le_max)*sqrt(three*dens))/
191 . (sqrt(twelve*(len**2)+(
min(le_min,le_max)**2)))))
192
193 IF (ndof>1) THEN
194 IF (ndof > 2) THEN
195 ALLOCATE(vpred(nel,ndof+1))
196 ndnod = ndof + 1
197 ELSE
198 ALLOCATE(vpred(nel,ndof))
199 ndnod = ndof
200 ENDIF
201 ENDIF
202
203 IF (.NOT.ALLOCATED(var_reg)) ALLOCATE(var_reg(nel,ndof))
204
205
206# include "vectorize.inc"
207 DO i = 1,nel
208
209 IF (nxref == 0) THEN
210 x1(i)=x(1,ixc(2,nft+i))
211 y1(i)=x(2,ixc(2,nft+i))
212 z1(i)=x(3,ixc(2,nft+i))
213 x2(i)=x(1,ixc(3,nft+i))
214 y2(i)=x(2,ixc(3,nft+i))
215 z2(i)=x(3,ixc(3,nft+i))
216 x3(i)=x(1,ixc(4,nft+i))
217 y3(i)=x(2,ixc(4,nft+i))
218 z3(i)=x(3,ixc(4,nft+i))
219 x4(i)=x(1,ixc(5,nft+i))
220 y4(i)=x(2,ixc(5,nft+i))
221 z4(i)=x(3,ixc(5,nft+i))
222 ELSE
223 x1(i)=xrefc(1,1,nft+i)
224 y1(i)=xrefc(1,2,nft+i)
225 z1(i)=xrefc(1,3,nft+i)
226 x2(i)=xrefc(2,1,nft+i)
227 y2(i)=xrefc(2,2,nft+i)
228 z2(i)=xrefc(2,3,nft+i)
229 x3(i)=xrefc(3,1,nft+i)
230 y3(i)=xrefc(3,2,nft+i)
231 z3(i)=xrefc(3,3,nft+i)
232 x4(i)=xrefc(4,1,nft+i)
233 y4(i)=xrefc(4,2,nft+i)
234 z4(i)=xrefc(4,3,nft+i)
235 ENDIF
236
237
238 n1 = nloc_dmg%IDXI(ixc(2,nft+i))
239 n2 = nloc_dmg%IDXI(ixc(3,nft+i))
240 n3 = nloc_dmg%IDXI(ixc(4,nft+i))
241 n4 = nloc_dmg%IDXI(ixc(5,nft+i))
242
243 pos1(i) = nloc_dmg%POSI(n1)
244 pos2(i) = nloc_dmg%POSI(n2)
245 pos3(i) = nloc_dmg%POSI(n3)
246 pos4(i) = nloc_dmg%POSI(n4)
247 ENDDO
248
249
250
251 DO k = 1,ndof
252
253 DO i = 1,nel
254 var_reg(i,k) = fourth*(dnl(pos1(i)+k-1) + dnl(pos2(i)+k-1)
255 . + dnl(pos3(i)+k-1) + dnl(pos4(i)+k-1))
256 ENDDO
257 ENDDO
258
260 . x1 ,x2 ,x3 ,x4 ,y1 ,y2 ,
261 . y3 ,y4 ,z1 ,z2 ,z3 ,z4 ,
262 . e1x ,e2x ,e3x ,e1y ,e2y ,e3y ,
263 . e1z ,e2z ,e3z )
264
265
267 . bufmat ,time ,var_reg ,
268 . failure )
269
270
271
272
273
274# include "vectorize.inc"
275 DO i=1,nel
276
277
278 x2l(i) = e1x(i)*(x2(i)-x1(i))+e1y(i)*(y2(i)-y1(i))+e1z(i)*(z2(i)-z1(i))
279 y2l(i) = e2x(i)*(x2(i)-x1(i))+e2y(i)*(y2(i)-y1(i))+e2z(i)*(z2(i)-z1(i))
280 x3l(i) = e1x(i)*(x3(i)-x1(i))+e1y(i)*(y3(i)-y1(i))+e1z(i)*(z3(i)-z1(i))
281 y3l(i) = e2x(i)*(x3(i)-x1(i))+e2y(i)*(y3(i)-y1(i))+e2z(i)*(z3(i)-z1(i))
282 x4l(i) = e1x(i)*(x4(i)-x1(i))+e1y(i)*(y4(i)-y1(i))+e1z(i)*(z4(i)-z1(i))
283 y4l(i) = e2x(i)*(x4(i)-x1(i))+e2y(i)*(y4(i)-y1(i))+e2z(i)*(z4(i)-z1(i))
284
285 px1(i) = half *(y2l(i)-y4l(i))
286 px2(i) = half * y3l(i)
287 py1(i) = -half *(x2l(i)-x4l(i))
288 py2(i) = -half * x3l(i)
289
290
291 btb11(i) = px1(i)**2 + py1(i)**2
292 btb12(i) = px1(i)*px2(i) + py1(i)*py2(i)
293 btb22(i) = px2(i)**2 + py2(i)**2
294
295
296 vols(i) =
area(nft+i)*thck(i)
297
298
299 offg(i) = elbuf_tab(ng)%GBUF%OFF(i)
300
301 ENDDO
302
303
304
305
306
307 IF ((ndof > 1).AND.(len>zero)) THEN
308
309
310 bufnl => elbuf_tab(ng)%NLOC(1,1)
311
312
313 massth => bufnl%MASSTH(1:nel,1:ndnod)
314 unlth => bufnl%UNLTH(1:nel ,1:ndnod)
315 vnlth => bufnl%VNLTH(1:nel ,1:ndnod)
316 fnlth => bufnl%FNLTH(1:nel ,1:ndnod)
317
318 DO k = 1,ndnod
319 DO i = 1,nel
320
321 vpred(i,k) = vnlth(i,k) - (fnlth(i,k)/massth(i,k))*(dt_nl/two)
322 ENDDO
323 ENDDO
324 DO k = 1,ndnod
325 DO i = 1,nel
326
327 fnlth(i,k) = zero
328 ENDDO
329 ENDDO
330
331
332 DO k = 1, ndof
333
334
335 IF ((ndof==2).AND.(k==2)) THEN
336 nth1 = (z01(k,ndof) - zn1(k,ndof))/(zn1(k-1,ndof) - zn1(k,ndof))
337 nth2 = (z01(k,ndof) - zn1(k-1,ndof)) /(zn1(k,ndof) - zn1(k-1,ndof))
338 ELSE
339 nth1 = (z01(k,ndof) - zn1(k+1,ndof))/(zn1(k,ndof) - zn1(k+1,ndof))
340 nth2 = (z01(k,ndof) - zn1(k,ndof)) /(zn1(k+1,ndof) - zn1(k,ndof))
341 ENDIF
342
343
344 DO i = 1,nel
345
346 IF ((ndof==2).AND.(k==2)) THEN
347 bth1 = (one/(zn1(k-1,ndof) - zn1(k,ndof)))*(one/thck(i))
348 bth2 = (one/(zn1(k,ndof) - zn1(k-1,ndof)))*(one/thck(i))
349 ELSE
350 bth1 = (one/(zn1(k,ndof) - zn1(k+1,ndof)))*(one/thck(i))
351 bth2 = (one/(zn1(k+1,ndof) - zn1(k,ndof)))*(one/thck(i))
352 ENDIF
353
354
355 k1 = (len**2)*(bth1**2) + nth1**2
356 k12 = (len**2)*(bth1*bth2)+ (nth1*nth2)
357 k2 = (len**2)*(bth2**2) + nth2**2
358
359
360 IF ((ndof==2).AND.(k==2)) THEN
361 fnlth(i,k-1) = fnlth(i,k-1) + (k1*unlth(i,k-1) + k12*unlth(i,k)
362 . + damp*((nth1**2)*vpred
363 . + (nth1*nth2)*vpred(i,k))
364 . - (nth1*var_reg(i,k)))*vols(i)*wf1(k,ndof)
365 fnlth(i,k) = fnlth(i,k) + (k12*unlth(i,k-1) + k2*unlth(i,k)
366 . + damp*(nth1*nth2*vpred(i,k-1)
367 . + (nth2**2)*vpred(i,k))
368 . - nth2*var_reg(i,k))*vols(i)*wf1(k,ndof)
369 ELSE
370 fnlth(i,k) = fnlth(i,k) + (k1*unlth(i,k) + k12*unlth(i,k+1)
371 . + damp*((nth1**2)*vpred(i,k)
372 . + (nth1*nth2)*vpred(i,k+1))
373 . - (nth1*var_reg(i,k)))*vols(i)*wf1(k,ndof)
374 fnlth(i,k+1) = fnlth(i,k+1) + (k12*unlth(i,k) + k2*unlth(i,k+1)
375 . + damp*(nth1*nth2*vpred(i,k)
376 . + (nth2**2)*vpred(i,k+1))
377 . - nth2*var_reg(i,k))*vols(i)*wf1(k,ndof)
378 ENDIF
379 ENDDO
380 ENDDO
381
382 DO k = 1,ndnod
383 DO i = 1,nel
384
385 vnlth(i,k) = vnlth(i,k) - (fnlth(i,k)/massth(i,k))*dt_nl
386 ENDDO
387 ENDDO
388
389 DO k = 1,ndnod
390 DO i = 1,nel
391
392 unlth(i,k) = unlth(i,k) + vnlth(i,k)*dt_nl
393 ENDDO
394 ENDDO
395
396
397 DO k = 1, ndof
398
399 IF ((ndof==2).AND.(k==2)) THEN
400 nth1 = (z01(k,ndof) - zn1(k,ndof))/(zn1(k-1,ndof) - zn1(k,ndof))
401 nth2 = (z01(k,ndof) - zn1(k-1,ndof)) /(zn1(k,ndof) - zn1(k-1,ndof))
402 ELSE
403 nth1 = (z01(k,ndof) - zn1(k+1,ndof))/(zn1(k,ndof) - zn1(k+1,ndof))
404 nth2 = (z01(k,ndof) - zn1(k,ndof)) /(zn1(k+1,ndof) - zn1(k,ndof))
405 ENDIF
406
407 DO i = 1,nel
408
409 IF ((ndof==2).AND.(k==2)) THEN
410 var_reg(i,k) = nth1*unlth(i,k-1) + nth2*unlth(i,k)
411 ELSE
412 var_reg(i,k) = nth1*unlth(i,k) + nth2*unlth(i,k+1)
413 ENDIF
414 ENDDO
415 ENDDO
416 ENDIF
417
418
419 ! computation of
the elementary non-local forces
420
421
422 DO k = 1,ndof
423
424
425# include "vectorize.inc"
426 DO i = 1,nel
427
428
429 IF (offg(i) > zero) THEN
430
431 b1 = ((len**2)/vols(i))*wf1(k,ndof)*(btb11(i)*unl(pos1(i)+k-1) + btb12(i)*unl(pos2(i)+k-1)
432 . - btb11(i)*unl
433 b2 = ((len**2)/vols(i))*wf1(k,ndof)*(btb12(i)*unl
434 . - btb12(i)*unl(pos3(i)+k-1) - btb22(i)*unl(pos4(i)+k-1))
435 b3 = ((len**2)/vols(i))*wf1(k,ndof)*(-btb11(i)*unl(pos1(i)+k-1) - btb12(i)*unl(pos2(i)+k-1)
436 . + btb11(i)*unl(pos3(i)+k-1) + btb12(i)*unl(pos4(i)+k-1))
437 b4 = ((len**2)/vols(i))*wf1(k,ndof)*(-btb12(i)*unl(pos1(i)+k-1) - btb22(i)*unl(pos2(i)+k-1)
438 . + btb12(i)*unl(pos3(i)+k-1) + btb22(i)*unl(pos4(i)+k-1))
439
440 ntn_unl = ((unl(pos1(i)+k-1) + unl(pos2(i)+k-1) + unl(pos3(i)+k-1) +
441 . unl(pos4(i)+k-1))*fourth*fourth)*vols(i)*wf1(k,ndof)
442
443 ntn_vnl = ((vnl(pos1(i)+k-1) + vnl(pos2(i)+k-1) + vnl(pos3(i)+k-1) +
444 . vnl(pos4(i)+k-1))*fourth*fourth)*damp*vols(i)*wf1(k,ndof)
445
446 ntvar = var_reg(i,k)*fourth*vols(i)*wf1(k,ndof)
447
448 fnl(pos1(i)+k-1) = fnl(pos1(i)+k-1) - (ntn_unl + ntn_vnl - ntvar + b1)
449 fnl(pos2(i)+k-1) = fnl(pos2(i)+k-1) - (ntn_unl + ntn_vnl - ntvar + b2)
450 fnl(pos3(i)+k-1) = fnl(pos3(i)+k-1) - (ntn_unl + ntn_vnl - ntvar + b3)
451 fnl(pos4(i)+k-1) = fnl(pos4(i)+k-1) - (ntn_unl + ntn_vnl - ntvar + b4)
452
453 ELSE
454
455 fnl(pos1(i)+k-1) = fnl(pos1(i)+k-1) - wf1(k,ndof)*dens*sspnl*vnl(pos1(i)+k-1)*le_max*thck(i)
456 fnl(pos2(i)+k-1) = fnl(pos2(i)+k-1) - wf1(k,ndof)*dens*sspnl*vnl(pos2(i)+k-1)*le_max*thck(i)
457 fnl(pos3(i)+k-1) = fnl(pos3(i)+k-1) - wf1(k,ndof)*dens*sspnl*vnl(pos3(i)+k-1)*le_max*thck(i)
458 fnl(pos4(i)+k-1) = fnl(pos4(i)+k-1) - wf1(k,ndof)*dens*sspnl*vnl(pos4(i)+k-1)*le_max*thck(i)
459 ENDIF
460 ENDDO
461 ENDDO
462
463 IF (ALLOCATED(var_reg)) DEALLOCATE(var_reg)
464 IF (ALLOCATED(vpred)) DEALLOCATE(vpred)
subroutine ceveci(jft, jlt, area, x1, x2, x3, x4, y1, y2, y3, y4, z1, z2, z3, z4, e1x, e2x, e3x, e1y, e2y, e3y, e1z, e2z, e3z)
subroutine cnloc_matini(elbuf_str, nel, ipm, bufmat, time, varnl, failure)
end diagonal values have been computed in the(sparse) matrix id.SOL
subroutine area(d1, x, x2, y, y2, eint, stif0)