OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
cfint_reg.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!|| cfint_reg ../engine/source/elements/shell/coque/cfint_reg.F
25!||--- called by ------------------------------------------------------
26!|| cforc3 ../engine/source/elements/shell/coque/cforc3.F
27!|| czforc3 ../engine/source/elements/shell/coquez/czforc3.F
28!||--- uses -----------------------------------------------------
29!|| elbufdef_mod ../common_source/modules/mat_elem/elbufdef_mod.F90
30!|| nlocal_reg_mod ../common_source/modules/nlocal_reg_mod.F
31!||====================================================================
32 SUBROUTINE cfint_reg(
33 1 NLOC_DMG,VAR_REG, THK, NEL,
34 2 OFF, AREA, NC1, NC2,
35 3 NC3, NC4, PX1, PY1,
36 4 PX2, PY2, BUFNL, IMAT,
37 5 NDDL, ITASK, DT2T, LE,
38 6 THK0, AREA0, NFT)
39C-----------------------------------------------
40C M o d u l e s
41C-----------------------------------------------
42 USE nlocal_reg_mod
43 USE elbufdef_mod
44C-----------------------------------------------
45C I m p l i c i t T y p e s
46C-----------------------------------------------
47#include "implicit_f.inc"
48C-----------------------------------------------
49C G l o b a l P a r a m e t e r s
50C-----------------------------------------------
51#include "parit_c.inc"
52#include "com20_c.inc"
53#include "com08_c.inc"
54#include "scr02_c.inc"
55#include "scr18_c.inc"
56C-----------------------------------------------
57C D u m m y A r g u m e n t s
58C-----------------------------------------------
59 INTEGER, INTENT(IN) :: NFT
60 INTEGER :: NEL,IMAT,NDDL,ITASK
61 INTEGER, DIMENSION(NEL) :: NC1,NC2,NC3,NC4
62 my_real, DIMENSION(NEL,NDDL), INTENT(INOUT)::
63 . VAR_REG
64 my_real, DIMENSION(NEL), INTENT(IN) ::
65 . area,off,px1,py1,px2,py2,thk,le,thk0,area0
66 my_real, INTENT(INOUT) ::
67 . dt2t
68 TYPE(nlocal_str_),TARGET :: NLOC_DMG
69 TYPE(BUF_NLOC_) ,TARGET :: BUFNL
70C-----------------------------------------------
71C L o c a l V a r i a b l e s
72C-----------------------------------------------
73 INTEGER I,II,K,N1,N2,N3,N4,L_NLOC,J,NDNOD
74 my_real
75 . dx, dy, dz, l2,ntn,ntn_unl,ntn_vnl,xi,ntvar,a,zeta,
76 . b1,b2,b3,b4,sspnl,nth1,nth2,bth1,bth2,k1,k2,k12,
77 . dtnl_th,dtnl,le_max,maxstif,dtnod,dt2p
78 my_real, DIMENSION(:,:), ALLOCATABLE ::
79 . f1,f2,f3,f4,sti1,sti2,sti3,sti4
80 my_real, DIMENSION(:) ,ALLOCATABLE ::
81 . btb11,btb12,btb22,vol
82 INTEGER, DIMENSION(:), ALLOCATABLE ::
83 . POS1,POS2,POS3,POS4
84 my_real, POINTER, DIMENSION(:) ::
85 . vnl,fnl,unl,stifnl,mass,mass0,vnl0
86 my_real, POINTER, DIMENSION(:,:) ::
87 . massth,unlth,vnlth,fnlth
88 my_real, DIMENSION(:,:), ALLOCATABLE ::
89 . stifnlth,dtn
90 ! Safety coefficient for non-local stability vs mechanical stability
91 ! (it has been slightly increased vs nloc_dmg_init.F)
92 my_real, PARAMETER :: csta = 40.0d0
93 ! Coefficient for non-local stability to take into account damping
94 my_real, PARAMETER :: cdamp = 0.7d0
95c-----------------------------------------------------------------------
96c VAR_REG : variable a regulariser (local cumulated plastic strain)
97c NTVAR = NT * VAR_REG
98C=======================================================================
99 ! Recovering non-local parameters
100 l2 = nloc_dmg%LEN(imat)**2 ! Non-local internal length ** 2
101 xi = nloc_dmg%DAMP(imat) ! Non-local damping parameter
102 zeta = nloc_dmg%DENS(imat) ! Non-local density
103 sspnl = nloc_dmg%SSPNL(imat) ! Non-local sound speed
104 l_nloc = nloc_dmg%L_NLOC ! Length of non-local tables
105 le_max = nloc_dmg%LE_MAX(imat) ! Maximal length of convergence
106 ! Allocation of elementary forces vectors
107 ALLOCATE(f1(nel,nddl),f2(nel,nddl),f3(nel,nddl),f4(nel,nddl))
108 ! Only for nodal timestep
109 IF (nodadt > 0) THEN
110 ! Nodal stiffness
111 ALLOCATE(sti1(nel,nddl),sti2(nel,nddl),sti3(nel,nddl),sti4(nel,nddl))
112 ! Non-local mass
113 mass => nloc_dmg%MASS(1:l_nloc)
114 ! Initial non-local mass
115 mass0 => nloc_dmg%MASS0(1:l_nloc)
116 ENDIF
117 ALLOCATE(btb11(nel),btb12(nel),btb22(nel),vol(nel),
118 . pos1(nel),pos2(nel),pos3(nel),pos4(nel))
119 ! Recovering non-local data
120 vnl => nloc_dmg%VNL(1:l_nloc) ! Non-local variable velocities
121 vnl0 => nloc_dmg%VNL_OLD(1:l_nloc) ! Non-local variable velocities
122 unl => nloc_dmg%UNL(1:l_nloc) ! Non-local cumulated variable
123 ! Constant for quad elements
124 ntn = four*four
125c
126 !-----------------------------------------------------------------------
127 ! Computation of the element volume and the BtB matrix product
128 !-----------------------------------------------------------------------
129 ! Loop over elements
130# include "vectorize.inc"
131 DO i=1,nel
132c
133 ! Recovering the nodes of the triangle element
134 n1 = nloc_dmg%IDXI(nc1(i))
135 n2 = nloc_dmg%IDXI(nc2(i))
136 n3 = nloc_dmg%IDXI(nc3(i))
137 n4 = nloc_dmg%IDXI(nc4(i))
138c
139 ! Recovering the positions of the first d.o.fs of each nodes
140 pos1(i) = nloc_dmg%POSI(n1)
141 pos2(i) = nloc_dmg%POSI(n2)
142 pos3(i) = nloc_dmg%POSI(n3)
143 pos4(i) = nloc_dmg%POSI(n4)
144c
145 ! Computation of the element volume
146 vol(i) = area(i)*thk(i)
147c
148 ! Computation of the product BtxB
149 btb11(i) = px1(i)**2 + py1(i)**2
150 btb12(i) = px1(i)*px2(i) + py1(i)*py2(i)
151 btb22(i) = px2(i)**2 + py2(i)**2
152c
153 ENDDO
154c
155 !-----------------------------------------------------------------------
156 ! Pre-treatment non-local regularization in the thickness
157 !-----------------------------------------------------------------------
158 ! Only if NDDL > 1
159 IF ((nddl > 1).AND.(l2>zero)) THEN
160c
161 ! Allocation of the velocities predictor
162 IF (nddl > 2) THEN
163 IF (nodadt > 0) THEN
164 ALLOCATE(stifnlth(nel,nddl+1))
165 ALLOCATE(dtn(nel,nddl+1))
166 ENDIF
167 ndnod = nddl+1
168 ELSE
169 IF (nodadt > 0) THEN
170 ALLOCATE(stifnlth(nel,nddl))
171 ALLOCATE(dtn(nel,nddl))
172 ENDIF
173 ndnod = nddl
174 ENDIF
175c
176 ! Pointing the non-local values in the thickness of the corresponding element
177 massth => bufnl%MASSTH(1:nel,1:ndnod)
178 unlth => bufnl%UNLTH(1:nel,1:ndnod)
179 vnlth => bufnl%VNLTH(1:nel,1:ndnod)
180 fnlth => bufnl%FNLTH(1:nel,1:ndnod)
181c
182 DO k = 1,ndnod
183 DO i = 1,nel
184 ! Resetting non-local forces
185 fnlth(i,k) = zero
186 ! Resetting non-local nodal stiffness
187 IF (nodadt > 0) THEN
188 stifnlth(i,k) = em20
189 ENDIF
190 ENDDO
191 ENDDO
192c
193 ! Computation of non-local forces in the shell thickness
194 DO k = 1, nddl
195c
196 ! Computation of shape functions value
197 IF ((nddl==2).AND.(k==2)) THEN
198 nth1 = (z0(k,nddl) - zth(k,nddl)) / (zth(k-1,nddl) - zth(k,nddl))
199 nth2 = (z0(k,nddl) - zth(k-1,nddl)) / (zth(k,nddl) - zth(k-1,nddl))
200 ELSE
201 nth1 = (z0(k,nddl) - zth(k+1,nddl)) / (zth(k,nddl) - zth(k+1,nddl))
202 nth2 = (z0(k,nddl) - zth(k,nddl)) / (zth(k+1,nddl) - zth(k,nddl))
203 ENDIF
204c
205 ! Loop over elements
206 DO i = 1,nel
207 ! Computation of B-matrix values
208 IF ((nddl==2).AND.(k==2)) THEN
209 bth1 = (one/(zth(k-1,nddl) - zth(k,nddl)))*(one/thk(i))
210 bth2 = (one/(zth(k,nddl) - zth(k-1,nddl)))*(one/thk(i))
211 ELSE
212 bth1 = (one/(zth(k,nddl) - zth(k+1,nddl)))*(one/thk(i))
213 bth2 = (one/(zth(k+1,nddl) - zth(k,nddl)))*(one/thk(i))
214 ENDIF
215c
216 ! Computation of the non-local K matrix
217 k1 = l2*(bth1**2) + nth1**2
218 k12 = l2*(bth1*bth2)+ (nth1*nth2)
219 k2 = l2*(bth2**2) + nth2**2
220c
221 ! Computation of the non-local forces
222 IF ((nddl==2).AND.(k==2)) THEN
223 fnlth(i,k-1) = fnlth(i,k-1) + (k1*unlth(i,k-1) + k12*unlth(i,k)
224 . + xi*((nth1**2)*vnlth(i,k-1)
225 . + (nth1*nth2)*vnlth(i,k))
226 . - (nth1*var_reg(i,k)))*vol(i)*wf(k,nddl)
227 fnlth(i,k) = fnlth(i,k) + (k12*unlth(i,k-1) + k2*unlth(i,k)
228 . + xi*(nth1*nth2*vnlth(i,k-1)
229 . + (nth2**2)*vnlth(i,k))
230 . - nth2*var_reg(i,k))*vol(i)*wf(k,nddl)
231 ELSE
232 fnlth(i,k) = fnlth(i,k) + (k1*unlth(i,k) + k12*unlth(i,k+1)
233 . + xi*((nth1**2)*vnlth(i,k)
234 . + (nth1*nth2)*vnlth(i,k+1))
235 . - (nth1*var_reg(i,k)))*vol(i)*wf(k,nddl)
236 fnlth(i,k+1) = fnlth(i,k+1) + (k12*unlth(i,k) + k2*unlth(i,k+1)
237 . + xi*(nth1*nth2*vnlth(i,k)
238 . + (nth2**2)*vnlth(i,k+1))
239 . - nth2*var_reg(i,k))*vol(i)*wf(k,nddl)
240 ENDIF
241c
242 ! Computation of non-local nodal stiffness
243 IF (nodadt > 0) THEN
244 IF ((nddl==2).AND.(k==2)) THEN
245 stifnlth(i,k-1) = stifnlth(i,k-1) + max(abs(k1)+abs(k12),abs(k12)+abs(k2))*vol(i)*wf(k,nddl)
246 stifnlth(i,k) = stifnlth(i,k) + max(abs(k1)+abs(k12),abs(k12)+abs(k2))*vol(i)*wf(k,nddl)
247 ELSE
248 stifnlth(i,k) = stifnlth(i,k) + max(abs(k1)+abs(k12),abs(k12)+abs(k2))*vol(i)*wf(k,nddl)
249 stifnlth(i,k+1) = stifnlth(i,k+1) + max(abs(k1)+abs(k12),abs(k12)+abs(k2))*vol(i)*wf(k,nddl)
250 ENDIF
251 ENDIF
252c
253 ENDDO
254 ENDDO
255c
256 ! Updating non-local mass with /DT/NODA
257 IF (nodadt > 0) THEN
258C
259 ! Initial computation of the nodal timestep
260 dtnod = ep20
261 DO k = 1,ndnod
262 DO i = 1,nel
263 dtn(i,k) = dtfac1(11)*cdamp*sqrt(two*massth(i,k)/max(stifnlth(i,k),em20))
264 dtnod = min(dtn(i,k),dtnod)
265 ENDDO
266 ENDDO
267C
268 ! /DT/NODA/CSTX - Constant timestep with added mass
269 IF ((idtmin(11)==3).OR.(idtmin(11)==4).OR.(idtmin(11)==8)) THEN
270 ! Added mass computation if necessary
271 IF (dtnod < dtmin1(11)*(sqrt(csta))) THEN
272 DO k = 1,ndnod
273 DO i = 1,nel
274 IF (dtn(i,k) < dtmin1(11)) THEN
275 dt2p = dtmin1(11)/(dtfac1(11)*cdamp)
276 massth(i,k) = max(massth(i,k),csta*half*stifnlth(i,k)*dt2p*dt2p*onep00001)
277 ENDIF
278 ENDDO
279 ENDDO
280 ENDIF
281 dtnod = dtmin1(11)*(sqrt(csta))
282 ENDIF
283C
284 ! Classical nodal timestep check
285 IF (dtnod < dt2t) THEN
286 dt2t = min(dt2t,dtnod)
287 ENDIF
288 ENDIF
289c
290 DO k = 1,ndnod
291 DO i = 1,nel
292 ! Updating the non-local in-thickness velocities
293 vnlth(i,k) = vnlth(i,k) - (fnlth(i,k)/massth(i,k))*dt12
294 ENDDO
295 ENDDO
296c
297 DO k = 1,ndnod
298 DO i = 1,nel
299 ! Computing the non-local in-thickness cumulated values
300 unlth(i,k) = unlth(i,k) + vnlth(i,k)*dt1
301 ENDDO
302 ENDDO
303c
304 ! Transfert at the integration point
305 DO k = 1, nddl
306 !Computation of shape functions value
307 IF ((nddl==2).AND.(k==2)) THEN
308 nth1 = (z0(k,nddl) - zth(k,nddl))/(zth(k-1,nddl) - zth(k,nddl))
309 nth2 = (z0(k,nddl) - zth(k-1,nddl)) /(zth(k,nddl) - zth(k-1,nddl))
310 ELSE
311 nth1 = (z0(k,nddl) - zth(k+1,nddl))/(zth(k,nddl) - zth(k+1,nddl))
312 nth2 = (z0(k,nddl) - zth(k,nddl)) /(zth(k+1,nddl) - zth(k,nddl))
313 ENDIF
314 ! Loop over elements
315 DO i = 1,nel
316 !Integration points non-local variables
317 IF ((nddl==2).AND.(k==2)) THEN
318 var_reg(i,k) = nth1*unlth(i,k-1) + nth2*unlth(i,k)
319 ELSE
320 var_reg(i,k) = nth1*unlth(i,k) + nth2*unlth(i,k+1)
321 ENDIF
322 ENDDO
323 ENDDO
324 ENDIF
325c
326 !-----------------------------------------------------------------------
327 ! Computation of the elementary non-local forces
328 !-----------------------------------------------------------------------
329 ! Loop over additional degrees of freedom
330 DO k = 1, nddl
331c
332 ! Loop over elements
333# include "vectorize.inc"
334 DO i = 1, nel
335c
336 ! If the element is not broken, normal computation
337 IF (off(i) /= zero) THEN
338 ! Computation of the product LEN**2 * BtxB
339 ! Warning: the derivatives of the shape function does not take into account the volume of the element.
340 ! That is why a factor (1/VOL**2) is added in B1, B2, B3 and B4.
341 b1 = (l2 / vol(i))*wf(k,nddl)*(btb11(i)*unl(pos1(i)+k-1) + btb12(i)*unl(pos2(i)+k-1)
342 . - btb11(i)*unl(pos3(i)+k-1) - btb12(i)*unl(pos4(i)+k-1))
343c
344 b2 = (l2 / vol(i))*wf(k,nddl)*(btb12(i)*unl(pos1(i)+k-1) + btb22(i)*unl(pos2(i)+k-1)
345 . - btb12(i)*unl(pos3(i)+k-1) - btb22(i)*unl(pos4(i)+k-1))
346c
347 b3 = (l2 / vol(i))*wf(k,nddl)*(-btb11(i)*unl(pos1(i)+k-1) - btb12(i)*unl(pos2(i)+k-1)
348 . + btb11(i)*unl(pos3(i)+k-1) + btb12(i)*unl(pos4(i)+k-1))
349c
350 b4 = (l2 / vol(i))*wf(k,nddl)*(-btb12(i)*unl(pos1(i)+k-1) - btb22(i)*unl(pos2(i)+k-1)
351 . + btb12(i)*unl(pos3(i)+k-1) + btb22(i)*unl(pos4(i)+k-1))
352c
353 ! Computing the product NtN*UNL
354 ntn_unl = (unl(pos1(i)+k-1) + unl(pos2(i)+k-1) + unl(pos3(i)+k-1) + unl(pos4(i)+k-1))/ntn
355c
356 ! Computing the product XDAMP*NtN*VNL
357 ntn_vnl = (vnl(pos1(i)+k-1) + vnl(pos2(i)+k-1) + vnl(pos3(i)+k-1) + vnl(pos4(i)+k-1))/ntn
358 IF (nodadt > 0) THEN
359 ntn_vnl = min(sqrt(mass(pos1(i)+k-1)/mass0(pos1(i)+k-1)),
360 . sqrt(mass(pos2(i)+k-1)/mass0(pos2(i)+k-1)),
361 . sqrt(mass(pos3(i)+k-1)/mass0(pos3(i)+k-1)),
362 . sqrt(mass(pos4(i)+k-1)/mass0(pos4(i)+k-1)))*ntn_vnl
363 ENDIF
364c
365 ! Multiplication by the volume of the element
366 ntn_unl = ntn_unl * vol(i) * wf(k,nddl)
367 ntn_vnl = ntn_vnl * xi * vol(i) * wf(k,nddl)
368c
369 ! Introducing the internal variable to be regularized
370 ntvar = var_reg(i,k)*fourth*vol(i)*wf(k,nddl)
371c
372 ! Computing the elementary non-local forces
373 a = ntn_unl + ntn_vnl - ntvar
374 f1(i,k) = a + b1
375 f2(i,k) = a + b2
376 f3(i,k) = a + b3
377 f4(i,k) = a + b4
378c
379 ! Computing nodal equivalent stiffness
380 IF (nodadt > 0) THEN
381 sti1(i,k) = wf(k,nddl)*(abs((l2/vol(i))*btb11(i) + one/ntn*vol(i)) +
382 . abs((l2/vol(i))*btb12(i) + one/ntn*vol(i)) +
383 . abs(-(l2/vol(i))*btb11(i) + one/ntn*vol(i)) +
384 . abs(-(l2/vol(i))*btb12(i) + one/ntn*vol(i)))
385 sti2(i,k) = wf(k,nddl)*(abs((l2/vol(i))*btb12(i) + one/ntn*vol(i)) +
386 . abs((l2/vol(i))*btb22(i) + one/ntn*vol(i)) +
387 . abs(-(l2/vol(i))*btb12(i) + one/ntn*vol(i)) +
388 . abs(-(l2/vol(i))*btb22(i) + one/ntn*vol(i)))
389 sti3(i,k) = wf(k,nddl)*(abs(-(l2/vol(i))*btb11(i) + one/ntn*vol(i)) +
390 . abs(-(l2/vol(i))*btb12(i) + one/ntn*vol(i)) +
391 . abs((l2/vol(i))*btb11(i) + one/ntn*vol(i)) +
392 . abs((l2/vol(i))*btb12(i) + one/ntn*vol(i)))
393 sti4(i,k) = wf(k,nddl)*(abs(-(l2/vol(i))*btb12(i) + one/ntn*vol(i)) +
394 . abs(-(l2/vol(i))*btb22(i) + one/ntn*vol(i)) +
395 . abs((l2/vol(i))*btb12(i) + one/ntn*vol(i)) +
396 . abs((l2/vol(i))*btb22(i) + one/ntn*vol(i)))
397 ENDIF
398c
399 ! If the element is broken, the non-local wave is absorbed
400 ELSE
401 IF (nodadt > 0) THEN
402 ! Non-local absorbing forces
403 f1(i,k) = wf(k,nddl)*sqrt(mass(pos1(i)+k-1)/mass0(pos1(i)+k-1))*zeta*sspnl*
404 . half*(vnl(pos1(i)+k-1)+vnl0(pos1(i)+k-1))*sqrt(area0(i))*thk0(i)
405 f2(i,k) = wf(k,nddl)*sqrt(mass(pos2(i)+k-1)/mass0(pos2(i)+k-1))*zeta*sspnl*
406 . half*(vnl(pos2(i)+k-1)+vnl0(pos2(i)+k-1))*sqrt(area0(i))*thk0(i)
407 f3(i,k) = wf(k,nddl)*sqrt(mass(pos3(i)+k-1)/mass0(pos3(i)+k-1))*zeta*sspnl*
408 . half*(vnl(pos3(i)+k-1)+vnl0(pos3(i)+k-1))*sqrt(area0(i))*thk0(i)
409 f4(i,k) = wf(k,nddl)*sqrt(mass(pos4(i)+k-1)/mass0(pos4(i)+k-1))*zeta*sspnl*
410 . half*(vnl(pos4(i)+k-1)+vnl0(pos4(i)+k-1))*sqrt(area0(i))*thk0(i)
411 ! Computing nodal equivalent stiffness
412 sti1(i,k) = em20
413 sti2(i,k) = em20
414 sti3(i,k) = em20
415 sti4(i,k) = em20
416 ELSE
417 ! Non-local absorbing forces
418 f1(i,k) = wf(k,nddl)*zeta*sspnl*half*(vnl(pos1(i)+k-1)+vnl0(pos1(i)+k-1))*sqrt(area0(i))*thk0(i)
419 f2(i,k) = wf(k,nddl)*zeta*sspnl*half*(vnl(pos2(i)+k-1)+vnl0(pos2(i)+k-1))*sqrt(area0(i))*thk0(i)
420 f3(i,k) = wf(k,nddl)*zeta*sspnl*half*(vnl(pos3(i)+k-1)+vnl0(pos3(i)+k-1))*sqrt(area0(i))*thk0(i)
421 f4(i,k) = wf(k,nddl)*zeta*sspnl*half*(vnl(pos4(i)+k-1)+vnl0(pos4(i)+k-1))*sqrt(area0(i))*thk0(i)
422 ENDIF
423 ENDIF
424 ENDDO
425 ENDDO
426c
427 !-----------------------------------------------------------------------
428 ! Assembling of the non-local forces
429 !-----------------------------------------------------------------------
430c
431 ! If PARITH/OFF
432 IF (iparit == 0) THEN
433 ! Recovering non-local internal forces
434 fnl => nloc_dmg%FNL(1:l_nloc,itask+1) ! Non-local forces
435 IF (nodadt > 0) stifnl => nloc_dmg%STIFNL(1:l_nloc,itask+1) ! Non-local equivalent nodal stiffness
436 ! Loop over elements
437 DO i=1,nel
438 ! Loop over non-local degrees of freedom
439# include "vectorize.inc"
440 DO k = 1,nddl
441 ! Assembling non-local forces
442 fnl(pos1(i)+k-1) = fnl(pos1(i)+k-1) - f1(i,k)
443 fnl(pos2(i)+k-1) = fnl(pos2(i)+k-1) - f2(i,k)
444 fnl(pos3(i)+k-1) = fnl(pos3(i)+k-1) - f3(i,k)
445 fnl(pos4(i)+k-1) = fnl(pos4(i)+k-1) - f4(i,k)
446 IF (nodadt > 0) THEN
447 ! Spectral radius of stiffness matrix
448 maxstif = max(sti1(i,k),sti2(i,k),sti3(i,k),sti4(i,k))
449 ! Computing nodal stiffness
450 stifnl(pos1(i)+k-1) = stifnl(pos1(i)+k-1) + maxstif
451 stifnl(pos2(i)+k-1) = stifnl(pos2(i)+k-1) + maxstif
452 stifnl(pos3(i)+k-1) = stifnl(pos3(i)+k-1) + maxstif
453 stifnl(pos4(i)+k-1) = stifnl(pos4(i)+k-1) + maxstif
454 ENDIF
455 ENDDO
456 ENDDO
457c
458 ! If PARITH/ON
459 ELSE
460 ! Loop over additional d.o.fs
461 DO j = 1,nddl
462c
463 ! Loop over elements
464 DO i=1,nel
465 ii = i + nft
466c
467 ! Spectral radius of stiffness matrix
468 IF (nodadt > 0) THEN
469 maxstif = max(sti1(i,j),sti2(i,j),sti3(i,j),sti4(i,j))
470 ENDIF
471c
472 k = nloc_dmg%IADC(1,ii)
473 nloc_dmg%FSKY(k,j) = -f1(i,j)
474 IF (nodadt > 0) nloc_dmg%STSKY(k,j) = maxstif
475c
476 k = nloc_dmg%IADC(2,ii)
477 nloc_dmg%FSKY(k,j) = -f2(i,j)
478 IF (nodadt > 0) nloc_dmg%STSKY(k,j) = maxstif
479c
480 k = nloc_dmg%IADC(3,ii)
481 nloc_dmg%FSKY(k,j) = -f3(i,j)
482 IF (nodadt > 0) nloc_dmg%STSKY(k,j) = maxstif
483c
484 k = nloc_dmg%IADC(4,ii)
485 nloc_dmg%FSKY(k,j) = -f4(i,j)
486 IF (nodadt > 0) nloc_dmg%STSKY(k,j) = maxstif
487c
488 ENDDO
489c
490 ENDDO
491c
492 ENDIF
493c
494 !-----------------------------------------------------------------------
495 ! Computing non-local timestep
496 !-----------------------------------------------------------------------
497 IF (nodadt == 0) THEN
498 DO i = 1,nel
499 ! If the element is not broken, normal computation
500 IF (off(i) /= zero) THEN
501 ! Non-local critical time-step in the plane
502 dtnl = (two*(min(le(i),le_max))*sqrt(three*zeta))/
503 . sqrt(twelve*l2 + (min(le(i),le_max))**2)
504 ! Non-local critical time-step in the thickness
505 IF (nddl>1) THEN
506 IF (nddl > 2) THEN
507 dtnl_th = (two*(min(thk(i)/nddl,le_max))*sqrt(three*zeta))/
508 . sqrt(twelve*l2 + (min(thk(i)/nddl,le_max))**2)
509 ELSE
510 dtnl_th = (two*(min(thk(i),le_max))*sqrt(three*zeta))/
511 . sqrt(twelve*l2 + (min(thk(i),le_max))**2)
512 ENDIF
513 ELSE
514 dtnl_th = ep20
515 ENDIF
516 ! Retaining the minimal value
517 dt2t = min(dt2t,dtfac1(1)*cdamp*dtnl_th,dtfac1(1)*cdamp*dtnl)
518 ENDIF
519 ENDDO
520 ENDIF
521c
522 ! Deallocation of tables
523 IF (ALLOCATED(f1)) DEALLOCATE(f1)
524 IF (ALLOCATED(f2)) DEALLOCATE(f2)
525 IF (ALLOCATED(f3)) DEALLOCATE(f3)
526 IF (ALLOCATED(f4)) DEALLOCATE(f4)
527 IF (ALLOCATED(sti1)) DEALLOCATE(sti1)
528 IF (ALLOCATED(sti2)) DEALLOCATE(sti2)
529 IF (ALLOCATED(sti3)) DEALLOCATE(sti3)
530 IF (ALLOCATED(sti4)) DEALLOCATE(sti4)
531 IF (ALLOCATED(stifnlth)) DEALLOCATE(stifnlth)
532 IF (ALLOCATED(btb11)) DEALLOCATE(btb11)
533 IF (ALLOCATED(btb12)) DEALLOCATE(btb12)
534 IF (ALLOCATED(btb22)) DEALLOCATE(btb22)
535 IF (ALLOCATED(pos1)) DEALLOCATE(pos1)
536 IF (ALLOCATED(pos2)) DEALLOCATE(pos2)
537 IF (ALLOCATED(pos3)) DEALLOCATE(pos3)
538 IF (ALLOCATED(pos4)) DEALLOCATE(pos4)
539 IF (ALLOCATED(vol)) DEALLOCATE(vol)
540c
541 END
cfint_reg
subroutine cfint_reg(nloc_dmg, var_reg, thk, nel, off, area, nc1, nc2, nc3, nc4, px1, py1, px2, py2, bufnl, imat, nddl, itask, dt2t, le, thk0, area0, nft)
Definition cfint_reg.F:39
area
subroutine area(d1, x, x2, y, y2, eint, stif0)
Definition hm_read_prop32.F:567
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21
nlocal_reg_mod
Definition nlocal_reg_mod.F:117
nlocal_reg_mod::nlocal_str_
Definition nlocal_reg_mod.F:124