OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
chvis3.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!|| chvis3 ../engine/source/elements/shell/coque/chvis3.f
25!||--- called by ------------------------------------------------------
26!|| cforc3 ../engine/source/elements/shell/coque/cforc3.F
27!|| cforc3_crk ../engine/source/elements/xfem/cforc3_crk.F
28!||--- uses -----------------------------------------------------
29!|| element_mod ../common_source/modules/elements/element_mod.f90
30!||====================================================================
31 SUBROUTINE chvis3(
32 2 JFT ,JLT ,THK ,HOUR,OFF ,PX1 ,PX2 ,PY1 ,PY2 ,
33 3 IXC ,DT1C,SSP ,RHO ,STI ,VX1 ,VX2 ,VX3 ,VX4 ,VY1 ,
34 4 VY2 ,VY3 ,VY4 ,VZ1 ,VZ2 ,VZ3 ,VZ4 ,AREA,THK0,VHX ,
35 5 VHY ,SHF ,Z2 ,EANI,STIR,VISCMX,G ,A11 ,
36 6 H1 ,H2 ,H3 ,YM ,NU ,THK02 ,ALPE ,H11 ,
37 7 H12 ,H13 ,H21 ,H22 ,H23 ,H31 ,H32 ,H33 ,
38 8 B11 ,B12 ,B13 ,B14 ,B21 ,B22 ,B23 ,B24 ,
39 9 RX1 ,RX2 ,RX3 ,RX4 ,RY1 ,RY2 ,RY3 ,RY4 ,
40 A IPARTC,PARTSAV,
41 B IHBE ,NFT ,ISMSTR,KFTS ,
42 C SRH1 ,SRH2 ,SRH3 ,IGTYP ,
43 D IGMAT ,A11R,NEL)
44 use element_mod , only : nixc
45C-----------------------------------------------
46C I m p l i c i t T y p e s
47C-----------------------------------------------
48#include "implicit_f.inc"
49#include "comlock.inc"
50C-----------------------------------------------
51C G l o b a l P a r a m e t e r s
52C-----------------------------------------------
53#include "mvsiz_p.inc"
54C-----------------------------------------------
55C C o m m o n B l o c k s
56C-----------------------------------------------
57#include "com04_c.inc"
58#include "com06_c.inc"
59#include "param_c.inc"
60#include "remesh_c.inc"
61#include "scr02_c.inc"
62#include "scr06_c.inc"
63#include "scr14_c.inc"
64#include "scr16_c.inc"
65#include "scr17_c.inc"
66#include "sms_c.inc"
67C-----------------------------------------------
68C D u m m y A r g u m e n t s
69C-----------------------------------------------
70 INTEGER IXC(NIXC,*),IPARTC(*), JFT, JLT,
71 . IHBE ,NFT ,ISMSTR,KFTS,IGTYP,
72 . IGMAT,NEL
73C REAL
74 my_real
75 . THK(*), HOUR(NEL,5), OFF(*),PARTSAV(NPSAV,*),
76 . PX1(*), PX2(*), PY1(*), PY2(*),DT1C(*),EANI(*),
77 . SSP(MVSIZ), RHO(MVSIZ),STI(MVSIZ),STIR(MVSIZ),
78 . VX1(MVSIZ), VX2(MVSIZ), VX3(MVSIZ), VX4(MVSIZ),
79 . VY1(MVSIZ), VY2(MVSIZ), VY3(MVSIZ), VY4(MVSIZ),
80 . VZ1(MVSIZ), VZ2(MVSIZ), VZ3(MVSIZ), VZ4(MVSIZ),
81 . AREA(MVSIZ),THK0(MVSIZ),VHX(MVSIZ), VHY(MVSIZ),
82 . shf(mvsiz),z2(mvsiz),viscmx(mvsiz),g(mvsiz),
83 . h1(mvsiz), h2(mvsiz), h3(mvsiz), a11(mvsiz),
84 . b11(mvsiz), b12(mvsiz), b13(mvsiz), b14(mvsiz),
85 . b21(mvsiz), b22(mvsiz), b23(mvsiz), b24(mvsiz),
86 . h11(mvsiz), h12(mvsiz), h13(mvsiz),
87 . h21(mvsiz), h22(mvsiz), h23(mvsiz),
88 . h31(mvsiz), h32(mvsiz), h33(mvsiz),
89 . rx1(mvsiz), rx2(mvsiz), rx3(mvsiz), rx4(mvsiz),
90 . ry1(mvsiz), ry2(mvsiz), ry3(mvsiz), ry4(mvsiz),
91 . thk02(mvsiz),ym(mvsiz), nu(mvsiz), alpe(mvsiz),
92 . srh1(*) ,srh2(*) ,srh3(*),a11r(*)
93C-----------------------------------------------
94C L o c a l V a r i a b l e s
95C-----------------------------------------------
96 INTEGER I,MX,J, II, IC, JST(MVSIZ+1)
97C REAL
98 my_real
99 . H1L(MVSIZ), H2L(MVSIZ), H3L(MVSIZ), H1Q(MVSIZ), H2Q(MVSIZ),
100 . H3Q(MVSIZ), HG1(MVSIZ), HG2(MVSIZ),
101 . HH1(MVSIZ), HH2(MVSIZ), HH3(MVSIZ),
102 . A1(MVSIZ), A2(MVSIZ), A3(MVSIZ), A4(MVSIZ),
103 . A5(MVSIZ), A6(MVSIZ), A7(MVSIZ), A8(MVSIZ),
104 . B1(MVSIZ), B2(MVSIZ), EHOU(MVSIZ),
105 .
106 . GAMA1(MVSIZ), GAMA2(MVSIZ), GAMA3(MVSIZ), GAMA4(MVSIZ)
107 my_real
108 . px1v, px2v, py1v, py2v, hvish1, hvish2,
109 . shfpr3, fac, ehourt,
110 . r0,r1,t2a,tsa,td,vv,
111 . hour1a,hour2a,hour3a,
112 .
113 . sr2d2,srshfpr3, inv9, inv12,scale(mvsiz)
114C=======================================================================
115 ehourt = zero
116 sr2d2 = sqrt(two)* half
117 inv9 = one_over_9
118 inv12 = one_over_12
119
120!! IGTYP = IGEO(11,MAT(1))
121C
122 IF(ismstr/=1.AND.ismstr/=11.AND.ihbe>=1)THEN
123 DO i=jft,jlt
124 px1v = px1(i)*vhx(i)
125 px2v = px2(i)*vhx(i)
126 py1v = py1(i)*vhy(i)
127 py2v = py2(i)*vhy(i)
128 gama1(i)= off(i)*( one - px1v-py1v)
129 gama3(i)= off(i)*( one + px1v+py1v)
130 gama2(i)= off(i)*(-one - px2v-py2v)
131 gama4(i)= off(i)*(-one + px2v+py2v)
132 ENDDO
133 ELSE
134 DO i=jft,jlt
135 gama1(i)= off(i)
136 gama3(i)= off(i)
137 gama2(i)= -off(i)
138 gama4(i)= -off(i)
139 ENDDO
140 ENDIF
141C----------------------------------------
142 DO i=jft,jlt
143 shfpr3 = shf(i)/(three*(one + nu(i)))
144 hvish1 = hvisc*h1(i)
145 hvish2 = hvisc*h2(i)
146 r0 = fourth*rho(i)
147 r1 = r0*hundred
148 r0 = r0*hvlin
149 a1(i) = r1*hvish1
150 a2(i) = r0*sr2d2*srh1(i)
151 srshfpr3 = sqrt(shfpr3)
152 a3(i) = r1*hvish2*srshfpr3
153 a4(i) = r0*sr2d2*srh2(i)*srshfpr3
154 hh3(i) = helas*h3(i)
155 a5(i) = hh3(i)*r1*zep072169
156 hh3(i) = sr2d2*srh3(i)
157 a6(i) = hh3(i)*r0*zep072169
158 r0 = fourth*ym(i)*helas
159 a7(i) = h1(i)*r0
160 a8(i) = h2(i)*r0*shfpr3
161 ENDDO
162C-------------------------------------
163 DO i=jft,jlt
164 t2a = thk02(i)*area(i)
165 tsa = sqrt(t2a)
166 h1q(i) = a1(i)*tsa
167 h1l(i) = a2(i)*ssp(i)*tsa
168 h2q(i) = a3(i)*thk02(i)
169 h2l(i) = a4(i)*ssp(i)*thk02(i)
170 h3q(i) = a5(i)*t2a
171 h3l(i) = a6(i)*ssp(i)*t2a
172 td = thk0(i)*dt1c(i)
173 hh1(i) = a7(i)*td
174 b1(i) = px1(i)*px1(i)+py1(i)*py1(i)
175 b2(i) = px2(i)*px2(i)+py2(i)*py2(i)
176 hh2(i) = a8(i)*thk02(i)*td/(b1(i)+b2(i))
177 ENDDO
178C-----------------------------------
179C TRIANGLES
180C-----------------------------------
181 DO i=jft,jlt
182 IF(ixc(4,i)==ixc(5,i))THEN
183 h1q(i)=zero
184 h1l(i)=zero
185 h2q(i)=zero
186 h2l(i)=zero
187 h3q(i)=zero
188 h3l(i)=zero
189 hh1(i)=zero
190 hh2(i)=zero
191 END IF
192 END DO
193C-----------------------------------
194C STIFFNESS - DT
195C-----------------------------------
196 DO i=jft,jlt
197 scale(i) = zero
198 ENDDO
199 IF(nodadt/=0.OR.idt1sh==1.OR.idtmins==2)THEN
200 DO i=jft,jlt
201 scale(i) = max(gama1(i)*gama1(i),gama2(i)*gama2(i),
202 . gama3(i)*gama3(i),gama4(i)*gama4(i)) *
203 . dt1c(i)*max(hh1(i)+h1l(i),hh2(i)+h2l(i),h3l(i))
204 . /max(dt1c(i)*dt1c(i),em20)
205 sti(i) = sti(i) + scale(i)
206 ENDDO
207 !! IF(IDTMINS==2 )IGMAT = 1
208C
209 IF(igtyp == 52 .OR.
210 . ((igtyp == 11 .OR. igtyp == 17 .OR. igtyp == 51)
211 . .AND. igmat > 0 )) THEN
212 IF(nadmesh==0)THEN
213 DO i=jft,jlt
214 IF (off(i)==zero) THEN
215 sti(i) = zero
216 stir(i) = zero
217 ELSE
218 vv = viscmx(i) * viscmx(i) * alpe(i)
219 fac = max(b1(i),b2(i))/(area(i) * vv)
220 sti(i) = sti(i) + fac * thk0(i) * a11(i)
221 stir(i) = fac * a11r(i)*one_over_12*thk0(i)**3 +
222 . fac * a11(i)*thk0(i)*area(i)*inv9 +
223 . fac*scale(i)*(one_over_12*thk0(i)**2 + area(i)*inv9)
224 ENDIF
225 ENDDO
226 ELSE
227 DO i=jft,jlt
228 IF (off(i)==zero) THEN
229 sti(i) = zero
230 stir(i) = zero
231 ELSE
232 vv = viscmx(i) * viscmx(i) * alpe(i)
233 fac = max(b1(i),b2(i))/(area(i) * vv)
234 sti(i) = sti(i) + fac * thk0(i) * a11(i)
235 stir(i) = fac * a11r(i)*one_over_12*thk0(i)**3 +
236 . fac*scale(i)*one_over_12*thk0(i)**2
237 ENDIF
238 ENDDO
239 END IF
240 ELSE
241 IF(nadmesh==0)THEN
242 DO i=jft,jlt
243 IF (off(i)==zero) THEN
244 sti(i) = zero
245 stir(i) = zero
246 ELSE
247 vv = viscmx(i) * viscmx(i) * alpe(i)
248 sti(i) = sti(i) + max(b1(i),b2(i))
249 . * thk0(i) * a11(i) / (area(i) * vv)
250 stir(i) = sti(i) * (thk02(i)*inv12 + area(i)*inv9)
251c .. + 0.5 * SHF(I) * AREA(I) * G(I)/A11(I))
252 ENDIF
253 ENDDO
254 ELSE
255 DO i=jft,jlt
256 IF (off(i)==zero) THEN
257 sti(i) = zero
258 stir(i) = zero
259 ELSE
260 vv = viscmx(i) * viscmx(i) * alpe(i)
261 sti(i) = sti(i) + max(b1(i),b2(i))
262 . * thk0(i) * a11(i) / (area(i) * vv)
263 stir(i) = sti(i) * thk02(i)*inv12
264 ENDIF
265 ENDDO
266 END IF
267 endif! IGTYP
268 ENDIF
269C-----------------------------------
270C ANTI-HOURGLASS MEMBRANE FORCES
271C-----------------------------------
272 IF(ismstr==1.OR.ismstr==11.OR.ihbe<1)THEN
273 DO i=jft,jlt
274 hg1(i)=(vx1(i)-vx2(i)+vx3(i)-vx4(i))*off(i)
275 hg2(i)=(vy1(i)-vy2(i)+vy3(i)-vy4(i))*off(i)
276 ENDDO
277 ELSE
278 DO i=jft,jlt
279 hg1(i)=vx1(i)*gama1(i)+vx2(i)*gama2(i)
280 . +vx3(i)*gama3(i)+vx4(i)*gama4(i)
281 hg2(i)=vy1(i)*gama1(i)+vy2(i)*gama2(i)
282 . +vy3(i)*gama3(i)+vy4(i)*gama4(i)
283 ENDDO
284 ENDIF
285 DO i=jft,jlt
286 hour(i,1)=hour(i,1)+hg1(i)*hh1(i)
287 hour(i,2)=hour(i,2)+hg2(i)*hh1(i)
288 hour1a =hour(i,1)+hg1(i)*(h1l(i)+h1q(i)*abs(hg1(i)))
289 h11(i)=hour1a*gama1(i)
290 h12(i)=hour1a*gama2(i)
291 h13(i)=hour1a*gama3(i)
292C H14(I)=HOUR1A*GAMA4(I)
293 hour2a =hour(i,2)+hg2(i)*(h1l(i)+h1q(i)*abs(hg2(i)))
294 h21(i)=hour2a*gama1(i)
295 h22(i)=hour2a*gama2(i)
296 h23(i)=hour2a*gama3(i)
297C H24(I)=HOUR2A*GAMA4(I)
298 ehou(i) = hour1a*hg1(i) + hour2a*hg2(i)
299 ENDDO
300C
301C----------------------------------
302C ANTI-HOURGLASS BENDING FORCES
303C----------------------------------
304 IF(ismstr==1.OR.ismstr==11.OR.ihbe<1)THEN
305 DO i=jft,jlt
306 hg1(i)=(vz1(i)-vz2(i)+vz3(i)-vz4(i))*off(i)
307 ENDDO
308 ELSE
309 DO i=jft,jlt
310 hg1(i)=vz1(i)*gama1(i)+vz2(i)*gama2(i)
311 . +vz3(i)*gama3(i)+vz4(i)*gama4(i)
312 ENDDO
313 ENDIF
314 DO i=jft,jlt
315 hour(i,3)=hour(i,3)+hg1(i)*hh2(i)
316 hour3a =hour(i,3)+hg1(i)*(h2l(i)+h2q(i)*abs(hg1(i)))
317 h31(i)=hour3a*gama1(i)
318 h32(i)=hour3a*gama2(i)
319 h33(i)=hour3a*gama3(i)
320C H34(I)=HOUR3(I)*GAMA4(I)
321 ehou(i) = ehou(i) + hour3a*hg1(i)
322 ENDDO
323C
324C------------------------
325C ANTI-HOURGLASS MOMENTS
326C------------------------
327 DO i=jft,jlt
328 hg1(i)=rx1(i)-rx2(i)+rx3(i)-rx4(i)
329 hg2(i)=ry1(i)-ry2(i)+ry3(i)-ry4(i)
330 ENDDO
331C
332 DO i=jft,jlt
333 hour(i,4)=hg1(i)*(h3l(i)+h3q(i)*abs(hg1(i)))
334 hour(i,5)=hg2(i)*(h3l(i)+h3q(i)*abs(hg2(i)))
335 ehou(i) = ehou(i) +
336 . hour(i,4)*hg1(i) + hour(i,5)*hg2(i)
337 ehou(i) = dt1c(i) * ehou(i) * off(i)
338C EHOURT = EHOURT + EHOU(I)
339 b11(i)= hour(i,4)*off(i)
340 b12(i)=-hour(i,4)*off(i)
341 b13(i)= hour(i,4)*off(i)
342 b14(i)=-hour(i,4)*off(i)
343 b21(i)= hour(i,5)*off(i)
344 b22(i)=-hour(i,5)*off(i)
345 b23(i)= hour(i,5)*off(i)
346 b24(i)=-hour(i,5)*off(i)
347 ENDDO
348C separate loop that does not vectorize w/o simd directive
349 DO i=jft,jlt
350 ehourt = ehourt + ehou(i)
351 ENDDO
352C
353 ic=1
354 jst(ic)=jft
355 DO j=jft+1,jlt
356 IF (ipartc(j)/=ipartc(j-1)) THEN
357 ic=ic+1
358 jst(ic)=j
359 ENDIF
360 ENDDO
361
362 jst(ic+1)=jlt+1
363 IF(ic==1) THEN
364 mx = ipartc(jft)
365 DO i=jft,jlt
366 partsav(8,mx)=partsav(8,mx) + ehou(i)
367 ENDDO
368
369 ELSEIF(ic==2.AND.kfts>0)THEN
370 mx = ipartc(jft)
371 DO i=jft, kfts-1
372 partsav(8,mx)=partsav(8,mx) + ehou(i)
373 ENDDO
374 mx = ipartc(jlt)
375 DO i=kfts,jlt
376 partsav(8,mx)=partsav(8,mx) + ehou(i)
377 ENDDO
378
379 ELSE
380C
381 DO ii=1,ic
382 mx=ipartc(jst(ii))
383 IF (jst(ii+1)-jst(ii)>15) THEN
384 DO j=jst(ii),jst(ii+1)-1
385 partsav(8,mx)=partsav(8,mx)+ehou(j)
386 ENDDO
387 ELSE
388 DO j=jst(ii),jst(ii+1)-1
389 partsav(8,mx)=partsav(8,mx)+ehou(j)
390 ENDDO
391 ENDIF
392 ENDDO
393 ENDIF
394C
395!$OMP ATOMIC
396 ehour = ehour + ehourt
397C
398 DO i=jft,jlt
399 eani(nft+numels+i) = eani(nft+numels+i) + ehou(i)
400 ENDDO
401C
402 RETURN
403 END
chvis3
subroutine chvis3(jft, jlt, thk, hour, off, px1, px2, py1, py2, ixc, dt1c, ssp, rho, sti, vx1, vx2, vx3, vx4, vy1, vy2, vy3, vy4, vz1, vz2, vz3, vz4, area, thk0, vhx, vhy, shf, z2, eani, stir, viscmx, g, a11, h1, h2, h3, ym, nu, thk02, alpe, h11, h12, h13, h21, h22, h23, h31, h32, h33, b11, b12, b13, b14, b21, b22, b23, b24, rx1, rx2, rx3, rx4, ry1, ry2, ry3, ry4, ipartc, partsav, ihbe, nft, ismstr, kfts, srh1, srh2, srh3, igtyp, igmat, a11r, nel)
Definition chvis3.F:44
max
#define max(a, b)
Definition macros.h:21