OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
cdkfint_reg.F File Reference
#include "implicit_f.inc"
#include "parit_c.inc"
#include "com20_c.inc"
#include "com08_c.inc"
#include "scr02_c.inc"
#include "scr18_c.inc"
#include "vectorize.inc"

Go to the source code of this file.

Functions/Subroutines

subroutine cdkfint_reg (nloc_dmg, var_reg, thk, nel, off, area, nc1, nc2, nc3, px2, py2, px3, py3, ksi, eta, bufnl, imat, nddl, itask, ng, dt2t, thk0, area0, nft)

Function/Subroutine Documentation

◆ cdkfint_reg()

subroutine cdkfint_reg ( type(nlocal_str_), target nloc_dmg,
intent(inout) var_reg,
intent(in) thk,
integer nel,
intent(in) off,
intent(in) area,
integer, dimension(nel) nc1,
integer, dimension(nel) nc2,
integer, dimension(nel) nc3,
intent(in) px2,
intent(in) py2,
intent(in) px3,
intent(in) py3,
ksi,
eta,
type(buf_nloc_), target bufnl,
integer imat,
integer nddl,
integer itask,
integer ng,
intent(inout) dt2t,
intent(in) thk0,
intent(in) area0,
integer, intent(in) nft )

Definition at line 31 of file cdkfint_reg.F.

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