38
39
40
42 USE elbufdef_mod
43
44
45
46#include "implicit_f.inc"
47
48
49
50#include "parit_c.inc"
51#include "com20_c.inc"
52#include "com08_c.inc"
53#include "scr02_c.inc"
54#include "scr18_c.inc"
55
56
57
58 INTEGER, INTENT(IN) :: NFT
59 INTEGER :: NEL,IMAT,NDDL,ITASK
60 INTEGER, DIMENSION(NEL) :: NC1,NC2,NC3
61 my_real,
DIMENSION(NEL,NDDL),
INTENT(INOUT)::
62 . var_reg
63 my_real,
DIMENSION(NEL),
INTENT(IN) ::
64 .
area,off,px2,py2,px3,py3,thk,le,thk0,area0
66 . dt2t
67 TYPE(NLOCAL_STR_), TARGET :: NLOC_DMG
68 TYPE(BUF_NLOC_) , TARGET :: BUFNL
69
70
71
72 INTEGER I,K,N1,N2,N3,L_NLOC,II,J,NDNOD
74 . l2,ntn,ntn_unl,ntn_vnl,xi,ntvar,a,
75 . b1,b2,b3,zeta,sspnl,
76 . nth1,nth2,bth1,bth2,k1,k2,k12,
77 . dtnl_th,dtnl,le_max,maxstif,
78 . dtnod,dt2p
79 my_real,
DIMENSION(:,:),
ALLOCATABLE ::
80 . f1,f2,f3,sti1,sti2,sti3
81 my_real,
DIMENSION(:) ,
ALLOCATABLE ::
82 . btb11,btb12,btb13,btb22,btb23,btb33,vol
83 INTEGER, DIMENSION(:), ALLOCATABLE ::
84 . POS1,,POS3
85 my_real,
POINTER,
DIMENSION(:) ::
86 . vnl,fnl,unl,stifnl,mass,mass0,vnl0
87 my_real,
POINTER,
DIMENSION(:,:) ::
88 . massth,unlth,vnlth,fnlth
89 my_real,
DIMENSION(:,:),
ALLOCATABLE ::
90 . stifnlth,dtn
91
92
93
94
95
96 l2 = nloc_dmg%LEN(imat)**2
97 xi = nloc_dmg%DAMP(imat)
98 zeta = nloc_dmg%DENS(imat)
99 sspnl = nloc_dmg%SSPNL(imat)
100 l_nloc = nloc_dmg%L_NLOC
101 le_max = nloc_dmg%LE_MAX(imat)
102
103 ALLOCATE(f1(nel,nddl),f2(nel,nddl),f3(nel,nddl))
104 ! only
for nodal timestep
105 IF (nodadt > 0) THEN
106
107 ALLOCATE(sti1(nel,nddl),sti2(nel,nddl),sti3(nel,nddl))
108
109 mass => nloc_dmg%MASS(1:l_nloc)
110
111 mass0 => nloc_dmg%MASS0(1:l_nloc)
112 ENDIF
113 ALLOCATE(btb11(nel),btb12(nel),btb13(nel),btb22(nel),
114 . btb23(nel),btb33(nel),vol(nel),pos1(nel),
115 . pos2(nel),pos3(nel))
116
117 vnl => nloc_dmg%VNL(1:l_nloc)
118 vnl0 => nloc_dmg%VNL_OLD(1:l_nloc)
119 unl => nloc_dmg%UNL(1:l_nloc)
120
121 ntn = three*three
122
123
124
125
126
127# include "vectorize.inc"
128 DO i=1,nel
129
130
131 n1 = nloc_dmg%IDXI(nc1(i))
132 n2 = nloc_dmg%IDXI(nc2(i))
133 n3 = nloc_dmg%IDXI(nc3(i))
134
135
136 pos1(i) = nloc_dmg%POSI(n1)
137 pos2(i) = nloc_dmg%POSI(n2)
138 pos3(i) = nloc_dmg%POSI(n3)
139
140
141 vol(i) =
area(i)*thk(i)
142
143
144 btb11(i) = (-px3(i)-px2(i))**2 + (-py2(i)-py3(i))**2
145 btb12(i) = (-px3(i)-px2(i))*px2(i) + (-py2(i)-py3(i))*py2(i)
146 btb13(i) = (-px3(i)-px2(i))*px3(i) + (-py2(i)-py3(i))*py3(i)
147 btb22(i) = px2(i)**2 + py2(i)**2
148 btb23(i) = px2(i)*px3(i) + py2(i)*py3(i)
149 btb33(i) = px3(i)**2 + py3(i)**2
150
151 ENDDO
152
153
154
155
156
157 IF ((nddl > 1).AND.(l2>zero)) THEN
158
159
160 IF (nddl > 2) THEN
161 IF (nodadt > 0) THEN
162 ALLOCATE(stifnlth(nel,nddl+1))
163 ALLOCATE(dtn(nel,nddl+1))
164 ENDIF
165 ndnod = nddl+1
166 ELSE
167 IF (nodadt > 0) THEN
168 ALLOCATE(stifnlth(nel,nddl))
169 ALLOCATE(dtn(nel,nddl))
170 ENDIF
171 ndnod = nddl
172 ENDIF
173
174
175 massth => bufnl%MASSTH(1:nel,1:ndnod)
176 unlth => bufnl%UNLTH(1:nel,1:ndnod)
177 vnlth => bufnl%VNLTH(1:nel,1:ndnod)
178 fnlth => bufnl%FNLTH(1:nel,1:ndnod)
179
180 DO k = 1,ndnod
181 DO i = 1,nel
182
183 fnlth(i,k) = zero
184
185 IF (nodadt > 0) THEN
186 stifnlth(i,k) = em20
187 ENDIF
188 ENDDO
189 ENDDO
190
191
192 DO k = 1, nddl
193
194
195 IF ((nddl==2).AND.(k==2)) THEN
196 nth1 = (z0(k,nddl) - zth(k,nddl)) / (zth(k-1,nddl) - zth(k,nddl))
197 nth2 = (z0(k,nddl) - zth(k-1,nddl)) / (zth(k,nddl) - zth(k-1,nddl))
198 ELSE
199 nth1 = (z0(k,nddl) - zth(k+1,nddl)) / (zth(k,nddl) - zth(k+1,nddl))
200 nth2 = (z0(k,nddl) - zth(k,nddl)) / (zth(k+1,nddl) - zth(k,nddl))
201 ENDIF
202
203
204 DO i = 1,nel
205
206 IF ((nddl==2).AND.(k==2)) THEN
207 bth1 = (one/(zth(k-1,nddl) - zth(k,nddl)))*(one/thk(i))
208 bth2 = (one/(zth(k,nddl) - zth(k-1,nddl)))*(one/thk(i))
209 ELSE
210 bth1 = (one/(zth(k,nddl) - zth(k+1,nddl)))*(one/thk(i))
211 bth2 = (one/(zth(k+1,nddl) - zth(k,nddl)))*(one/thk(i))
212 ENDIF
213
214
215 k1 = l2*(bth1**2) + nth1**2
216 k12 = l2*(bth1*bth2)+ (nth1*nth2)
217 k2 = l2*(bth2**2) + nth2**2
218
219
220 IF ((nddl==2).AND.(k==2)) THEN
221 fnlth(i,k-1) = fnlth(i,k-1) + (k1*unlth(i,k-1) + k12*unlth(i,k)
222 . + xi*((nth1**2)*vnlth(i,k-1)
223 . + (nth1*nth2)*vnlth(i,k))
224 . - (nth1*var_reg(i,k)))*vol(i)*wf(k,nddl)
225 fnlth(i,k) = fnlth(i,k) + (k12*unlth(i,k-1) + k2*unlth(i,k)
226 . + xi*(nth1*nth2*vnlth(i,k-1)
227 . + (nth2**2)*vnlth(i,k))
228 . - nth2*var_reg(i,k))*vol(i)*wf(k,nddl)
229 ELSE
230 fnlth(i,k) = fnlth(i,k) + (k1*unlth(i,k) + k12*unlth(i,k+1)
231 . + xi*((nth1**2)*vnlth(i,k)
232 . + (nth1*nth2)*vnlth(i,k+1))
233 . - (nth1*var_reg(i,k)))*vol(i)*wf(k,nddl)
234 fnlth(i,k+1) = fnlth(i,k+1) + (k12*unlth(i,k) + k2*unlth(i,k+1)
235 . + xi*(nth1*nth2*vnlth(i,k)
236 . + (nth2**2)*vnlth(i,k+1))
237 . - nth2*var_reg(i,k))*vol(i)*wf(k,nddl)
238 ENDIF
239
240
241 IF (nodadt > 0) THEN
242 IF ((nddl==2).AND.(k==2)) THEN
243 stifnlth(i,k-1) = stifnlth(i,k-1) +
max(abs(k1)+abs(k12),abs(k12)+abs(k2))*vol(i)*wf(k,nddl)
244 stifnlth(i,k) = stifnlth(i,k) +
max(abs(k1)+abs(k12),abs(k12)+abs(k2))*vol(i)*wf(k,nddl)
245 ELSE
246 stifnlth(i,k) = stifnlth(i,k) +
max(abs(k1)+abs(k12),abs(k12)+abs(k2))*vol(i)*wf(k,nddl)
247 stifnlth(i,k+1) = stifnlth(i,k+1) +
max(abs(k1)+abs(k12),abs(k12)+abs(k2))*vol(i
248 ENDIF
249 ENDIF
250
251 ENDDO
252 ENDDO
253
254
255 IF (nodadt > 0) THEN
256
257
258 dtnod = ep20
259 DO k = 1,ndnod
260 DO i = 1,nel
261 dtn(i,k) = dtfac1(11)*cdamp*sqrt(two*massth(i,k)/
max(stifnlth(i,k),em20))
262 dtnod =
min(dtn(i,k),dtnod)
263 ENDDO
264 ENDDO
265
266
267 IF ((idtmin(11)==3).OR.(idtmin(11)==4).OR.(idtmin(11)==8)) THEN
268
269 IF (dtnod < dtmin1(11)*(sqrt(csta))) THEN
270 DO k = 1,ndnod
271 DO i = 1,nel
272 IF (dtn(i,k) < dtmin1(11)) THEN
273 dt2p = dtmin1(11)/(dtfac1(11)*cdamp)
274 massth(i,k) =
max(massth(i,k),csta*half*stifnlth(i,k)*dt2p*dt2p*onep00001)
275 ENDIF
276 ENDDO
277 ENDDO
278 dtnod = dtmin1(11)*(sqrt(csta))
279 ENDIF
280 ENDIF
281
282
283 IF (dtnod < dt2t) THEN
284 dt2t =
min(dt2t,dtnod)
285 ENDIF
286 ENDIF
287
288 DO k = 1,ndnod
289 DO i = 1,nel
290
291 vnlth(i,k) = vnlth(i,k) - (fnlth(i,k)/massth(i,k))*dt12
292 ENDDO
293 ENDDO
294
295 DO k = 1,ndnod
296 DO i = 1,nel
297
298 unlth(i,k) = unlth(i,k) + vnlth(i,k)*dt1
299 ENDDO
300 ENDDO
301
302
303 DO k = 1, nddl
304
305 IF ((nddl==2).AND.(k==2)) THEN
306 nth1 = (z0(k,nddl) - zth(k,nddl))/(zth(k-1,nddl) - zth(k,nddl))
307 nth2 = (z0(k,nddl) - zth(k-1,nddl)) /(zth(k,nddl) - zth(k-1,nddl))
308 ELSE
309 nth1 = (z0(k,nddl) - zth(k+1,nddl))/(zth(k,nddl) - zth(k+1,nddl))
310 nth2 = (z0(k,nddl) - zth(k,nddl)) /(zth(k+1,nddl) - zth(k,nddl))
311 ENDIF
312
313 DO i = 1,nel
314
315 IF ((nddl==2).AND.(k==2)) THEN
316 var_reg(i,k) = nth1*unlth(i,k-1) + nth2*unlth(i,k)
317 ELSE
318 var_reg(i,k) = nth1*unlth(i,k) + nth2*unlth(i,k+1)
319 ENDIF
320 ENDDO
321 ENDDO
322 ENDIF
323
324
325
326
327
328 DO k = 1, nddl
329
330
331# include "vectorize.inc"
332 DO i=1,nel
333
334
335 IF (off(i)/=zero) THEN
336
337 b1 = (l2 * vol(i)) * wf(k,nddl)*(btb11(i)*unl(pos1(i)+k-1) + btb12(i)*unl(pos2(i)+k-1)
338 . + btb13(i)*unl(pos3(i)+k-1))
339
340 b2 = (l2 * vol(i)) * wf(k,nddl)*(btb12(i)*unl(pos1(i)+k-1) + btb22(i)*unl(pos2(i)+k-1)
341 . + btb23(i)*unl(pos3(i)+k-1))
342
343 b3 = (l2 * vol(i)) * wf(k,nddl)*(btb13(i)*unl(pos1(i)+k-1) + btb23(i)*unl(pos2(i)+k-1)
344 . + btb33(i)*unl(pos3(i)+k-1))
345
346
347 ntn_unl = (unl(pos1(i)+k-1) + unl(pos2(i)+k-1) + unl(pos3(i)+k-1))/ntn
348
349
350 ntn_vnl = (vnl(pos1(i)+k-1) + vnl(pos2(i)+k-1) + vnl(pos3(i)+k-1))/ntn
351 IF (nodadt > 0) THEN
352 ntn_vnl =
min(sqrt(mass(pos1(i)+k-1)/mass0(pos1(i)+k-1)),
353 . sqrt(mass(pos2(i)+k-1)/mass0(pos2(i)+k-1)),
354 . sqrt(mass(pos3(i)+k-1)/mass0(pos3(i)+k-1)))*ntn_vnl
355 ENDIF
356
357
358 ntn_unl = ntn_unl * vol(i) * wf(k,nddl)
359 ntn_vnl = ntn_vnl * xi * vol(i) * wf(k,nddl)
360
361
362 ntvar = var_reg(i,k)*third*vol(i)*wf(k,nddl)
363
364
365 a = ntn_unl + ntn_vnl - ntvar
366 f1(i,k) = a + b1
367 f2(i,k) = a + b2
368 f3(i,k) = a + b3
369
370
371 IF (nodadt > 0) THEN
372 sti1(i,k) = (abs(l2*btb11(i) + one/ntn) + abs(l2*btb12(i) + one/ntn) +
373 . abs(l2*btb13(i) + one/ntn))* vol(i) * wf(k,nddl)
374 sti2(i,k) = (abs(l2*btb12(i) + one/ntn) + abs(l2*btb22(i) + one/ntn) +
375 . abs(l2*btb23(i) + one/ntn))* vol(i) * wf(k,nddl)
376 sti3(i,k) = (abs(l2*btb13(i) + one/ntn) + abs(l2*btb23(i) + one/ntn) +
377 . abs(l2*btb33(i) + one/ntn))* vol(i) * wf(k,nddl)
378 ENDIF
379
380
381 ELSE
382 IF (nodadt > 0) THEN
383
384 f1(i,k) = wf(k,nddl)*sqrt(mass(pos1(i)+k-1)/mass0(pos1(i)+k-1))*zeta*sspnl*
385 . half*(vnl(pos1(i)+k-1)+vnl0(pos1(i)+k-1))*sqrt((four/sqrt(three))*(area0(i)))*thk0(i)
386 f2(i,k) = wf(k,nddl)*sqrt(mass(pos2(i)+k-1)/mass0(pos2(i)+k-1))*zeta*sspnl*
387 . half*(vnl(pos2(i)+k-1)+vnl0(pos2(i)+k-1))*sqrt((four/sqrt(three))*(area0(i)))*thk0(i)
388 f3(i,k) = wf(k,nddl)*sqrt(mass(pos3(i)+k-1)/mass0(pos3(i)+k-1))*zeta*sspnl*
389 . half*(vnl(pos3(i)+k-1)+vnl0(pos3(i)+k-1))*sqrt((four/sqrt(three))*(area0(i)))*thk0(i)
390
391 sti1(i,k) = em20
392 sti2(i,k) = em20
393 sti3(i,k) = em20
394 ELSE
395
396 f1(i,k) = wf(k,nddl)*zeta*sspnl*half*(vnl(pos1(i)+k-1)+vnl0(pos1(i)+k-1))*
397 . sqrt((four/sqrt(three))*(area0(i)))*thk0(i)
398 f2(i,k) = wf(k,nddl)*zeta*sspnl*half*(vnl(pos2(i)+k-1)+vnl0(pos2(i)+k-1))*
399 . sqrt((four/sqrt(three))*(area0(i)))*thk0(i)
400 f3(i,k) = wf(k,nddl)*zeta*sspnl*half*(vnl(pos3(i)+k-1)+vnl0(pos3(i)+k-1))*
401 . sqrt((four/sqrt(three))*(area0(i)))*thk0(i)
402 ENDIF
403 ENDIF
404 ENDDO
405 ENDDO
406
407
408
409
410
411
412 IF (iparit == 0) THEN
413
414 fnl => nloc_dmg%FNL(1:l_nloc,itask+1)
415 IF (nodadt > 0) stifnl => nloc_dmg%STIFNL(1:l_nloc,itask+1)
416
417 DO i=1,nel
418
419#include "vectorize.inc"
420 DO k = 1,nddl
421
422 fnl(pos1(i)+k-1) = fnl(pos1(i)+k-1) - f1(i,k)
423 fnl(pos2(i)+k-1) = fnl(pos2(i)+k-1) - f2(i,k)
424 fnl(pos3(i)+k-1) = fnl(pos3(i)+k-1) - f3(i,k)
425 IF (nodadt > 0) THEN
426
427 maxstif =
max(sti1(i,k),sti2(i,k),sti3(i,k))
428
429 stifnl(pos1(i)+k-1) = stifnl(pos1(i)+k-1) + maxstif
430 stifnl(pos2(i)+k-1) = stifnl(pos2(i)+k-1) + maxstif
431 stifnl(pos3(i)+k-1) = stifnl(pos3(i)+k-1) + maxstif
432 ENDIF
433 ENDDO
434 ENDDO
435
436
437 ELSE
438
439 DO j = 1,nddl
440
441
442 DO i=1,nel
443 ii = i + nft
444
445
446 IF (nodadt > 0) THEN
447 maxstif =
max(sti1(i,j),sti2(i,j),sti3(i,j))
448 ENDIF
449
450 k = nloc_dmg%IADTG(1,ii)
451 nloc_dmg%FSKY(k,j) = -f1(i,j)
452 IF (nodadt > 0) nloc_dmg%STSKY(k,j) = maxstif
453
454 k = nloc_dmg%IADTG(2,ii)
455 nloc_dmg%FSKY(k,j) = -f2(i,j)
456 IF (nodadt > 0) nloc_dmg%STSKY(k,j) = maxstif
457
458 k = nloc_dmg%IADTG(3,ii)
459 nloc_dmg%FSKY(k,j) = -f3(i,j)
460 IF (nodadt > 0) nloc_dmg%STSKY(k,j) = maxstif
461
462 ENDDO
463
464 ENDDO
465 ENDIF
466
467
468
469
470 IF (nodadt == 0) THEN
471 DO i = 1,nel
472
473 IF (off(i)/=zero) THEN
474
475 dtnl = (two*(
min(le(i),le_max))*sqrt(three*zeta))/
476 . sqrt(twelve*l2 + (
min(le(i),le_max))**2)
477
478 IF (nddl>1) THEN
479 IF (nddl > 2) THEN
480 dtnl_th = (two*(
min(thk(i)/nddl,le_max))*sqrt(three*zeta))/
481 . sqrt(twelve*l2 + (
min(thk(i)/nddl,le_max))**2)
482 ELSE
483 dtnl_th = (two*(
min(thk(i),le_max))*sqrt(three*zeta))/
484 . sqrt(twelve*l2 + (
min(thk(i),le_max))**2)
485 ENDIF
486 ELSE
487 dtnl_th = ep20
488 ENDIF
489
490 dt2t =
min(dt2t,dtfac1(1)*cdamp*dtnl_th,dtfac1(1)*cdamp*dtnl)
491 ENDIF
492 ENDDO
493 ENDIF
494
495
496 IF (ALLOCATED(f1)) DEALLOCATE(f1)
497 IF (ALLOCATED(f2)) DEALLOCATE(f2)
498 IF (ALLOCATED(f3)) DEALLOCATE(f3)
499 IF (ALLOCATED(sti1)) DEALLOCATE(sti1)
500 IF (ALLOCATED(sti2)) DEALLOCATE(sti2)
501 IF (ALLOCATED(sti3)) DEALLOCATE(sti3)
502 IF (ALLOCATED(btb11)) DEALLOCATE(btb11)
503 IF (ALLOCATED(btb12)) DEALLOCATE(btb12)
504 IF (ALLOCATED(btb13)) DEALLOCATE(btb13)
505 IF (ALLOCATED(btb22)) DEALLOCATE(btb22)
506 IF (ALLOCATED(btb23)) DEALLOCATE(btb23)
507 IF (ALLOCATED(btb33)) DEALLOCATE(btb33)
508 IF (ALLOCATED(stifnlth)) DEALLOCATE(stifnlth)
509 IF (ALLOCATED(pos1)) DEALLOCATE(pos1)
510 IF (ALLOCATED(pos2)) DEALLOCATE(pos2)
511 IF (ALLOCATED(pos3)) DEALLOCATE(pos3)
512 IF (ALLOCATED(vol)) DEALLOCATE(vol)
subroutine area(d1, x, x2, y, y2, eint, stif0)
for(i8=*sizetab-1;i8 >=0;i8--)