34 . DT_NL ,X ,XREFC ,NFT ,
35 . NEL ,NG ,IPM ,BUFMAT ,
46#include "implicit_f.inc"
57 TYPE(elbuf_struct_),
TARGET,
DIMENSION(NGROUP) :: ELBUF_TAB
59INTEGER IXC(NIXC,*),NFT,NEL,NG,IPM(NPROPMI,*)
60 my_real ,
DIMENSION(NUMELC+NUMELTG)INTENT(IN)
63 . x(3,*),xrefc(4,3,*),dt_nl,bufmat(*),time
68 INTEGER :: IMAT,NDOF,L_NLOC,N1,N2,N3,N4,
69 . K,I,NPTR,NPTS,IR,IS,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 :: vpred,var_reg
77 .
DIMENSION(:),
POINTER :: fnl,unl,vnl,dnl,mnl,thck
78 my_real,
DIMENSION(NEL) :: x1,x2,x3,x4,
80 . e1y,e2y,e3y,e1z,e2z,e3z,x2l,y2l,x3l,y3l,
81 . x4l,y4l,z1,z2,z3,z4,surf,offg,vols,btb11,
83 TYPE(buf_nloc_),
POINTER :: BUFNL
84 my_real,
DIMENSION(:,:),
POINTER ::
85 . MASSTH,FNLTH,VNLTH,UNLTH
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. ,
94 4 -.5 ,-.1666667,0.1666667,0.5 ,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. ,
104 9 -.5 ,-.375 ,-.25 ,-.125 ,0.0 ,
105 9 0.125 ,0.25 ,0.375 ,0.5 ,0. ,0.
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 /
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. ,
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 /
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.0
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 /
160 l_nloc = nloc_dmg%L_NLOC
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)
172 ndof = elbuf_tab(ng)%BUFLY(1)%NPTT
174 thck => elbuf_tab(ng)%GBUF%THK(1:nel)
177 le_min =
min(le_min,minval(thck(1:nel))/ndof)
180 len = nloc_dmg%LEN(imat)
182 le_max = nloc_dmg%LE_MAX(imat)
184 damp = nloc_dmg%DAMP(imat)
186 dens = nloc_dmg%DENS(imat)
188 sspnl = nloc_dmg%SSPNL(imat)
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)))))
195 ALLOCATE(vpred(nel,ndof+1))
198 ALLOCATE(vpred(nel,ndof))
203 nptr = elbuf_tab(ng)%NPTR
204 npts = elbuf_tab(ng)%NPTS
206 IF (.NOT.
ALLOCATED(var_reg))
ALLOCATE(var_reg(nel,ndof))
209# include "vectorize.inc"
213 x1(i)=x(1,ixc(2,nft+i))
214 y1(i)=x(2,ixc(2,nft+i))
215 z1(i)=x(3,ixc(2,nft+i))
216 x2(i)=x(1,ixc(3,nft+i))
217 y2(i)=x(2,ixc(3,nft+i))
218 z2(i)=x(3,ixc(3,nft+i))
219 x3(i)=x(1,ixc(4,nft+i))
220 y3(i)=x(2,ixc(4,nft+i))
221 z3(i)=x(3,ixc(4,nft+i))
222 x4(i)=x(1,ixc(5,nft+i))
223 y4(i)=x(2,ixc(5,nft+i))
224 z4(i)=x(3,ixc(5,nft+i))
226 x1(i)=xrefc(1,1,nft+i)
227 y1(i)=xrefc(1,2,nft+i)
228 z1(i)=xrefc(1,3,nft+i)
229 x2(i)=xrefc(2,1,nft+i)
230 y2(i)=xrefc(2,2,nft+i)
231 z2(i)=xrefc(2,3,nft+i)
232 x3(i)=xrefc(3,1,nft+i)
233 y3(i)=xrefc(3,2,nft+i)
234 z3(i)=xrefc(3,3,nft+i)
235 x4(i)=xrefc(4,1,nft+i)
236 y4(i)=xrefc(4,2,nft+i)
237 z4(i)=xrefc(4,3,nft+i)
241 n1 = nloc_dmg%IDXI(ixc(2,nft+i))
242 n2 = nloc_dmg%IDXI(ixc(3,nft+i))
243 n3 = nloc_dmg%IDXI(ixc(4,nft+i))
244 n4 = nloc_dmg%IDXI(ixc(5,nft+i))
246 pos1(i) = nloc_dmg%POSI(n1)
247 pos2(i) = nloc_dmg%POSI(n2)
248 pos3(i) = nloc_dmg%POSI(n3)
249 pos4(i) = nloc_dmg%POSI(n4)
256# include "vectorize.inc"
258 var_reg(i,k) = fourth*(dnl(pos1(i)+k-1) + dnl(pos2(i)+k-1)
259 . + dnl(pos3(i)+k-1) + dnl(pos4(i)+k-1))
264 . x1 ,x2 ,x3 ,x4 ,y1 ,y2 ,
265 . y3 ,y4 ,z1 ,z2 ,z3 ,z4 ,
266 . e1x ,e2x ,e3x ,e1y ,e2y ,e3y ,
279# include "vectorize.inc"
283 x2l(i)=e1x(i)*(x2(i)-x1(i))+e1y(i)*(y2(i)-y1(i))+e1z(i)*(z2(i)-z1(i))
284 y2l(i)=e2x(i)*(x2(i)-x1(i))+e2y(i)*(y2(i)-y1(i))+e2z(i)*(z2(i)-z1(i))
285 x3l(i)=e1x(i)*(x3(i)-x1(i))+e1y(i)*(y3(i)-y1(i))+e1z(i)*(z3(i
286 y3l(i)=e2x(i)*(x3(i)-x1(i))+e2y(i)*(y3(i)-y1(i))+e2z(i)*(z3(i)-z1(i))
287 x4l(i)=e1x(i)*(x4(i)-x1(i))+e1y(i)*(y4
288 y4l(i)=e2x(i)*(x4(i)-x1(i))+e2y(i)*(y4(i)-y1(i))+e2z(i)*(z4(i)-z1(i))
290 px1(i) = half *(y2l(i)-y4l(i))
291 px2(i) = half * y3l(i)
292 py1(i) = -half *(x2l(i)-x4l(i))
293 py2(i) = -half * x3l(i)
296 btb11(i) = px1(i)**2 + py1(i)**2
297 btb12(i) = px1(i)*px2(i) + py1(i)*py2(i)
298 btb22(i) = px2(i)**2 + py2(i)**2
301 vols(i) =
area(nft+i)*thck(i)
304 offg(i) = elbuf_tab(ng)%GBUF%OFF(i)
312 IF ((ndof > 1).AND.(len>zero))
THEN
315 bufnl => elbuf_tab(ng)%NLOC(1,1)
318 massth => bufnl%MASSTH(1:nel,1:ndnod)
319 unlth => bufnl%UNLTH(1:nel ,1:ndnod)
320 vnlth => bufnl%VNLTH(1:nel ,1:ndnod)
321 fnlth => bufnl%FNLTH(1:nel ,1:ndnod)
325 ! prediction of
the velocities
326 vpred(i,k) = vnlth(i,k) - (fnlth(i,k)/massth(i,k))*(dt_nl/two)
340 IF ((ndof==2).AND.(k==2))
THEN
341 nth1 = (z01(k,ndof) - zn1(k,ndof))/(zn1(k-1,ndof) - zn1(k,ndof))
342 nth2 = (z01(k,ndof) - zn1(k-1,ndof)) /(zn1(k,ndof) - zn1(k-1,ndof))
344 nth1 = (z01(k,ndof) - zn1(k+1,ndof))/(zn1(k,ndof) - zn1(k+1,ndof))
345 nth2 = (z01(k,ndof) - zn1(k,ndof)) /(zn1(k+1,ndof) - zn1(k,ndof))
351 IF ((ndof==2).AND.(k==2))
THEN
352 bth1 = (one/(zn1(k-1,ndof) - zn1(k,ndof)))*(one/thck(i))
353 bth2 = (one/(zn1(k,ndof) - zn1(k-1,ndof)))*(one/thck(i))
355 bth1 = (one/(zn1(k,ndof) - zn1(k+1,ndof)))*(one/thck(i))
356 bth2 = (one/(zn1(k+1,ndof) - zn1(k,ndof)))*(one/thck(i))
360 k1 = (len**2)*(bth1**2) + nth1**2
361 k12 = (len**2)*(bth1*bth2)+ (nth1*nth2)
362 k2 = (len**2)*(bth2**2) + nth2**2
365 IF ((ndof==2).AND.(k==2))
THEN
366 fnlth(i,k-1) = fnlth(i,k-1) + (k1
368 . + (nth1*nth2)*vpred(i,k))
369 . - (nth1*var_reg(i,k)))*vols(i)*wf1(k,ndof)
370 fnlth(i,k) = fnlth(i,k) + (k12*unlth(i,k-1) + k2*unlth(i,k)
371 . + damp*(nth1*nth2*vpred(i,k-1)
372 . + (nth2**2)*vpred(i,k))
373 . - nth2*var_reg(i,k))*vols(i)*wf1(k,ndof)
375 fnlth(i,k) = fnlth(i,k) + (k1*unlth(i,k) + k12*unlth(i,k+1)
376 . + damp*((nth1**2)*vpred(i,k)
377 . + (nth1*nth2)*vpred(i,k+1))
378 . - (nth1*var_reg(i,k)))*vols(i)*wf1(k,ndof)
379 fnlth(i,k+1) = fnlth(i,k+1) + (k12*unlth(i,k) + k2*unlth(i,k+1)
380 . + damp*(nth1*nth2*vpred(i,k)
381 . + (nth2**2)*vpred(i,k+1))
382 . - nth2*var_reg(i,k))*vols(i)*wf1(k,ndof)
390 vnlth(i,k) = vnlth(i,k) - (fnlth(i,k)/massth(i,k))*dt_nl
397 unlth(i,k) = unlth(i,k) + vnlth(i,k)*dt_nl
404 IF ((ndof==2).AND.(k==2))
THEN
405 nth1 = (z01(k,ndof) - zn1(k,ndof))/(zn1(k-1,ndof) - zn1(k,ndof))
406 nth2 = (z01(k,ndof) - zn1(k-1,ndof)) /(zn1(k,ndof) - zn1(k-1,ndof))
408 nth1 = (z01(k,ndof) - zn1(k+1,ndof))/(zn1(k,ndof) - zn1(k+1,ndof))
409 nth2 = (z01(k,ndof) - zn1(k,ndof)) /(zn1(k+1,ndof) - zn1(k,ndof))
414 IF ((ndof==2).AND.(k==2))
THEN
415 var_reg(i,k) = nth1*unlth(i,k-1) + nth2*unlth(i,k)
417 var_reg(i,k) = nth1*unlth(i,k) + nth2*unlth(i,k+1)
422 IF ((nptr>1).OR.(npts>1))
THEN
425 bufnl => elbuf_tab(ng)%NLOC(ir,is)
426 bufnl%MASSTH(1:nel,1:ndnod) = elbuf_tab(ng)%NLOC(1,1)%MASSTH(1:nel,1:ndnod)
427 bufnl%MASSTH(1:nel,1:ndnod) = elbuf_tab(ng)%NLOC(1,1)%UNLTH(1:nel,1:ndnod)
428 bufnl%MASSTH(1:nel,1:ndnod) = elbuf_tab(ng)%NLOC
429 bufnl%MASSTH(1:nel,1:ndnod) = elbuf_tab(ng)%NLOC(1,1)%FNLTH(1:nel,1:ndnod)
442# include "vectorize.inc"
446 IF (offg(i) > zero)
THEN
448 b1 = ((len**2)/vols(i))*wf1(k,ndof)*(btb11(i)*unl(pos1(i)+k-1) + btb12(i)*unl(pos2(i)+k-1)
449 . - btb11(i)*unl(pos3(i)+k-1) - btb12(i)*unl(pos4(i)+k-1))
450 b2 = ((len**2)/vols(i))*wf1(k,ndof)*(btb12(i)*unl(pos1(i)+k-1) + btb22(i)*unl(pos2(i)+k-1)
451 . - btb12(i)*unl(pos3(i)+k-1) - btb22
452 b3 = ((len**2)/vols(i))*wf1(k,ndof)*(-btb11(i)*unl(pos1(i)+k-1) - btb12(i)*unl(pos2(i)+k-1)
453 . + btb11(i)*unl(pos3(i)+k-1) + btb12(i)*unl(pos4(i)+k-1))
454 b4 = ((len**2)/vols(i))*wf1(k,ndof)*(-btb12(i)*unl(pos1(i)+k-1) - btb22(i)*unl(pos2(i)+k-1)
455 . + btb12(i)*unl(pos3(i)+k-1) + btb22(i)*unl(pos4(i)+k-1))
457 ntn_unl = ((unl(pos1(i)+k-1) + unl(pos2(i)+k-1) + unl(pos3(i)+k-1) +
458 . unl(pos4(i)+k-1))*fourth*fourth)*vols(i)*wf1(k,ndof)
460 ntn_vnl = ((vnl(pos1(i)+k-1) + vnl(pos2(i)+k-1) + vnl(pos3(i)+k-1) +
461 . vnl(pos4(i)+k-1))*fourth*fourth)*damp*vols(i)*wf1(k,ndof
463 ntvar = var_reg(i,k)*fourth*vols(i)*wf1(k,ndof)
465 fnl(pos1(i)+k-1) = fnl(pos1(i)+k-1) - (ntn_unl + ntn_vnl - ntvar + b1)
466 fnl(pos2(i)+k-1) = fnl(pos2(i)+k-1) - (ntn_unl + ntn_vnl - ntvar + b2)
467 fnl(pos3(i)+k-1) = fnl(pos3(i)+k-1) - (ntn_unl + ntn_vnl - ntvar + b3
468 fnl(pos4(i)+k-1) = fnl(pos4(i)+k-1) - (ntn_unl + ntn_vnl - ntvar + b4)
472 fnl(pos1(i)+k-1) = fnl(pos1(i)+k-1) - wf1(k,ndof)*dens*sspnl*vnl(pos1(i)+k-1)*le_max*thck
473 fnl(pos2(i)+k-1) = fnl(pos2(i)+k-1) - wf1(k,ndof)*dens*sspnl*vnl(pos2(i)+k
474 fnl(pos3(i)+k-1) = fnl(pos3(i)+k-1) - wf1(k,ndof)*dens*sspnl*vnl(pos3(i)+k-1)*le_max*thck(i)
475 fnl(pos4(i)+k-1) = fnl(pos4(i)+k-1) - wf1(k,ndof)*dens*sspnl*vnl(pos4(i)+k-1)*le_max*thck(i)
480 IF (
ALLOCATED(var_reg))
DEALLOCATE(var_reg)
481 IF (
ALLOCATED(vpred))
DEALLOCATE(vpred)