OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
cbafint_reg_ini.F
Go to the documentation of this file.
1Copyright> OpenRadioss
2Copyright> Copyright (C) 1986-2025 Altair Engineering Inc.
3Copyright>
4Copyright> This program is free software: you can redistribute it and/or modify
5Copyright> it under the terms of the GNU Affero General Public License as published by
6Copyright> the Free Software Foundation, either version 3 of the License, or
7Copyright> (at your option) any later version.
8Copyright>
9Copyright> This program is distributed in the hope that it will be useful,
10Copyright> but WITHOUT ANY WARRANTY; without even the implied warranty of
11Copyright> MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12Copyright> GNU Affero General Public License for more details.
13Copyright>
14Copyright> You should have received a copy of the GNU Affero General Public License
15Copyright> along with this program. If not, see <https://www.gnu.org/licenses/>.
16Copyright>
17Copyright>
18Copyright> Commercial Alternative: Altair Radioss Software
19Copyright>
20Copyright> As an alternative to this open-source version, Altair also offers Altair Radioss
21Copyright> software under a commercial license. Contact Altair to discuss further if the
22Copyright> commercial version may interest you: https://www.altair.com/radioss/.
23!||====================================================================
24!|| cbafint_reg_ini ../starter/source/elements/shell/coqueba/cbafint_reg_ini.F
25!||--- called by ------------------------------------------------------
26!|| nlocal_init_sta ../starter/source/materials/fail/nlocal_init_sta.F
27!||--- calls -----------------------------------------------------
28!|| cneveci ../starter/source/elements/shell/coqueba/cneveci.F
29!|| cnloc_matini ../starter/source/materials/mat_share/cnloc_matini.F
30!||--- uses -----------------------------------------------------
31!|| message_mod ../starter/share/message_module/message_mod.F
32!||====================================================================
33 SUBROUTINE cbafint_reg_ini(ELBUF_TAB,NLOC_DMG ,AREA ,IXC ,
34 . DT_NL ,X ,XREFC ,NFT ,
35 . NEL ,NG ,IPM ,BUFMAT ,
36 . TIME ,FAILURE )
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)
483 END
subroutine cbafint_reg_ini(elbuf_tab, nloc_dmg, area, ixc, dt_nl, x, xrefc, nft, nel, ng, ipm, bufmat, time, failure)
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