OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
cbafint_reg_ini.F File Reference
#include "implicit_f.inc"
#include "param_c.inc"
#include "com01_c.inc"
#include "com04_c.inc"
#include "scr03_c.inc"
#include "vectorize.inc"

Go to the source code of this file.

Functions/Subroutines

subroutine cbafint_reg_ini (elbuf_tab, nloc_dmg, area, ixc, dt_nl, x, xrefc, nft, nel, ng, ipm, bufmat, time, failure)

Function/Subroutine Documentation

◆ cbafint_reg_ini()

subroutine cbafint_reg_ini ( type(elbuf_struct_), dimension(ngroup), target elbuf_tab,
type (nlocal_str_), target nloc_dmg,
intent(in) area,
integer, dimension(nixc,*) ixc,
dt_nl,
x,
xrefc,
integer nft,
integer nel,
integer ng,
integer, dimension(npropmi,*) ipm,
bufmat,
time,
logical failure )

Definition at line 33 of file cbafint_reg_ini.F.

37C-----------------------------------------------
38C M o d u l e s
39C-----------------------------------------------
40 USE elbufdef_mod
41 USE message_mod
43 use element_mod , only : nixc
44C-----------------------------------------------
45C I m p l i c i t T y p e s
46C-----------------------------------------------
47#include "implicit_f.inc"
48C-----------------------------------------------
49C C o m m o n B l o c k s
50C-----------------------------------------------
51#include "param_c.inc"
52#include "com01_c.inc"
53#include "com04_c.inc"
54#include "scr03_c.inc"
55C-----------------------------------------------
56C D u m m y A r g u m e n t s
57C-----------------------------------------------
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) ::
62 . area
64 . x(3,*),xrefc(4,3,*),dt_nl,bufmat(*),time
65 LOGICAL :: FAILURE
66C-----------------------------------------------
67C L o c a l V a r i a b l e s
68C-----------------------------------------------
69 INTEGER :: IMAT,NDOF,L_NLOC,N1,N2,N3,N4,
70 . K,I,NPTR,NPTS,IR,IS,NDNOD
71 INTEGER, DIMENSION(NEL) :: POS1,POS2,POS3,POS4
73 . le_min,len,damp,dens,ntn_unl,ntn_vnl,
74 . ntvar,z01(11,11),wf1(11,11),zn1(12,11),b1,b2,
75 . b3,b4,nth1,nth2,bth1,bth2,k1,k12,k2,sspnl,le_max
76 my_real, DIMENSION(:,:), ALLOCATABLE :: vpred,var_reg
77 my_real,
78 . DIMENSION(:), POINTER :: fnl,unl,vnl,dnl,mnl,thck
79 my_real, DIMENSION(NEL) :: x1,x2,x3,x4,
80 . y1,y2,y3,y4,px1,px2,py1,py2,e1x,e2x,e3x,
81 . e1y,e2y,e3y,e1z,e2z,e3z,x2l,y2l,x3l,y3l,
82 . x4l,y4l,z1,z2,z3,z4,surf,offg,vols,btb11,
83 . btb12,btb22
84 TYPE(BUF_NLOC_),POINTER :: BUFNL
85 my_real, DIMENSION(:,:), POINTER ::
86 . massth,fnlth,vnlth,unlth
87 ! Position of integration points in the thickness
88 DATA z01/
89 1 0. ,0. ,0. ,0. ,0. ,
90 1 0. ,0. ,0. ,0. ,0. ,0. ,
91 2 -.5 ,0.5 ,0. ,0. ,0. ,
92 2 0. ,0. ,0. ,0. ,0. ,0. ,
93 3 -.5 ,0. ,0.5 ,0. ,0. ,
94 3 0. ,0. ,0. ,0. ,0. ,0. ,
95 4 -.5 ,-.1666667,0.1666667,0.5 ,0. ,
96 4 0. ,0. ,0. ,0. ,0. ,0. ,
97 5 -.5 ,-.25 ,0. ,0.25 ,0.5 ,
98 5 0. ,0. ,0. ,0. ,0. ,0. ,
99 6 -.5 ,-.3 ,-.1 ,0.1 ,0.3 ,
100 6 0.5 ,0. ,0. ,0. ,0. ,0. ,
101 7 -.5 ,-.3333333,-.1666667,0.0 ,0.1666667,
102 7 0.3333333,0.5 ,0. ,0. ,0. ,0. ,
103 8 -.5 ,-.3571429,-.2142857,-.0714286,0.0714286,
104 8 0.2142857,0.3571429,0.5 ,0. ,0. ,0. ,
105 9 -.5 ,-.375 ,-.25 ,-.125 ,0.0 ,
106 9 0.125 ,0.25 ,0.375 ,0.5 ,0. ,0. ,
107 a -.5 ,-.3888889,-.2777778,-.1666667,-.0555555,
108 a 0.0555555,0.1666667,0.2777778,0.3888889,0.5 ,0. ,
109 b -.5 ,-.4 ,-.3 ,-.2 ,-.1 ,
110 b 0. ,0.1 ,0.2 ,0.3 ,0.4 ,0.5 /
111 ! Weight of integration in the thickness
112 DATA wf1/
113 1 1. ,0. ,0. ,0. ,0. ,
114 1 0. ,0. ,0. ,0. ,0. ,0. ,
115 2 0.5 ,0.5 ,0. ,0. ,0. ,
116 2 0. ,0. ,0. ,0. ,0. ,0. ,
117 3 0.25 ,0.5 ,0.25 ,0. ,0. ,
118 3 0. ,0. ,0. ,0. ,0. ,0. ,
119 4 0.1666667,0.3333333,0.3333333,0.1666667,0. ,
120 4 0. ,0. ,0. ,0. ,0. ,0. ,
121 5 0.125 ,0.25 ,0.25 ,0.25 ,0.125 ,
122 5 0. ,0. ,0. ,0. ,0. ,0. ,
123 6 0.1 ,0.2 ,0.2 ,0.2 ,0.2 ,
124 6 0.1 ,0. ,0. ,0. ,0. ,0. ,
125 7 0.0833333,0.1666667,0.1666667,0.1666667,0.1666667,
126 7 0.1666667,0.0833333,0. ,0. ,0. ,0. ,
127 8 0.0714286,0.1428571,0.1428571,0.1428571,0.1428571,
128 8 0.1428571,0.1428571,0.0714286,0. ,0. ,0. ,
129 9 0.0625 ,0.125 ,0.125 ,0.125 ,0.125 ,
130 9 0.125 ,0.125 ,0.125 ,0.0625 ,0. ,0. ,
131 a 0.0555556,0.1111111,0.1111111,0.1111111,0.1111111,
132 a 0.1111111,0.1111111,0.1111111,0.1111111,0.0555556,0. ,
133 b 0.05 ,0.1 ,0.1 ,0.1 ,0.1 ,
134 b 0.1 ,0.1 ,0.1 ,0.1 ,0.1 ,0.05 /
135 ! Position of nodes in the shell thickness
136 DATA zn1/
137 1 0. ,0. ,0. ,0. ,0. ,0. ,
138 1 0. ,0. ,0. ,0. ,0. ,0. ,
139 2 -.5 ,0.5 ,0. ,0. ,0. ,0. ,
140 2 0. ,0. ,0. ,0. ,0. ,0. ,
141 3 -.5 ,-.25 ,0.25 ,0.5 ,0. ,0. ,
142 3 0. ,0. ,0. ,0. ,0. ,0. ,
143 4 -.5 ,-.3333333,0. ,0.3333333,0.5 ,0. ,
144 4 0. ,0. ,0. ,0. ,0. ,0. ,
145 5 -.5 ,-.375 ,-0.125 ,0.125 ,0.375 ,0.5 ,
146 5 0. ,0. ,0. ,0. ,0. ,0. ,
147 6 -.5 ,-.4 ,-.2 ,0.0 ,0.2 ,0.4 ,
148 6 0.5 ,0. ,0. ,0. ,0. ,0. ,
149 7 -.5 ,-.4166667,-.25 ,-.0833333,0.0833333,0.25 ,
150 7 0.4166667,0.5 ,0. ,0. ,0. ,0. ,
151 8 -.5 ,-.4285715,-.2857143,-.1428572,0.0 ,0.1428572,
152 8 0.2857143,0.4285715,0.5 ,0. ,0. ,0. ,
153 9 -.5 ,-.4375 ,-.3125 ,-.1875 ,-.0625 ,0.0625 ,
154 9 0.1875 ,0.3125 ,0.4375 ,0.5 ,0. ,0. ,
155 a -.5 ,-.4444444,-.3333333,-.2222222,-.1111111,0. ,
156 a 0.1111111,0.2222222,0.3333333,0.4444444,0.5 ,0. ,
157 b -.5 ,-.45 ,-.35 ,-.25 ,-.15 ,-.05 ,
158 b 0.05 ,0.15 ,0.25 ,0.35 ,0.45 ,0.5 /
159C=======================================================================
160 ! Size of the non-local vectors
161 l_nloc = nloc_dmg%L_NLOC
162 ! Pointing the non-local forces vector
163 fnl => nloc_dmg%FNL(1:l_nloc,1)
164 vnl => nloc_dmg%VNL(1:l_nloc)
165 dnl => nloc_dmg%DNL(1:l_nloc)
166 unl => nloc_dmg%UNL(1:l_nloc)
167 mnl => nloc_dmg%MASS(1:l_nloc)
168 ! Number of the material law
169 imat = ixc(1,1+nft)
170 ! Minimal length in the surface
171 le_min = sqrt(minval(area(nft+1:nft+nel)))
172 ! Number of integration points
173 ndof = elbuf_tab(ng)%BUFLY(1)%NPTT
174 ! Thickness of the shell
175 thck => elbuf_tab(ng)%GBUF%THK(1:nel)
176 ! Global minimal length
177 IF (ndof>1) THEN
178 le_min = min(le_min,minval(thck(1:nel))/ndof)
179 ENDIF
180 ! Non-local internal length
181 len = nloc_dmg%LEN(imat)
182 ! Maximal length of convergence
183 le_max = nloc_dmg%LE_MAX(imat)
184 ! Non-local damping
185 damp = nloc_dmg%DAMP(imat)
186 ! Non-local density
187 dens = nloc_dmg%DENS(imat)
188 ! Non-local soundspeed
189 sspnl = nloc_dmg%SSPNL(imat)
190 ! Non-local timestep
191 dt_nl = min(dt_nl,0.5d0*((two*min(le_min,le_max)*sqrt(three*dens))/
192 . (sqrt(twelve*(len**2)+(min(le_min,le_max)**2)))))
193 ! Allocation of the velocities predictor
194 IF (ndof>1) THEN
195 IF (ndof > 2) THEN
196 ALLOCATE(vpred(nel,ndof+1))
197 ndnod = ndof+1
198 ELSE
199 ALLOCATE(vpred(nel,ndof))
200 ndnod = ndof
201 ENDIF
202 ENDIF
203 ! Number of integration points in the surface of the element
204 nptr = elbuf_tab(ng)%NPTR
205 npts = elbuf_tab(ng)%NPTS
206 ! Variable to regularize
207 IF (.NOT.ALLOCATED(var_reg)) ALLOCATE(var_reg(nel,ndof))
208c
209 ! Computation of kinematical data
210# include "vectorize.inc"
211 DO i = 1,nel
212 ! Coordinates of the nodes of the element
213 IF (nxref == 0) THEN
214 x1(i)=x(1,ixc(2,nft+i))
215 y1(i)=x(2,ixc(2,nft+i))
216 z1(i)=x(3,ixc(2,nft+i))
217 x2(i)=x(1,ixc(3,nft+i))
218 y2(i)=x(2,ixc(3,nft+i))
219 z2(i)=x(3,ixc(3,nft+i))
220 x3(i)=x(1,ixc(4,nft+i))
221 y3(i)=x(2,ixc(4,nft+i))
222 z3(i)=x(3,ixc(4,nft+i))
223 x4(i)=x(1,ixc(5,nft+i))
224 y4(i)=x(2,ixc(5,nft+i))
225 z4(i)=x(3,ixc(5,nft+i))
226 ELSE
227 x1(i)=xrefc(1,1,nft+i)
228 y1(i)=xrefc(1,2,nft+i)
229 z1(i)=xrefc(1,3,nft+i)
230 x2(i)=xrefc(2,1,nft+i)
231 y2(i)=xrefc(2,2,nft+i)
232 z2(i)=xrefc(2,3,nft+i)
233 x3(i)=xrefc(3,1,nft+i)
234 y3(i)=xrefc(3,2,nft+i)
235 z3(i)=xrefc(3,3,nft+i)
236 x4(i)=xrefc(4,1,nft+i)
237 y4(i)=xrefc(4,2,nft+i)
238 z4(i)=xrefc(4,3,nft+i)
239 ENDIF
240c
241 ! Recovering the nodes of the non-local element
242 n1 = nloc_dmg%IDXI(ixc(2,nft+i))
243 n2 = nloc_dmg%IDXI(ixc(3,nft+i))
244 n3 = nloc_dmg%IDXI(ixc(4,nft+i))
245 n4 = nloc_dmg%IDXI(ixc(5,nft+i))
246 ! Recovering the positions of the first d.o.fs of each nodes
247 pos1(i) = nloc_dmg%POSI(n1)
248 pos2(i) = nloc_dmg%POSI(n2)
249 pos3(i) = nloc_dmg%POSI(n3)
250 pos4(i) = nloc_dmg%POSI(n4)
251 ENDDO
252c
253 ! Non-local variable transfer at Gauss point
254 ! Loop over Gauss points in the thickness
255 DO k = 1,ndof
256 ! Loop over element
257# include "vectorize.inc"
258 DO i = 1,nel
259 var_reg(i,k) = fourth*(dnl(pos1(i)+k-1) + dnl(pos2(i)+k-1)
260 . + dnl(pos3(i)+k-1) + dnl(pos4(i)+k-1))
261 ENDDO
262 ENDDO
263c
264 CALL cneveci(1 ,nel ,surf,
265 . x1 ,x2 ,x3 ,x4 ,y1 ,y2 ,
266 . y3 ,y4 ,z1 ,z2 ,z3 ,z4 ,
267 . e1x ,e2x ,e3x ,e1y ,e2y ,e3y ,
268 . e1z ,e2z ,e3z )
269C
270 ! Filling internal variable data of the non-local material
271 CALL cnloc_matini(elbuf_tab(ng),nel ,ipm ,
272 . bufmat ,time ,var_reg ,
273 . failure )
274C
275C
276 !-----------------------------------------------------------------------
277 ! Computation of the element volume and the BtB matrix product
278 !-----------------------------------------------------------------------
279 ! Loop over elements
280# include "vectorize.inc"
281 DO i=1,nel
282c
283 ! Computation of shape functions derivatives
284 x2l(i)=e1x(i)*(x2(i)-x1(i))+e1y(i)*(y2(i)-y1(i))+e1z(i)*(z2(i)-z1(i))
285 y2l(i)=e2x(i)*(x2(i)-x1(i))+e2y(i)*(y2(i)-y1(i))+e2z(i)*(z2(i)-z1(i))
286 x3l(i)=e1x(i)*(x3(i)-x1(i))+e1y(i)*(y3(i)-y1(i))+e1z(i)*(z3(i)-z1(i))
287 y3l(i)=e2x(i)*(x3(i)-x1(i))+e2y(i)*(y3(i)-y1(i))+e2z(i)*(z3(i)-z1(i))
288 x4l(i)=e1x(i)*(x4(i)-x1(i))+e1y(i)*(y4(i)-y1(i))+e1z(i)*(z4(i)-z1(i))
289 y4l(i)=e2x(i)*(x4(i)-x1(i))+e2y(i)*(y4(i)-y1(i))+e2z(i)*(z4(i)-z1(i))
290c
291 px1(i) = half *(y2l(i)-y4l(i))
292 px2(i) = half * y3l(i)
293 py1(i) = -half *(x2l(i)-x4l(i))
294 py2(i) = -half * x3l(i)
295c
296 ! Computation of the product BtxB
297 btb11(i) = px1(i)**2 + py1(i)**2
298 btb12(i) = px1(i)*px2(i) + py1(i)*py2(i)
299 btb22(i) = px2(i)**2 + py2(i)**2
300c
301 ! Computation of the element volume
302 vols(i) = area(nft+i)*thck(i)
303c
304 ! To check if element is not broken
305 offg(i) = elbuf_tab(ng)%GBUF%OFF(i)
306c
307 ENDDO
308C
309 !-----------------------------------------------------------------------
310 ! Pre-treatment non-local regularization in the thickness
311 !-----------------------------------------------------------------------
312 ! Only if NDOF > 1
313 IF ((ndof > 1).AND.(len>zero)) THEN
314c
315 ! Pointing the non-local values in the thickness of the corresponding element
316 bufnl => elbuf_tab(ng)%NLOC(1,1)
317c
318 ! Pointing the non-local values in the thickness of the corresponding element
319 massth => bufnl%MASSTH(1:nel,1:ndnod)
320 unlth => bufnl%UNLTH(1:nel ,1:ndnod)
321 vnlth => bufnl%VNLTH(1:nel ,1:ndnod)
322 fnlth => bufnl%FNLTH(1:nel ,1:ndnod)
323c
324 DO k = 1,ndnod
325 DO i = 1,nel
326 ! Prediction of the velocities
327 vpred(i,k) = vnlth(i,k) - (fnlth(i,k)/massth(i,k))*(dt_nl/two)
328 ENDDO
329 ENDDO
330 DO k = 1,ndnod
331 DO i = 1,nel
332 ! Resetting non-local forces
333 fnlth(i,k) = zero
334 ENDDO
335 ENDDO
336c
337 ! Computation of non-local forces in the shell thickness
338 DO k = 1, ndof
339c
340 ! Computation of shape functions value
341 IF ((ndof==2).AND.(k==2)) THEN
342 nth1 = (z01(k,ndof) - zn1(k,ndof))/(zn1(k-1,ndof) - zn1(k,ndof))
343 nth2 = (z01(k,ndof) - zn1(k-1,ndof)) /(zn1(k,ndof) - zn1(k-1,ndof))
344 ELSE
345 nth1 = (z01(k,ndof) - zn1(k+1,ndof))/(zn1(k,ndof) - zn1(k+1,ndof))
346 nth2 = (z01(k,ndof) - zn1(k,ndof)) /(zn1(k+1,ndof) - zn1(k,ndof))
347 ENDIF
348c
349 ! Loop over elements
350 DO i = 1,nel
351 ! Computation of B-matrix values
352 IF ((ndof==2).AND.(k==2)) THEN
353 bth1 = (one/(zn1(k-1,ndof) - zn1(k,ndof)))*(one/thck(i))
354 bth2 = (one/(zn1(k,ndof) - zn1(k-1,ndof)))*(one/thck(i))
355 ELSE
356 bth1 = (one/(zn1(k,ndof) - zn1(k+1,ndof)))*(one/thck(i))
357 bth2 = (one/(zn1(k+1,ndof) - zn1(k,ndof)))*(one/thck(i))
358 ENDIF
359c
360 ! Computation of the non-local K matrix
361 k1 = (len**2)*(bth1**2) + nth1**2
362 k12 = (len**2)*(bth1*bth2)+ (nth1*nth2)
363 k2 = (len**2)*(bth2**2) + nth2**2
364c
365 ! Computation of the non-local forces
366 IF ((ndof==2).AND.(k==2)) THEN
367 fnlth(i,k-1) = fnlth(i,k-1) + (k1*unlth(i,k-1) + k12*unlth(i,k)
368 . + damp*((nth1**2)*vpred(i,k-1)
369 . + (nth1*nth2)*vpred(i,k))
370 . - (nth1*var_reg(i,k)))*vols(i)*wf1(k,ndof)
371 fnlth(i,k) = fnlth(i,k) + (k12*unlth(i,k-1) + k2*unlth(i,k)
372 . + damp*(nth1*nth2*vpred(i,k-1)
373 . + (nth2**2)*vpred(i,k))
374 . - nth2*var_reg(i,k))*vols(i)*wf1(k,ndof)
375 ELSE
376 fnlth(i,k) = fnlth(i,k) + (k1*unlth(i,k) + k12*unlth(i,k+1)
377 . + damp*((nth1**2)*vpred(i,k)
378 . + (nth1*nth2)*vpred(i,k+1))
379 . - (nth1*var_reg(i,k)))*vols(i)*wf1(k,ndof)
380 fnlth(i,k+1) = fnlth(i,k+1) + (k12*unlth(i,k) + k2*unlth(i,k+1)
381 . + damp*(nth1*nth2*vpred(i,k)
382 . + (nth2**2)*vpred(i,k+1))
383 . - nth2*var_reg(i,k))*vols(i)*wf1(k,ndof)
384 ENDIF
385 ENDDO
386 ENDDO
387c
388 DO k = 1,ndnod
389 DO i = 1,nel
390 ! Updating the non-local in-thickness velocities
391 vnlth(i,k) = vnlth(i,k) - (fnlth(i,k)/massth(i,k))*dt_nl
392 ENDDO
393 ENDDO
394c
395 DO k = 1,ndnod
396 DO i = 1,nel
397 ! Computing the non-local in-thickness cumulated values
398 unlth(i,k) = unlth(i,k) + vnlth(i,k)*dt_nl
399 ENDDO
400 ENDDO
401c
402 ! Transfer at the integration point
403 DO k = 1, ndof
404 !Computation of shape functions value
405 IF ((ndof==2).AND.(k==2)) THEN
406 nth1 = (z01(k,ndof) - zn1(k,ndof))/(zn1(k-1,ndof) - zn1(k,ndof))
407 nth2 = (z01(k,ndof) - zn1(k-1,ndof)) /(zn1(k,ndof) - zn1(k-1,ndof))
408 ELSE
409 nth1 = (z01(k,ndof) - zn1(k+1,ndof))/(zn1(k,ndof) - zn1(k+1,ndof))
410 nth2 = (z01(k,ndof) - zn1(k,ndof)) /(zn1(k+1,ndof) - zn1(k,ndof))
411 ENDIF
412 ! Loop over elements
413 DO i = 1,nel
414 !Integration points non-local variables
415 IF ((ndof==2).AND.(k==2)) THEN
416 var_reg(i,k) = nth1*unlth(i,k-1) + nth2*unlth(i,k)
417 ELSE
418 var_reg(i,k) = nth1*unlth(i,k) + nth2*unlth(i,k+1)
419 ENDIF
420 ENDDO
421 ENDDO
422c
423 IF ((nptr>1).OR.(npts>1)) THEN
424 DO ir = 1,nptr
425 DO is = 1,npts
426 bufnl => elbuf_tab(ng)%NLOC(ir,is)
427 bufnl%MASSTH(1:nel,1:ndnod) = elbuf_tab(ng)%NLOC(1,1)%MASSTH(1:nel,1:ndnod)
428 bufnl%MASSTH(1:nel,1:ndnod) = elbuf_tab(ng)%NLOC(1,1)%UNLTH(1:nel,1:ndnod)
429 bufnl%MASSTH(1:nel,1:ndnod) = elbuf_tab(ng)%NLOC(1,1)%VNLTH(1:nel,1:ndnod)
430 bufnl%MASSTH(1:nel,1:ndnod) = elbuf_tab(ng)%NLOC(1,1)%FNLTH(1:nel,1:ndnod)
431 ENDDO
432 ENDDO
433 ENDIF
434 ENDIF
435c
436 !-----------------------------------------------------------------------
437 ! Computation of the elementary non-local forces
438 !-----------------------------------------------------------------------
439 ! Loop over integration points in the thickness
440 DO k = 1,ndof
441c
442 ! Computation of non-local forces
443# include "vectorize.inc"
444 DO i = 1,nel
445c
446 ! Computing the elementary non-local forces
447 IF (offg(i) > zero) THEN
448 ! Computing the product BtB*UNL
449 b1 = ((len**2)/vols(i))*wf1(k,ndof)*(btb11(i)*unl(pos1(i)+k-1) + btb12(i)*unl(pos2(i)+k-1)
450 . - btb11(i)*unl(pos3(i)+k-1) - btb12(i)*unl(pos4(i)+k-1))
451 b2 = ((len**2)/vols(i))*wf1(k,ndof)*(btb12(i)*unl(pos1(i)+k-1) + btb22(i)*unl(pos2(i)+k-1)
452 . - btb12(i)*unl(pos3(i)+k-1) - btb22(i)*unl(pos4(i)+k-1))
453 b3 = ((len**2)/vols(i))*wf1(k,ndof)*(-btb11(i)*unl(pos1(i)+k-1) - btb12(i)*unl(pos2(i)+k-1)
454 . + btb11(i)*unl(pos3(i)+k-1) + btb12(i)*unl(pos4(i)+k-1))
455 b4 = ((len**2)/vols(i))*wf1(k,ndof)*(-btb12(i)*unl(pos1(i)+k-1) - btb22(i)*unl(pos2(i)+k-1)
456 . + btb12(i)*unl(pos3(i)+k-1) + btb22(i)*unl(pos4(i)+k-1))
457 ! Computing the product NtN*UNL
458 ntn_unl = ((unl(pos1(i)+k-1) + unl(pos2(i)+k-1) + unl(pos3(i)+k-1) +
459 . unl(pos4(i)+k-1))*fourth*fourth)*vols(i)*wf1(k,ndof)
460 ! Computing the product DAMP*NtN*VNL!
461 ntn_vnl = ((vnl(pos1(i)+k-1) + vnl(pos2(i)+k-1) + vnl(pos3(i)+k-1) +
462 . vnl(pos4(i)+k-1))*fourth*fourth)*damp*vols(i)*wf1(k,ndof)
463 ! Introducing the internal variable to be regularized
464 ntvar = var_reg(i,k)*fourth*vols(i)*wf1(k,ndof)
465 !Computing the elementary non-local forces
466 fnl(pos1(i)+k-1) = fnl(pos1(i)+k-1) - (ntn_unl + ntn_vnl - ntvar + b1)
467 fnl(pos2(i)+k-1) = fnl(pos2(i)+k-1) - (ntn_unl + ntn_vnl - ntvar + b2)
468 fnl(pos3(i)+k-1) = fnl(pos3(i)+k-1) - (ntn_unl + ntn_vnl - ntvar + b3)
469 fnl(pos4(i)+k-1) = fnl(pos4(i)+k-1) - (ntn_unl + ntn_vnl - ntvar + b4)
470 ! If the element is broken, the non-local wave is absorbed
471 ELSE
472 ! Non-local absorbing forces
473 fnl(pos1(i)+k-1) = fnl(pos1(i)+k-1) - wf1(k,ndof)*dens*sspnl*vnl(pos1(i)+k-1)*le_max*thck(i)
474 fnl(pos2(i)+k-1) = fnl(pos2(i)+k-1) - wf1(k,ndof)*dens*sspnl*vnl(pos2(i)+k-1)*le_max*thck(i)
475 fnl(pos3(i)+k-1) = fnl(pos3(i)+k-1) - wf1(k,ndof)*dens*sspnl*vnl(pos3(i)+k-1)*le_max*thck(i)
476 fnl(pos4(i)+k-1) = fnl(pos4(i)+k-1) - wf1(k,ndof)*dens*sspnl*vnl(pos4(i)+k-1)*le_max*thck(i)
477 ENDIF
478 ENDDO
479 ENDDO
480c -------------------
481 IF (ALLOCATED(var_reg)) DEALLOCATE(var_reg)
482 IF (ALLOCATED(vpred)) DEALLOCATE(vpred)
subroutine cneveci(jft, jlt, area, x1, x2, x3, x4, y1, y2, y3, y4, z1, z2, z3, z4, e1x, e2x, e3x, e1y, e2y, e3y, e1z, e2z, e3z)
Definition cneveci.F:36
subroutine cnloc_matini(elbuf_str, nel, ipm, bufmat, time, varnl, failure)
#define my_real
Definition cppsort.cpp:32
subroutine area(d1, x, x2, y, y2, eint, stif0)
#define min(a, b)
Definition macros.h:20