OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
sfint_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!|| sfint_reg ../engine/source/elements/solid/solide/sfint_reg.F
25!||--- called by ------------------------------------------------------
26!|| sforc3 ../engine/source/elements/solid/solide/sforc3.F
27!|| szforc3 ../engine/source/elements/solid/solidez/szforc3.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 sfint_reg(
33 1 NLOC_DMG,VAR_REG, NEL, OFF,
34 2 VOL, NC1, NC2, NC3,
35 3 NC4, NC5, NC6, NC7,
36 4 NC8, PX1, PX2, PX3,
37 5 PX4, PY1, PY2, PY3,
38 6 PY4, PZ1, PZ2, PZ3,
39 7 PZ4, IMAT, ITASK, DT2T,
40 8 VOL0, NFT)
41C-----------------------------------------------
42C M o d u l e s
43C-----------------------------------------------
45 USE elbufdef_mod
46C-----------------------------------------------
47C I m p l i c i t T y p e s
48C-----------------------------------------------
49#include "implicit_f.inc"
50C-----------------------------------------------
51C G l o b a l P a r a m e t e r s
52C-----------------------------------------------
53#include "parit_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 TYPE(nlocal_str_), INTENT(INOUT), TARGET :: NLOC_DMG !< Non-local data structure
60 my_real, DIMENSION(NEL), INTENT(IN) :: VAR_REG !< Variable to regularize
61 INTEGER, INTENT(IN) :: NEL !< Number of elements in the group
62 my_real, DIMENSION(NEL), INTENT(IN) :: OFF !< Flag for element deletion
63 my_real, DIMENSION(NEL), INTENT(IN) :: VOL !< Current volume of the element
64 INTEGER, DIMENSION(NEL) :: NC1 !< Node 1 of the elements
65 INTEGER, DIMENSION(NEL) :: NC2 !< Node 2 of the elements
66 INTEGER, DIMENSION(NEL) :: NC3 !< Node 3 of the elements
67 INTEGER, DIMENSION(NEL) :: NC4 !< Node 4 of the elements
68 INTEGER, DIMENSION(NEL) :: NC5 !< Node 5 of the elements
69 INTEGER, DIMENSION(NEL) :: NC6 !< Node 6 of the elements
70 INTEGER, DIMENSION(NEL) :: NC7 !< Node 7 of the elements
71 INTEGER, DIMENSION(NEL) :: NC8 !< Node 8 of the elements
72 my_real, DIMENSION(NEL), INTENT(IN) :: px1 !< Derivative of the first node shape function along x
73 my_real, DIMENSION(NEL), INTENT(IN) :: px2 !< Derivative of the second node shape function along x
74 my_real, DIMENSION(NEL), INTENT(IN) :: px3 !< Derivative of the third node shape function along x
75 my_real, DIMENSION(NEL), INTENT(IN) :: px4 !< Derivative of the fourth node shape function along x
76 my_real, DIMENSION(NEL), INTENT(IN) :: py1 !< Derivative of the first node shape function along y
77 my_real, DIMENSION(NEL), INTENT(IN) :: py2 !< Derivative of the second node shape function along y
78 my_real, DIMENSION(NEL), INTENT(IN) :: py3 !< Derivative of the third node shape function along y
79 my_real, DIMENSION(NEL), INTENT(IN) :: py4 !< Derivative of the fourth node shape function along y
80 my_real, DIMENSION(NEL), INTENT(IN) :: pz1 !< Derivative of the first node shape function along z
81 my_real, DIMENSION(NEL), INTENT(IN) :: pz2 !< Derivative of the second node shape function along z
82 my_real, DIMENSION(NEL), INTENT(IN) :: pz3 !< Derivative of the third node shape function along z
83 my_real, DIMENSION(NEL), INTENT(IN) :: pz4 !< Derivative of the fourth node shape function along z
84 INTEGER, INTENT(IN) :: IMAT !< Material internal Id
85 INTEGER, INTENT(IN) :: ITASK !< Number of the thread
86 my_real, INTENT(INOUT) :: dt2t !< Element critical timestep
87 my_real, DIMENSION(NEL), INTENT(IN) :: vol0 !< Initial volume of the element
88 INTEGER, INTENT(IN) :: NFT !< Address of the first element of the group
89C-----------------------------------------------
90C L o c a l V a r i a b l e s
91C-----------------------------------------------
92 INTEGER I,II,K,N1,N2,N3,N4,N5,N6,N7,N8,L_NLOC
93 INTEGER, DIMENSION(:), ALLOCATABLE ::
94 . pos1,pos2,pos3,pos4,pos5,pos6,pos7,pos8
95 my_real
96 . l2,ntn_unl,ntn_vnl,xi,ntvar,a,
97 . b1,b2,b3,b4,b5,b6,b7,b8,zeta,sspnl,dtnl,le_max,
98 . maxstif,ntn
99 my_real, DIMENSION(:) ,ALLOCATABLE ::
100 . btb11,btb12,btb13,btb14,btb22,btb23,btb24,
101 . btb33,btb34,btb44,sti1,sti2,sti3,sti4,sti5,
102 . sti6,sti7,sti8,f1,f2,f3,f4,f5,f6,f7,f8,lc
103 my_real, POINTER, DIMENSION(:) ::
104 . vnl,fnl,unl,stifnl,mass,mass0,vnl0
105 ! Coefficient for non-local stability to take into account damping
106 my_real, target :: nothing(1)
107c-----------------------------------------------------------------------
108c VAR_REG : variable a regulariser (local cumulated plastic strain)
109c NTVAR = NT * VAR_REG
110C=======================================================================
111c
112 nothing = zero
113 vnl => nothing
114 fnl => nothing
115 unl => nothing
116 stifnl => nothing
117 mass => nothing
118 mass0 => nothing
119 vnl0 => nothing
120
121 ! Recover non-local parameters
122 l2 = nloc_dmg%LEN(imat)**2 !< Squared internal length
123 xi = nloc_dmg%DAMP(imat) !< Damping parameter
124 l_nloc = nloc_dmg%L_NLOC !< Size of non-local tables
125 zeta = nloc_dmg%DENS(imat) !< Density parameter
126 sspnl = nloc_dmg%SSPNL(imat) !< Non-local "Sound speed"
127 le_max = nloc_dmg%LE_MAX(imat) !< Maximal length of convergence
128c
129 ! Allocation of variable table
130 ALLOCATE(btb11(nel),btb12(nel),btb13(nel),btb14(nel),btb22(nel),
131 . btb23(nel),btb24(nel),btb33(nel),btb34(nel),btb44(nel),pos1(nel),
132 . pos2(nel),pos3(nel),pos4(nel),pos5(nel),pos6(nel),pos7(nel),pos8(nel),
133 . f1(nel),f2(nel),f3(nel),f4(nel),f5(nel),f6(nel),f7(nel),f8(nel),lc(nel))
134c
135 ! Factor for non-local forces
136 ntn = eight*eight
137c
138 ! Local variables initialization
139 lc(1:nel) = zero
140 ! For nodal timestep
141 IF (nodadt > 0) THEN
142 ! Non-local nodal stifness
143 ALLOCATE(sti1(nel),sti2(nel),sti3(nel),sti4(nel),sti5(nel),sti6(nel),
144 . sti7(nel),sti8(nel))
145 ! Non-local mass
146 mass => nloc_dmg%MASS(1:l_nloc)
147 ! Initial non-local mass
148 mass0 => nloc_dmg%MASS0(1:l_nloc)
149 ENDIF
150 ! Current timestep non-local velocities
151 vnl => nloc_dmg%VNL(1:l_nloc)
152 ! Previous timestep non-local velocities
153 vnl0 => nloc_dmg%VNL_OLD(1:l_nloc)
154 ! Non-local displacements
155 unl => nloc_dmg%UNL(1:l_nloc)
156c
157 !--------------------------------------------------------------------------------
158 ! Computation of the position of the non-local d.o.fs and the BtB matrix product
159 !--------------------------------------------------------------------------------
160 ! Loop over elements
161 DO i=1,nel
162c
163 ! Recovering the nodes of the brick element
164 n1 = nloc_dmg%IDXI(nc1(i))
165 n2 = nloc_dmg%IDXI(nc2(i))
166 n3 = nloc_dmg%IDXI(nc3(i))
167 n4 = nloc_dmg%IDXI(nc4(i))
168 n5 = nloc_dmg%IDXI(nc5(i))
169 n6 = nloc_dmg%IDXI(nc6(i))
170 n7 = nloc_dmg%IDXI(nc7(i))
171 n8 = nloc_dmg%IDXI(nc8(i))
172c
173 ! Recovering the positions of the first d.o.fs of each nodes
174 pos1(i) = nloc_dmg%POSI(n1)
175 pos2(i) = nloc_dmg%POSI(n2)
176 pos3(i) = nloc_dmg%POSI(n3)
177 pos4(i) = nloc_dmg%POSI(n4)
178 pos5(i) = nloc_dmg%POSI(n5)
179 pos6(i) = nloc_dmg%POSI(n6)
180 pos7(i) = nloc_dmg%POSI(n7)
181 pos8(i) = nloc_dmg%POSI(n8)
182c
183 ! Computation of the product BtxB
184 btb11(i) = px1(i)**2 + py1(i)**2 + pz1(i)**2
185 btb12(i) = px1(i)*px2(i) + py1(i)*py2(i) + pz1(i)*pz2(i)
186 btb13(i) = px1(i)*px3(i) + py1(i)*py3(i) + pz1(i)*pz3(i)
187 btb14(i) = px1(i)*px4(i) + py1(i)*py4(i) + pz1(i)*pz4(i)
188 btb22(i) = px2(i)**2 + py2(i)**2 + pz2(i)**2
189 btb23(i) = px2(i)*px3(i) + py2(i)*py3(i) + pz2(i)*pz3(i)
190 btb24(i) = px2(i)*px4(i) + py2(i)*py4(i) + pz2(i)*pz4(i)
191 btb33(i) = px3(i)**2 + py3(i)**2 + pz3(i)**2
192 btb34(i) = px3(i)*px4(i) + py3(i)*py4(i) + pz3(i)*pz4(i)
193 btb44(i) = px4(i)**2 + py4(i)**2 + pz4(i)**2
194c
195 ENDDO
196c
197 !-----------------------------------------------------------------------
198 ! Computation of non-local forces
199 !-----------------------------------------------------------------------
200 ! Loop over elements
201 DO i = 1, nel
202c
203 ! If the element is not broken, normal computation
204 IF (off(i) /= zero) THEN
205c
206 ! Computing the product NtN*UNL
207 ntn_unl = (unl(pos1(i)) + unl(pos2(i)) + unl(pos3(i)) + unl(pos4(i))
208 . + unl(pos5(i)) + unl(pos6(i)) + unl(pos7(i)) + unl(pos8(i))) / ntn
209c
210 ! Computing the product NtN*VNL
211 ntn_vnl = (vnl(pos1(i)) + vnl(pos2(i)) + vnl(pos3(i)) + vnl(pos4(i))
212 . + vnl(pos5(i)) + vnl(pos6(i)) + vnl(pos7(i)) + vnl(pos8(i))) / ntn
213 IF (nodadt > 0) THEN
214 ntn_vnl = min(sqrt(mass(pos1(i))/mass0(pos1(i))),
215 . sqrt(mass(pos2(i))/mass0(pos2(i))),
216 . sqrt(mass(pos3(i))/mass0(pos3(i))),
217 . sqrt(mass(pos4(i))/mass0(pos4(i))),
218 . sqrt(mass(pos5(i))/mass0(pos5(i))),
219 . sqrt(mass(pos6(i))/mass0(pos6(i))),
220 . sqrt(mass(pos7(i))/mass0(pos7(i))),
221 . sqrt(mass(pos8(i))/mass0(pos8(i))))*ntn_vnl
222 ENDIF
223c
224 ! Computation of the product LEN**2 * BtxB
225 b1 = l2 * vol(i) * ( btb11(i)*unl(pos1(i)) + btb12(i)*unl(pos2(i))
226 . + btb13(i)*unl(pos3(i)) + btb14(i)*unl(pos4(i)) - btb13(i)*unl(pos5(i))
227 . - btb14(i)*unl(pos6(i)) - btb11(i)*unl(pos7(i)) - btb12(i)*unl(pos8(i)))
228c
229 b2 = l2 * vol(i) * ( btb12(i)*unl(pos1(i)) + btb22(i)*unl(pos2(i))
230 . + btb23(i)*unl(pos3(i)) + btb24(i)*unl(pos4(i)) - btb23(i)*unl(pos5(i))
231 . - btb24(i)*unl(pos6(i)) - btb12(i)*unl(pos7(i)) - btb22(i)*unl(pos8(i)))
232c
233 b3 = l2 * vol(i) * ( btb13(i)*unl(pos1(i)) + btb23(i)*unl(pos2(i))
234 . + btb33(i)*unl(pos3(i)) + btb34(i)*unl(pos4(i)) - btb33(i)*unl(pos5(i))
235 . - btb34(i)*unl(pos6(i)) - btb13(i)*unl(pos7(i)) - btb23(i)*unl(pos8(i)))
236c
237 b4 = l2 * vol(i) * ( btb14(i)*unl(pos1(i)) + btb24(i)*unl(pos2(i))
238 . + btb34(i)*unl(pos3(i)) + btb44(i)*unl(pos4(i)) - btb34(i)*unl(pos5(i))
239 . - btb44(i)*unl(pos6(i)) - btb14(i)*unl(pos7(i)) - btb24(i)*unl(pos8(i)))
240c
241 b5 = l2 * vol(i) * ( -btb13(i)*unl(pos1(i)) - btb23(i)*unl(pos2(i))
242 . - btb33(i)*unl(pos3(i)) - btb34(i)*unl(pos4(i)) + btb33(i)*unl(pos5(i))
243 . + btb34(i)*unl(pos6(i)) + btb13(i)*unl(pos7(i)) + btb23(i)*unl(pos8(i)))
244c
245 b6 = l2 * vol(i) * ( -btb14(i)*unl(pos1(i)) - btb24(i)*unl(pos2(i))
246 . - btb34(i)*unl(pos3(i)) - btb44(i)*unl(pos4(i)) + btb34(i)*unl(pos5(i))
247 . + btb44(i)*unl(pos6(i)) + btb14(i)*unl(pos7(i)) + btb24(i)*unl(pos8(i)))
248c
249 b7 = l2 * vol(i) * ( -btb11(i)*unl(pos1(i)) - btb12(i)*unl(pos2(i))
250 . - btb13(i)*unl(pos3(i)) - btb14(i)*unl(pos4(i)) + btb13(i)*unl(pos5(i))
251 . + btb14(i)*unl(pos6(i)) + btb11(i)*unl(pos7(i)) + btb12(i)*unl(pos8(i)))
252c
253 b8 = l2 * vol(i) * ( -btb12(i)*unl(pos1(i)) - btb22(i)*unl(pos2(i))
254 . - btb23(i)*unl(pos3(i)) - btb24(i)*unl(pos4(i)) + btb23(i)*unl(pos5(i))
255 . + btb24(i)*unl(pos6(i)) + btb12(i)*unl(pos7(i)) + btb22(i)*unl(pos8(i)))
256c
257 ! Multiplication by the volume of the element (and damping parameter XI)
258 ntn_unl = ntn_unl * vol(i)
259 ntn_vnl = ntn_vnl * xi * vol(i)
260c
261 ! Introducing the internal variable to be regularized
262 ntvar = var_reg(i)*one_over_8* vol(i)
263c
264 ! Computing the elementary non-local forces
265 a = ntn_unl + ntn_vnl - ntvar
266 f1(i) = a + b1
267 f2(i) = a + b2
268 f3(i) = a + b3
269 f4(i) = a + b4
270 f5(i) = a + b5
271 f6(i) = a + b6
272 f7(i) = a + b7
273 f8(i) = a + b8
274c
275 ! Computing nodal equivalent stiffness
276 IF (nodadt > 0) THEN
277 sti1(i) = (abs(l2*btb11(i) + one/ntn) + abs(l2*btb12(i) + one/ntn) + abs(l2*btb13(i) + one/ntn) +
278 . abs(l2*btb14(i) + one/ntn) + abs(-l2*btb13(i) + one/ntn) + abs(-l2*btb14(i) + one/ntn) +
279 . abs(-l2*btb11(i) + one/ntn) + abs(-l2*btb12(i) + one/ntn))*vol(i)
280 sti2(i) = (abs(l2*btb12(i) + one/ntn) + abs(l2*btb22(i) + one/ntn) + abs(l2*btb23(i) + one/ntn) +
281 . abs(l2*btb24(i) + one/ntn) + abs(-l2*btb23(i) + one/ntn) + abs(-l2*btb24(i) + one/ntn) +
282 . abs(-l2*btb12(i) + one/ntn) + abs(-l2*btb22(i) + one/ntn))*vol(i)
283 sti3(i) = (abs(l2*btb13(i) + one/ntn) + abs(l2*btb23(i) + one/ntn) + abs(l2*btb33(i) + one/ntn) +
284 . abs(l2*btb34(i) + one/ntn) + abs(-l2*btb33(i) + one/ntn) + abs(-l2*btb34(i) + one/ntn) +
285 . abs(-l2*btb13(i) + one/ntn) + abs(-l2*btb23(i) + one/ntn))*vol(i)
286 sti4(i) = (abs(l2*btb14(i) + one/ntn) + abs(l2*btb24(i) + one/ntn) + abs(l2*btb34(i) + one/ntn) +
287 . abs(l2*btb44(i) + one/ntn) + abs(-l2*btb34(i) + one/ntn) + abs(-l2*btb44(i) + one/ntn) +
288 . abs(-l2*btb14(i) + one/ntn) + abs(-l2*btb24(i) + one/ntn))*vol(i)
289 sti5(i) = (abs(-l2*btb13(i) + one/ntn) + abs(-l2*btb23(i) + one/ntn) + abs(-l2*btb33(i) + one/ntn) +
290 . abs(-l2*btb34(i) + one/ntn) + abs(l2*btb33(i) + one/ntn) + abs(l2*btb34(i) + one/ntn) +
291 . abs(l2*btb13(i) + one/ntn) + abs(l2*btb23(i) + one/ntn))*vol(i)
292 sti6(i) = (abs(-l2*btb14(i) + one/ntn) + abs(-l2*btb24(i) + one/ntn) + abs(-l2*btb34(i) + one/ntn) +
293 . abs(-l2*btb44(i) + one/ntn) + abs(l2*btb34(i) + one/ntn) + abs(l2*btb44(i) + one/ntn) +
294 . abs(l2*btb14(i) + one/ntn) + abs(l2*btb24(i) + one/ntn))*vol(i)
295 sti7(i) = (abs(-l2*btb11(i) + one/ntn) + abs(-l2*btb12(i) + one/ntn) + abs(-l2*btb13(i) + one/ntn) +
296 . abs(-l2*btb14(i) + one/ntn) + abs(l2*btb13(i) + one/ntn) + abs(l2*btb14(i) + one/ntn) +
297 . abs(l2*btb11(i) + one/ntn) + abs(l2*btb12(i) + one/ntn))*vol(i)
298 sti8(i) = (abs(-l2*btb12(i) + one/ntn) + abs(-l2*btb22(i) + one/ntn) + abs(-l2*btb23(i) + one/ntn) +
299 . abs(-l2*btb24(i) + one/ntn) + abs(l2*btb23(i) + one/ntn) + abs(l2*btb24(i) + one/ntn) +
300 . abs(l2*btb12(i) + one/ntn) + abs(l2*btb22(i) + one/ntn))*vol(i)
301 ENDIF
302c
303 ! If the element is broken, the non-local wave is absorbed
304 ELSE
305c
306 ! Initial element characteristic length
307 lc(i) = vol0(i)**third
308c
309 IF (nodadt > 0) THEN
310
311 ! Non-local absorbing forces
312 f1(i) = sqrt(mass(pos1(i))/mass0(pos1(i)))*zeta*sspnl*half*
313 . (vnl(pos1(i))+vnl0(pos1(i)))*(three/four)*(lc(i)**2)
314 f2(i) = sqrt(mass(pos2(i))/mass0(pos2(i)))*zeta*sspnl*half*
315 . (vnl(pos2(i))+vnl0(pos2(i)))*(three/four)*(lc(i)**2)
316 f3(i) = sqrt(mass(pos3(i))/mass0(pos3(i)))*zeta*sspnl*half*
317 . (vnl(pos3(i))+vnl0(pos3(i)))*(three/four)*(lc(i)**2)
318 f4(i) = sqrt(mass(pos4(i))/mass0(pos4(i)))*zeta*sspnl*half*
319 . (vnl(pos4(i))+vnl0(pos4(i)))*(three/four)*(lc(i)**2)
320 f5(i) = sqrt(mass(pos5(i))/mass0(pos5(i)))*zeta*sspnl*half*
321 . (vnl(pos5(i))+vnl0(pos5(i)))*(three/four)*(lc(i)**2)
322 f6(i) = sqrt(mass(pos6(i))/mass0(pos6(i)))*zeta*sspnl*half*
323 . (vnl(pos6(i))+vnl0(pos6(i)))*(three/four)*(lc(i)**2)
324 f7(i) = sqrt(mass(pos7(i))/mass0(pos7(i)))*zeta*sspnl*half*
325 . (vnl(pos7(i))+vnl0(pos7(i)))*(three/four)*(lc(i)**2)
326 f8(i) = sqrt(mass(pos8(i))/mass0(pos8(i)))*zeta*sspnl*half*
327 . (vnl(pos8(i))+vnl0(pos8(i)))*(three/four)*(lc(i)**2)
328 ! Computing nodal equivalent stiffness
329 sti1(i) = em20
330 sti2(i) = em20
331 sti3(i) = em20
332 sti4(i) = em20
333 sti5(i) = em20
334 sti6(i) = em20
335 sti7(i) = em20
336 sti8(i) = em20
337 ELSE
338 ! Non-local absorbing forces
339 f1(i) = zeta*sspnl*half*(vnl(pos1(i))+vnl0(pos1(i)))*(three/four)*(lc(i)**2)
340 f2(i) = zeta*sspnl*half*(vnl(pos2(i))+vnl0(pos2(i)))*(three/four)*(lc(i)**2)
341 f3(i) = zeta*sspnl*half*(vnl(pos3(i))+vnl0(pos3(i)))*(three/four)*(lc(i)**2)
342 f4(i) = zeta*sspnl*half*(vnl(pos4(i))+vnl0(pos4(i)))*(three/four)*(lc(i)**2)
343 f5(i) = zeta*sspnl*half*(vnl(pos5(i))+vnl0(pos5(i)))*(three/four)*(lc(i)**2)
344 f6(i) = zeta*sspnl*half*(vnl(pos6(i))+vnl0(pos6(i)))*(three/four)*(lc(i)**2)
345 f7(i) = zeta*sspnl*half*(vnl(pos7(i))+vnl0(pos7(i)))*(three/four)*(lc(i)**2)
346 f8(i) = zeta*sspnl*half*(vnl(pos8(i))+vnl0(pos8(i)))*(three/four)*(lc(i)**2)
347 ENDIF
348 ENDIF
349 ENDDO
350c-----------------------------------------------------------------------
351c Assemblage
352c-----------------------------------------------------------------------
353 ! If PARITH/OFF
354 IF (iparit == 0) THEN
355 fnl => nloc_dmg%FNL(1:l_nloc,itask+1)
356 IF (nodadt > 0) stifnl => nloc_dmg%STIFNL(1:l_nloc,itask+1) ! Non-local equivalent nodal stiffness
357 ! Loop over elements
358 DO i=1,nel
359 ! Assembling the forces in the classic way
360 fnl(pos1(i)) = fnl(pos1(i)) - f1(i)
361 fnl(pos2(i)) = fnl(pos2(i)) - f2(i)
362 fnl(pos3(i)) = fnl(pos3(i)) - f3(i)
363 fnl(pos4(i)) = fnl(pos4(i)) - f4(i)
364 fnl(pos5(i)) = fnl(pos5(i)) - f5(i)
365 fnl(pos6(i)) = fnl(pos6(i)) - f6(i)
366 fnl(pos7(i)) = fnl(pos7(i)) - f7(i)
367 fnl(pos8(i)) = fnl(pos8(i)) - f8(i)
368 IF (nodadt > 0) THEN
369 ! Spectral radius of stiffness matrix
370 maxstif = max(sti1(i),sti2(i),sti3(i),sti4(i),sti5(i),sti6(i),sti7(i),sti8(i))
371 ! Computing nodal stiffness
372 stifnl(pos1(i)) = stifnl(pos1(i)) + maxstif
373 stifnl(pos2(i)) = stifnl(pos2(i)) + maxstif
374 stifnl(pos3(i)) = stifnl(pos3(i)) + maxstif
375 stifnl(pos4(i)) = stifnl(pos4(i)) + maxstif
376 stifnl(pos5(i)) = stifnl(pos5(i)) + maxstif
377 stifnl(pos6(i)) = stifnl(pos6(i)) + maxstif
378 stifnl(pos7(i)) = stifnl(pos7(i)) + maxstif
379 stifnl(pos8(i)) = stifnl(pos8(i)) + maxstif
380 ENDIF
381 ENDDO
382c
383 ! If PARITH/ON
384 ELSE
385c
386 ! Loop over elements
387 DO i=1,nel
388 ii = i + nft
389c
390 ! Spectral radius of stiffness matrix
391 IF (nodadt > 0) THEN
392 maxstif = max(sti1(i),sti2(i),sti3(i),sti4(i),sti5(i),sti6(i),sti7(i),sti8(i))
393 ENDIF
394c
395 k = nloc_dmg%IADS(1,ii)
396 nloc_dmg%FSKY(k,1) = -f1(i)
397 IF (nodadt > 0) nloc_dmg%STSKY(k,1) = maxstif
398c
399 k = nloc_dmg%IADS(2,ii)
400 nloc_dmg%FSKY(k,1) = -f2(i)
401 IF (nodadt > 0) nloc_dmg%STSKY(k,1) = maxstif
402c
403 k = nloc_dmg%IADS(3,ii)
404 nloc_dmg%FSKY(k,1) = -f3(i)
405 IF (nodadt > 0) nloc_dmg%STSKY(k,1) = maxstif
406c
407 k = nloc_dmg%IADS(4,ii)
408 nloc_dmg%FSKY(k,1) = -f4(i)
409 IF (nodadt > 0) nloc_dmg%STSKY(k,1) = maxstif
410c
411 k = nloc_dmg%IADS(5,ii)
412 nloc_dmg%FSKY(k,1) = -f5(i)
413 IF (nodadt > 0) nloc_dmg%STSKY(k,1) = maxstif
414c
415 k = nloc_dmg%IADS(6,ii)
416 nloc_dmg%FSKY(k,1) = -f6(i)
417 IF (nodadt > 0) nloc_dmg%STSKY(k,1) = maxstif
418c
419 k = nloc_dmg%IADS(7,ii)
420 nloc_dmg%FSKY(k,1) = -f7(i)
421 IF (nodadt > 0) nloc_dmg%STSKY(k,1) = maxstif
422c
423 k = nloc_dmg%IADS(8,ii)
424 nloc_dmg%FSKY(k,1) = -f8(i)
425 IF (nodadt > 0) nloc_dmg%STSKY(k,1) = maxstif
426c
427 ENDDO
428 ENDIF
429c
430 !-----------------------------------------------------------------------
431 ! Computing non-local timestep
432 !-----------------------------------------------------------------------
433 IF (nodadt == 0) THEN
434 DO i = 1,nel
435 ! If the element is not broken, normal computation
436 IF (off(i)/=zero) THEN
437 ! Non-local critical time-step in the plane
438 dtnl = (two*(min(vol(i)**third,le_max))*sqrt(three*zeta))/
439 . sqrt(twelve*l2 + (min(vol(i)**third,le_max))**2)
440 ! Retaining the minimal value
441 dt2t = min(dt2t,dtfac1(1)*cdamp*dtnl)
442 ENDIF
443 ENDDO
444 ENDIF
445c
446 ! Deallocation of tables
447 IF (ALLOCATED(btb11)) DEALLOCATE(btb11)
448 IF (ALLOCATED(btb12)) DEALLOCATE(btb12)
449 IF (ALLOCATED(btb13)) DEALLOCATE(btb13)
450 IF (ALLOCATED(btb14)) DEALLOCATE(btb14)
451 IF (ALLOCATED(btb22)) DEALLOCATE(btb22)
452 IF (ALLOCATED(btb23)) DEALLOCATE(btb23)
453 IF (ALLOCATED(btb24)) DEALLOCATE(btb24)
454 IF (ALLOCATED(btb33)) DEALLOCATE(btb33)
455 IF (ALLOCATED(btb34)) DEALLOCATE(btb34)
456 IF (ALLOCATED(btb44)) DEALLOCATE(btb44)
457 IF (ALLOCATED(pos1)) DEALLOCATE(pos1)
458 IF (ALLOCATED(pos2)) DEALLOCATE(pos2)
459 IF (ALLOCATED(pos3)) DEALLOCATE(pos3)
460 IF (ALLOCATED(pos4)) DEALLOCATE(pos4)
461 IF (ALLOCATED(pos5)) DEALLOCATE(pos5)
462 IF (ALLOCATED(pos6)) DEALLOCATE(pos6)
463 IF (ALLOCATED(pos7)) DEALLOCATE(pos7)
464 IF (ALLOCATED(pos8)) DEALLOCATE(pos8)
465 IF (ALLOCATED(f1)) DEALLOCATE(f1)
466 IF (ALLOCATED(f2)) DEALLOCATE(f2)
467 IF (ALLOCATED(f3)) DEALLOCATE(f3)
468 IF (ALLOCATED(f4)) DEALLOCATE(f4)
469 IF (ALLOCATED(f5)) DEALLOCATE(f5)
470 IF (ALLOCATED(f6)) DEALLOCATE(f6)
471 IF (ALLOCATED(f7)) DEALLOCATE(f7)
472 IF (ALLOCATED(f8)) DEALLOCATE(f8)
473 IF (ALLOCATED(sti1)) DEALLOCATE(sti1)
474 IF (ALLOCATED(sti2)) DEALLOCATE(sti2)
475 IF (ALLOCATED(sti3)) DEALLOCATE(sti3)
476 IF (ALLOCATED(sti4)) DEALLOCATE(sti4)
477 IF (ALLOCATED(sti5)) DEALLOCATE(sti5)
478 IF (ALLOCATED(sti6)) DEALLOCATE(sti6)
479 IF (ALLOCATED(sti7)) DEALLOCATE(sti7)
480 IF (ALLOCATED(sti8)) DEALLOCATE(sti8)
481 IF (ALLOCATED(lc)) DEALLOCATE(lc)
482c-----------
483 END
#define min(a, b)
Definition macros.h:20
#define max(a, b)
Definition macros.h:21
subroutine sfint_reg(nloc_dmg, var_reg, nel, off, vol, nc1, nc2, nc3, nc4, nc5, nc6, nc7, nc8, px1, px2, px3, px4, py1, py2, py3, py4, pz1, pz2, pz3, pz4, imat, itask, dt2t, vol0, nft)
Definition sfint_reg.F:41