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
43C-----------------------------------------------
44C I m p l i c i t T y p e s
45C-----------------------------------------------
46#include "implicit_f.inc"
47C-----------------------------------------------
48C C o m m o n B l o c k s
49C-----------------------------------------------
50#include "param_c.inc"
51#include "com01_c.inc"
52#include "com04_c.inc"
53#include "scr03_c.inc"
54C-----------------------------------------------
55C D u m m y A r g u m e n t s
56C-----------------------------------------------
57 TYPE(elbuf_struct_), TARGET, DIMENSION(NGROUP) :: ELBUF_TAB
58 TYPE (NLOCAL_STR_) , TARGET :: NLOC_DMG
59 INTEGER IXC(NIXC,*),NFT,NEL,NG,IPM(NPROPMI,*)
60 my_real ,DIMENSION(NUMELC+NUMELTG),INTENT(IN) ::
61 . area
63 . x(3,*),xrefc(4,3,*),dt_nl,bufmat(*),time
64 LOGICAL :: FAILURE
65C-----------------------------------------------
66C L o c a l V a r i a b l e s
67C-----------------------------------------------
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
76 my_real,
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 ! Position of integration points in the thickness
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 ! Weight of integration in the thickness
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 ! Position of nodes in the shell thickness
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 /
158C=======================================================================
159 ! Size of the non-local vectors
160 l_nloc = nloc_dmg%L_NLOC
161 ! Pointing the non-local forces vector
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 ! Number of the material law
168 imat = ixc(1,1+nft)
169 ! Minimal length in the surface
170 le_min = sqrt(minval(area(nft+1:nft+nel)))
171 ! Number of integration points
172 ndof = elbuf_tab(ng)%BUFLY(1)%NPTT
173 ! Thickness of the shell
174 thck => elbuf_tab(ng)%GBUF%THK(1:nel)
175 ! Global minimal length
176 IF (ndof>1) THEN
177 le_min = min(le_min,minval(thck(1:nel))/ndof)
178 ENDIF
179 ! Non-local internal length
180 len = nloc_dmg%LEN(imat)
181 ! Maximal length of convergence
182 le_max = nloc_dmg%LE_MAX(imat)
183 ! Non-local damping
184 damp = nloc_dmg%DAMP(imat)
185 ! Non-local density
186 dens = nloc_dmg%DENS(imat)
187 ! Non-local soundspeed
188 sspnl = nloc_dmg%SSPNL(imat)
189 ! Non-local timestep
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 ! Allocation of the velocities predictor
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 ! Number of integration points in the surface of the element
203 nptr = elbuf_tab(ng)%NPTR
204 npts = elbuf_tab(ng)%NPTS
205 ! Variable to regularize
206 IF (.NOT.ALLOCATED(var_reg)) ALLOCATE(var_reg(nel,ndof))
207c
208 ! Computation of kinematical data
209# include "vectorize.inc"
210 DO i = 1,nel
211 ! Coordinates of the nodes of the element
212 IF (nxref == 0) THEN
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))
225 ELSE
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)
238 ENDIF
239c
240 ! Recovering the nodes of the non-local element
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))
245 ! Recovering the positions of the first d.o.fs of each nodes
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)
250 ENDDO
251c
252 ! Non-local variable transfer at Gauss point
253 ! Loop over Gauss points in the thickness
254 DO k = 1,ndof
255 ! Loop over element
256# include "vectorize.inc"
257 DO i = 1,nel
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))
260 ENDDO
261 ENDDO
262c
263 CALL cneveci(1 ,nel ,surf,
264 . x1 ,x2 ,x3 ,x4 ,y1 ,y2 ,
265 . y3 ,y4 ,z1 ,z2 ,z3 ,z4 ,
266 . e1x ,e2x ,e3x ,e1y ,e2y ,e3y ,
267 . e1z ,e2z ,e3z )
268C
269 ! Filling internal variable data of the non-local material
270 CALL cnloc_matini(elbuf_tab(ng),nel ,ipm ,
271 . bufmat ,time ,var_reg ,
272 . failure )
273C
274C
275 !-----------------------------------------------------------------------
276 ! Computation of the element volume and the BtB matrix product
277 !-----------------------------------------------------------------------
278 ! Loop over elements
279# include "vectorize.inc"
280 DO i=1,nel
281c
282 ! Computation of shape functions derivatives
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)-z1(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(i)-y1(i))+e1z(i)*(z4(i)-z1(i))
288 y4l(i)=e2x(i)*(x4(i)-x1(i))+e2y(i)*(y4(i)-y1(i))+e2z(i)*(z4(i)-z1(i))
289c
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)
294c
295 ! Computation of the product BtxB
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
299c
300 ! Computation of the element volume
301 vols(i) = area(nft+i)*thck(i)
302c
303 ! To check if element is not broken
304 offg(i) = elbuf_tab(ng)%GBUF%OFF(i)
305c
306 ENDDO
307C
308 !-----------------------------------------------------------------------
309 ! Pre-treatment non-local regularization in the thickness
310 !-----------------------------------------------------------------------
311 ! Only if NDOF > 1
312 IF ((ndof > 1).AND.(len>zero)) THEN
313c
314 ! Pointing the non-local values in the thickness of the corresponding element
315 bufnl => elbuf_tab(ng)%NLOC(1,1)
316c
317 ! Pointing the non-local values in the thickness of the corresponding element
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)
322c
323 DO k = 1,ndnod
324 DO i = 1,nel
325 ! prediction of the velocities
326 vpred(i,k) = vnlth(i,k) - (fnlth(i,k)/massth(i,k))*(dt_nl/two)
327 ENDDO
328 ENDDO
329 DO k = 1,ndnod
330 DO i = 1,nel
331 ! Resetting non-local forces
332 fnlth(i,k) = zero
333 ENDDO
334 ENDDO
335c
336 ! Computation of non-local forces in the shell thickness
337 DO k = 1, ndof
338c
339 ! Computation of shape functions value
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))
343 ELSE
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))
346 ENDIF
347c
348 ! Loop over elements
349 DO i = 1,nel
350 ! Computation of B-matrix values
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))
354 ELSE
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))
357 ENDIF
358c
359 ! Computation of the non-local K matrix
360 k1 = (len**2)*(bth1**2) + nth1**2
361 k12 = (len**2)*(bth1*bth2)+ (nth1*nth2)
362 k2 = (len**2)*(bth2**2) + nth2**2
363c
364 ! Computation of the non-local forces
365 IF ((ndof==2).AND.(k==2)) THEN
366 fnlth(i,k-1) = fnlth(i,k-1) + (k1*unlth(i,k-1) + k12*unlth(i,k)
367 . + damp*((nth1**2)*vpred(i,k-1)
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)
374 ELSE
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)
383 ENDIF
384 ENDDO
385 ENDDO
386c
387 DO k = 1,ndnod
388 DO i = 1,nel
389 ! Updating the non-local in-thickness velocities
390 vnlth(i,k) = vnlth(i,k) - (fnlth(i,k)/massth(i,k))*dt_nl
391 ENDDO
392 ENDDO
393c
394 DO k = 1,ndnod
395 DO i = 1,nel
396 ! Computing the non-local in-thickness cumulated values
397 unlth(i,k) = unlth(i,k) + vnlth(i,k)*dt_nl
398 ENDDO
399 ENDDO
400c
401 ! Transfer at the integration point
402 DO k = 1, ndof
403 !Computation of shape functions value
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))
407 ELSE
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))
410 ENDIF
411 ! Loop over elements
412 DO i = 1,nel
413 !Integration points non-local variables
414 IF ((ndof==2).AND.(k==2)) THEN
415 var_reg(i,k) = nth1*unlth(i,k-1) + nth2*unlth(i,k)
416 ELSE
417 var_reg(i,k) = nth1*unlth(i,k) + nth2*unlth(i,k+1)
418 ENDIF
419 ENDDO
420 ENDDO
421c
422 IF ((nptr>1).OR.(npts>1)) THEN
423 DO ir = 1,nptr
424 DO is = 1,npts
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(1,1)%VNLTH(1:nel,1:ndnod)
429 bufnl%MASSTH(1:nel,1:ndnod) = elbuf_tab(ng)%NLOC(1,1)%FNLTH(1:nel,1:ndnod)
430 ENDDO
431 ENDDO
432 ENDIF
433 ENDIF
434c
435 !-----------------------------------------------------------------------
436 ! Computation of the elementary non-local forces
437 !-----------------------------------------------------------------------
438 ! Loop over integration points in the thickness
439 DO k = 1,ndof
440c
441 ! Computation of non-local forces
442# include "vectorize.inc"
443 DO i = 1,nel
444c
445 ! Computing the elementary non-local forces
446 IF (offg(i) > zero) THEN
447 ! Computing the product BtB*UNL
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(i)*unl(pos4(i)+k-1))
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))
456 ! Computing the product NtN*UNL
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)
459 ! Computing the product DAMP*NtN*VNL!
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)
462 ! Introducing the internal variable to be regularized
463 ntvar = var_reg(i,k)*fourth*vols(i)*wf1(k,ndof)
464 !Computing the elementary non-local forces
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)
469 ! If the element is broken, the non-local wave is absorbed
470 ELSE
471 ! Non-local absorbing forces
472 fnl(pos1(i)+k-1) = fnl(pos1(i)+k-1) - wf1(k,ndof)*dens*sspnl*vnl(pos1(i)+k-1)*le_max*thck(i)
473 fnl(pos2(i)+k-1) = fnl(pos2(i)+k-1) - wf1(k,ndof)*dens*sspnl*vnl(pos2(i)+k-1)*le_max*thck(i)
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)
476 ENDIF
477 ENDDO
478 ENDDO
479c -------------------
480 IF (ALLOCATED(var_reg)) DEALLOCATE(var_reg)
481 IF (ALLOCATED(vpred)) DEALLOCATE(vpred)
482 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
end diagonal values have been computed in the(sparse) matrix id.SOL
subroutine area(d1, x, x2, y, y2, eint, stif0)
#define min(a, b)
Definition macros.h:20