OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
joint_block_stiffness.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!|| joint_block_stiffness ../engine/source/elements/joint/joint_block_stiffness.F
25!||--- called by ------------------------------------------------------
26!|| resol ../engine/source/engine/resol.F
27!||--- calls -----------------------------------------------------
28!|| arret ../engine/source/system/arret.F
29!||--- uses -----------------------------------------------------
30!|| elbufdef_mod ../common_source/modules/mat_elem/elbufdef_mod.F90
31!|| element_mod ../common_source/modules/elements/element_mod.F90
32!||====================================================================
33 SUBROUTINE joint_block_stiffness(ITAB,MS,IN,STIFN,STIFR,
34 1 WEIGHT ,IXR,IPART,X,
35 2 IPARTR,IGEO,GEO,NPBY,IPARG,ELBUF_TAB,DMAS,DINER)
36C-----------------------------------------------
37C M o d u l e s
38C-----------------------------------------------
39 USE elbufdef_mod
40 use element_mod , only : nixr
41C----6---------------------------------------------------------------7---------8
42C I m p l i c i t T y p e s
43C-----------------------------------------------
44#include "implicit_f.inc"
45C-----------------------------------------------
46C G l o b a l P a r a m e t e r s
47C-----------------------------------------------
48#include "com08_c.inc"
49#include "com01_c.inc"
50#include "param_c.inc"
51#include "scr02_c.inc"
52#include "scr17_c.inc"
53#include "scr18_c.inc"
54#include "cong2_c.inc"
55#include "units_c.inc"
56C-----------------------------------------------------------------
57C D u m m y A r g u m e n t s
58C-----------------------------------------------
59 INTEGER ITAB(*),WEIGHT(*),IXR(NIXR,*),
60 . IPART(LIPART1,*),IPARTR(*),IGEO(NPROPGI,*),NPBY(NNPBY,*),
61 . iparg(nparg,*)
62C REAL
63 my_real stifn(*), stifr(*),ms(*) ,in(*),x(3,*),
64 . dmas,diner,geo(npropg,*)
65C
66 TYPE(elbuf_struct_), TARGET, DIMENSION(NGROUP) :: ELBUF_TAB
67C-----------------------------------------------
68C L o c a l V a r i a b l e s
69C-----------------------------------------------
70 INTEGER I,J,M1,M2,IG,IGTYP,N1,N2,ITYP,NG,JFT,JLT,NEL,
71 . FLAG,NFT,NUVAR,JNTYP,IRB1,IRB2,FLAG_S,FLAG_PR,NV
72 my_real iner1,iner2,km1,krm1,km2,krm2,xx,kx1,kx2,kr1,dtc,alpha,
73 . xl,kxmax,krmax,kx,kr,scf,get_u_geo,
74 . xx1,xx2,dta,mass1,mass2
75 CHARACTER(LEN=30) :: INFO_MESSAGE
76C
77 TYPE(g_bufel_),POINTER :: GBUF
78C-----------------------------------------------
79 EXTERNAL get_u_geo
80C----------------------------------------------------------
81
82 flag_pr = 0
83
84 dtc = dt2
85 dta = dtmin1(11)
86
87 IF (tt==0) THEN
88C----------------------------------------------------------
89 IF (nodadt>0) THEN
90 dtc = max(dtc,dtmin1(11)/dtfac1(11))
91 IF (dtmx>em20) dtc = min(dtc,dtmx)
92 ENDIF
93
94 IF (dtc>90000) THEN
95 IF (dta==0) THEN
96 print *,"ERROR NO TARGET TIME STEP DT=",dtc
97 print *,"STIFFNESS CAN NOT BE COMPUTED"
98 CALL arret(2)
99 ELSE
100 dtc = dta
101 ENDIF
102 ENDIF
103
104 alpha = 4.1
105
106 DO ng=1,ngroup
107 ityp = iparg(5,ng)
108 nel = iparg(2,ng)
109 nft = iparg(3,ng)
110 jft = 1
111 jlt = min(nvsiz,nel)
112 gbuf => elbuf_tab(ng)%GBUF
113 IF (ityp /= 6) GOTO 250
114C--------> Loop over spring elements -------
115 DO i=jft,jlt
116 j = i + nft
117 ig = ipart(2,ipartr(j))
118 igtyp = igeo(11,ig)
119 nuvar = nint(geo(25,ig))
120 nv = nuvar*(i-1) + 1
121 IF (igtyp==45) THEN
122 scf = get_u_geo(11,ig)
123 jntyp = nint(get_u_geo(1,ig))
124 flag = nint(get_u_geo(10,ig))
125 flag_s = 0
126 IF (flag==0) THEN
127C
128 n1 = ixr(2,j)
129 n2 = ixr(3,j)
130C
131 irb1 = nint(gbuf%VAR(nv + 37))
132 IF (irb1 > 0) THEN
133 m1 = npby(1,irb1)
134 ELSE
135 m1 = n1
136 ENDIF
137C
138 irb2 = nint(gbuf%VAR(nv + 38))
139 IF (irb2 > 0) THEN
140 m2 = npby(1,irb2)
141 ELSE
142 m2 = n2
143 ENDIF
144C
145 xl = ((x(1,n1)-x(1,n2))**2)+((x(2,n1)-x(2,n2))**2)
146 . +((x(3,n1)-x(3,n2))**2)
147 xx1 = ((x(1,n1)-x(1,m1))**2)+((x(2,n1)-x(2,m1))**2)
148 . +((x(3,n1)-x(3,m1))**2)
149 xx2 = ((x(1,n2)-x(1,m2))**2)+((x(2,n2)-x(2,m2))**2)
150 . +((x(3,n2)-x(3,m2))**2)
151 xx = max(xx1,xx2)
152C
153C--------> Stiffness calculation node 1 side -------
154 mass1 = ms(m1)
155 iner1 = in(m1)
156 km1 = stifn(m1)
157 krm1 = stifr(m1)
158C
159 kx1 = (2*mass1/(alpha*dtc*dtc)) - km1
160 IF (iner1 > zero) THEN
161 kx2 = 0.8*(iner1/(alpha*dtc*dtc)- krm1)/(max(em20,(xx+xl)))
162 kr = iner1/(alpha*dtc*dtc)- krm1
163 ELSE
164 kx2 = ep30
165 kr = zero
166 ENDIF
167 kxmax = min(kx1,kx2)
168C
169C--------> Stiffness calculation node 2 side -------
170 mass2 = ms(m2)
171 iner2 = in(m2)
172 km2 = stifn(m2)
173 krm2 = stifr(m2)
174C
175 kx1 = (2*mass2/(alpha*dtc*dtc)) - km2
176 IF (iner2 > zero) THEN
177 kx2 = 0.8*(iner2/(alpha*dtc*dtc)- krm2)/(max(em20,(xx+xl)))
178 kr1 = iner2/(alpha*dtc*dtc)- krm2
179 ELSE
180 kx2 = ep30
181 kr1 = zero
182 ENDIF
183C
184 kxmax = min(kx1,kx2,kxmax)
185 kr = min(kr,kr1)
186
187C--------> Final stiffness calculation and assignment-------
188 kx = max(kxmax,2*km1,2*km2)
189 IF ((kx - kxmax)>1e-8) flag_s = 1
190 IF ((iner1 == zero).OR.(iner2 == zero)) THEN
191C-- Rotational stiffness set to zero if kjoint attached to a solid node --
192 kr = zero
193 krmax = kr
194 flag_s = 0
195 ELSE
196 krmax = kr
197 kr = max(kr,2*krm1,2*krm2)
198 ENDIF
199C
200 info_message = ''
201 IF ((kx-kxmax)>1e-8 .OR. (kr-krmax)>1e-8) THEN
202 IF ((irb1 > 0).AND.(irb2 > 0)) THEN
203 info_message = '(FROM RIGID BODIES)'
204 ELSE
205 info_message = '(FROM CONNECTED STRUCTURES)'
206 ENDIF
207 ENDIF
208C
209 IF (flag_s==1) THEN
210 kr = max(kr,1.3*kx*(xx+xl))
211 ENDIF
212C
213 kx = scf*kx
214 kr = scf*kr
215C
216 WRITE(iout,'(i10,i10,8x,1pe11.4,8x,1pe11.4,1x,a)')
217 . IXR(NIXR,J),JNTYP,KX,KR,INFO_MESSAGE
218C
219 GBUF%VAR(NV + 16) = KX
220 GBUF%VAR(NV + 17) = KR
221
222C--------> Spherical Joint-------
223 IF (JNTYP==1) THEN
224 GBUF%VAR(NV + 18) = KX
225 GBUF%VAR(NV + 19) = KX
226 GBUF%VAR(NV + 20) = KX
227C--------> Revolutional Joint-------
228 ELSEIF (JNTYP==2) THEN
229 GBUF%VAR(NV + 18) = KX
230 GBUF%VAR(NV + 19) = KX
231 GBUF%VAR(NV + 20) = KX
232 GBUF%VAR(NV + 31) = KR
233 GBUF%VAR(NV + 32) = KR
234C--------> Cylindrical Joint-------
235 ELSEIF (JNTYP==3) THEN
236 GBUF%VAR(NV + 19) = KX
237 GBUF%VAR(NV + 20) = KX
238 GBUF%VAR(NV + 31) = KR
239 GBUF%VAR(NV + 32) = KR
240C--------> Planar Joint-------
241 ELSEIF (JNTYP==4) THEN
242 GBUF%VAR(NV + 18) = KX
243 GBUF%VAR(NV + 31) = KR
244 GBUF%VAR(NV + 32) = KR
245C--------> Universal Joint-------
246 ELSEIF (JNTYP==5) THEN
247 GBUF%VAR(NV + 18) = KX
248 GBUF%VAR(NV + 19) = KX
249 GBUF%VAR(NV + 20) = KX
250 GBUF%VAR(NV + 30) = KR
251C--------> Translational Joint-------
252 ELSEIF (JNTYP==6) THEN
253 GBUF%VAR(NV + 19) = KX
254 GBUF%VAR(NV + 20) = KX
255 GBUF%VAR(NV + 30) = KR
256 GBUF%VAR(NV + 31) = KR
257 GBUF%VAR(NV + 32) = KR
258C--------> Oldham Joint-------
259 ELSEIF (JNTYP==7) THEN
260 GBUF%VAR(NV + 18) = KX
261 GBUF%VAR(NV + 30) = KR
262 GBUF%VAR(NV + 31) = KR
263 GBUF%VAR(NV + 32) = KR
264C--------> Rigid Joint-------
265 ELSEIF (JNTYP==8) THEN
266 GBUF%VAR(NV + 18) = KX
267 GBUF%VAR(NV + 19) = KX
268 GBUF%VAR(NV + 20) = KX
269 GBUF%VAR(NV + 30) = KR
270 GBUF%VAR(NV + 31) = KR
271 GBUF%VAR(NV + 32) = KR
272 ENDIF ! IF (JNTYP)
273C
274 ENDIF ! IF (FLAG==0)
275 ENDIF ! IF (IGTYP==45)
276C
277 ENDDO
278250 CONTINUE
279 ENDDO
280
281C----------------------------------------------------------
282
283 ENDIF
284C
285 RETURN
286C
287 END
my_real
#define my_real
Definition cppsort.cpp:32
alpha
#define alpha
Definition eval.h:35
joint_block_stiffness
subroutine joint_block_stiffness(itab, ms, in, stifn, stifr, weight, ixr, ipart, x, ipartr, igeo, geo, npby, iparg, elbuf_tab, dmas, diner)
Definition joint_block_stiffness.F:36
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21
arret
subroutine arret(nn)
Definition arret.F:86