OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
static.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!|| static ../engine/source/general_controls/damping/static.F
25!||--- called by ------------------------------------------------------
26!|| resol ../engine/source/engine/resol.F
27!||--- calls -----------------------------------------------------
28!|| ngr2usr ../engine/source/input/freform.F
29!||--- uses -----------------------------------------------------
30!|| groupdef_mod ../common_source/modules/groupdef_mod.F
31!||====================================================================
32 SUBROUTINE static(V,VR,A,AR,MS,IN,IGRNOD,WEIGHT_MD,WFEXT)
33C-----------------------------------------------
34C M o d u l e s
35C-----------------------------------------------
36 USE groupdef_mod
37C-----------------------------------------------
38C I m p l i c i t T y p e s
39C-----------------------------------------------
40#include "implicit_f.inc"
41C-----------------------------------------------
42C C o m m o n B l o c k s
43C-----------------------------------------------
44#include "com01_c.inc"
45#include "com04_c.inc"
46#include "com08_c.inc"
47#include "scr11_c.inc"
48#include "statr_c.inc"
49#include "stati_c.inc"
50#include "task_c.inc"
51#include "sms_c.inc"
52C-----------------------------------------------
53C D u m m y A r g u m e n t s
54C-----------------------------------------------
55 INTEGER WEIGHT_MD(*)
56 my_real v(3,*), vr(3,*), a(3,*), ar(3,*),ms(*),in(*)
57 DOUBLE PRECISION,INTENT(INOUT) :: WFEXT
58C-----------------------------------------------
59 TYPE (GROUP_) , DIMENSION(NGRNOD) :: IGRNOD
60C-----------------------------------------------
61C L o c a l V a r i a b l e s
62C-----------------------------------------------
63 INTEGER J, I,LAST,N,NGR2USR,ISTAT2,ISTAT3
64 my_real encino, omega, uomega, domega, omega2, encint, encinn,encgrp,encgrpn, encgrpr,encgrprn
65 EXTERNAL ngr2usr
66C-----------------------------------------------
67 DATA encino /0.0/
68 DATA last /0/
69 SAVE encino,last
70C
71 IF(istatg<0)THEN
72 istatg=ngr2usr(-istatg,igrnod,ngrnod)
73 ENDIF
74 istat2=0
75 istat3=0
76 IF (ncycle==0) THEN
77 IF ((istat==2.OR.istat==3).AND.tst_stop==zero) tst_stop=tstop
78 END IF
79 IF (istat==2.AND.tt>=tst_start.AND.tt<=tst_stop) istat2=1
80 IF (istat==3.AND.tt>=tst_start.AND.tt<=tst_stop) istat3=1
81 IF((istat==1.OR.istat3==1).AND.istatg==0)THEN
82C
83 omega = betate * dt12
84 uomega = one - omega
85 domega = two*betate
86 omega2 = (one-two*omega)**2
87 encint = encin
88C ENCINT = ENCIN+ENROT
89 encinn = encint * omega2
90C
91 IF(idtmins==0) THEN
92 DO i = 1, numnod
93 a(1,i) = -domega*v(1,i) + uomega*a(1,i)
94 a(2,i) = -domega*v(2,i) + uomega*a(2,i)
95 a(3,i) = -domega*v(3,i) + uomega*a(3,i)
96 ENDDO
97 IF (ispmd==0) wfext = wfext - encint + encinn
98 END IF !(IDTMINS==0) THEN
99C
100 IF(iroddl/=0)THEN
101 encint = enrot
102 encinn = encint * omega2
103 DO i = 1, numnod
104 ar(1,i) = -domega*vr(1,i) + uomega*ar(1,i)
105 ar(2,i) = -domega*vr(2,i) + uomega*ar(2,i)
106 ar(3,i) = -domega*vr(3,i) + uomega*ar(3,i)
107 ENDDO
108 IF (ispmd==0) wfext = wfext - encint + encinn
109 ENDIF
110C
111 ELSEIF(istat2==1.AND.istatg==0)THEN
112C
113 last = last+1
114 encint = encin+enrot
115 IF(encint<encino)THEN
116 IF (ispmd==0)wfext = wfext - encint
117 encino=zero
118 last=0
119C
120 DO j=1,3
121 DO i = 1, numnod
122 v(j,i) = zero
123 ENDDO
124 IF(iroddl/=0)THEN
125 DO i = 1, numnod
126 vr(j,i) = zero
127 ENDDO
128 ENDIF
129 ENDDO
130 ELSE
131 encino = encint
132 ENDIF
133 ELSEIF((istat==1.OR.istat3==1).AND.istatg/=0)THEN
134C
135 omega = betate * dt12
136 uomega = one - omega
137 domega = two*betate
138 omega2 = (one-two*omega)**2
139 encint = encin+enrot
140 encinn = encint * omega2
141 encgrp = zero
142 encgrpr = zero
143C
144 IF(idtmins==0) THEN
145 DO j=1,3
146#include "vectorize.inc"
147 DO n=1,igrnod(istatg)%NENTITY
148 i=igrnod(istatg)%ENTITY(n)
149 encgrp = encgrp + ms(i)*weight_md(i)*
150 . v(j,i)*v(j,i)
151 a(j,i) = -domega*v(j,i) + uomega*a(j,i)
152 ENDDO
153 ENDDO
154 END IF !(IDTMINS==0) THEN
155C
156 IF(iroddl/=0)THEN
157 DO j=1,3
158#include "vectorize.inc"
159 DO n=1,igrnod(istatg)%NENTITY
160 i=igrnod(istatg)%ENTITY(n)
161 encgrpr = encgrpr + ms(i)*weight_md(i)*
162 . vr(j,i)*vr(j,i)
163 ar(j,i) = -domega*vr(j,i) + uomega*ar(j,i)
164 ENDDO
165 ENDDO
166 ENDIF
167C
168 encgrp = half*encgrp
169 encgrpr = half*encgrpr
170 encgrpn = encgrp * omega2
171 encgrprn = encgrpr * omega2
172 IF (ispmd==0)wfext = wfext + encgrpn - encgrp + encgrprn - encgrpr
173C
174 ELSEIF(istat2==1.AND.istatg/=0)THEN
175C
176 last = last+1
177 encint = encin+enrot
178 IF(encint<encino)THEN
179 IF (ispmd==0)wfext = wfext - encint
180 encino=zero
181 last=0
182C
183 DO j=1,3
184#include "vectorize.inc"
185 DO n=1,igrnod(istatg)%NENTITY
186 i=igrnod(istatg)%ENTITY(n)
187 v(j,i) = zero
188 ENDDO
189 IF(iroddl/=0)THEN
190#include "vectorize.inc"
191 DO n=1,igrnod(istatg)%NENTITY
192 i=igrnod(istatg)%ENTITY(n)
193 vr(j,i) = zero
194 ENDDO
195 ENDIF
196 ENDDO
197 ELSE
198 encino = encint
199 ENDIF
200 ENDIF
201C
202 RETURN
203 END
204!||====================================================================
205!|| e_period ../engine/source/general_controls/damping/static.F
206!||--- called by ------------------------------------------------------
207!|| ener_w0 ../engine/source/general_controls/damping/static.F
208!||--- calls -----------------------------------------------------
209!|| butterworth ../engine/source/tools/univ/butterworth.F
210!||====================================================================
211 SUBROUTINE e_period(IPI,IPC,F_FIL)
212C-----------------------------------------------
213C I m p l i c i t T y p e s
214C-----------------------------------------------
215#include "implicit_f.inc"
216C-----------------------------------------------
217C D u m m y A r g u m e n t s
218C-----------------------------------------------
219 INTEGER IPI,IPC
220C
221 my_real
222 . f_fil
223C-----------------------------------------------
224C C o m m o n B l o c k s
225C-----------------------------------------------
226#include "com08_c.inc"
227#include "scr11_c.inc"
228#include "statr_c.inc"
229C-----------------------------------------------
230C L o c a l V a r i a b l e s
231C-----------------------------------------------
232 INTEGER J, I , IDF
233C REAL
234 my_real
235 . encint, eint,vv1,vv2,fv1,fv2,fv,ff
236C-----------------------------------------------
237 IF (tt == zero) THEN
238 ff= encin+enrot
239 fil_ke(1) = ff
240 fil_ke(2) = ff
241 fil_ke(3) = ff
242 fil_ke(4) = ff
243 ff= enint
244 fil_ie(1) = ff
245 fil_ie(2) = ff
246 fil_ie(3) = ff
247 fil_ie(4) = ff
248 END IF
249 IF (f_fil>zero) THEN
250 ff = encin+enrot
251 vv1 = fil_ke(1)
252 vv2 = fil_ke(2)
253 fv1 = fil_ke(3)
254 fv2 = fil_ke(4)
255 CALL butterworth(dt2,f_fil,vv2,vv1,ff,fv2,fv1,fv)
256 fil_ke(1)= ff
257 fil_ke(2)= vv1
258 fil_ke(3)= fv
259 fil_ke(4)= fv1
260 encint = fv
261C
262 ff = enint
263 vv1 = fil_ie(1)
264 vv2 = fil_ie(2)
265 fv1 = fil_ie(3)
266 fv2 = fil_ie(4)
267 CALL butterworth(dt2,f_fil,vv2,vv1,ff,fv2,fv1,fv)
268 fil_ie(1) = ff
269 fil_ie(2) = vv1
270 fil_ie(3) = fv
271 fil_ie(4) = fv1
272 eint = fv
273C
274 ELSE
275 encint = encin+enrot
276 eint = enint
277 END IF
278 pcin = pcin +dt1
279 pint = pint +dt1
280 pimax = max(pimax,pint)
281 pcmax = max(pcmax,pcin)
282 IF(encint<encin_0.AND.encint>=zero)THEN
283 encin_0=zero
284 pcin = zero
285 ipc = 1
286 ELSE
287 ipc = 0
288 encin_0 = encint
289 ENDIF
290 IF(eint<zero) THEN
291 eint_0=zero
292 pint = zero
293 ipi = -2
294 ELSEIF(eint<eint_0.AND.eint>=zero)THEN
295 eint_0=zero
296 pint = zero
297 ipi = 1
298 ELSE
299 eint_0 = eint
300 ipi = 0
301 ENDIF
302C
303 RETURN
304 END
305!||====================================================================
306!|| ener_w0 ../engine/source/general_controls/damping/static.F
307!||--- called by ------------------------------------------------------
308!|| resol ../engine/source/engine/resol.F
309!||--- calls -----------------------------------------------------
310!|| e_period ../engine/source/general_controls/damping/static.F
311!||====================================================================
312 SUBROUTINE ener_w0
313C-----------------------------------------------
314C I m p l i c i t T y p e s
315C-----------------------------------------------
316#include "implicit_f.inc"
317C-----------------------------------------------
318C D u m m y A r g u m e n t s
319C-----------------------------------------------
320C
321C-----------------------------------------------
322C C o m m o n B l o c k s
323C-----------------------------------------------
324#include "com01_c.inc"
325#include "com08_c.inc"
326#include "statr_c.inc"
327#include "task_c.inc"
328#include "sms_c.inc"
329#include "warn_c.inc"
330#include "scr05_c.inc"
331#include "scr07_c.inc"
332C-----------------------------------------------
333C L o c a l V a r i a b l e s
334C-----------------------------------------------
335 INTEGER J, I,FREQ,IDF,IPI,IPC,NC_ACT,IFIRST
336 parameter(nc_act = 200)
337C REAL
338 my_real
339 . q_es,betate_n,fc,fi,betate_m
340 my_real
341 . f_max,f_min,f_0,ei_tol
342C REAL
343 CHARACTER*8 FILNAME
344C-----------------------------------------------
345 CALL e_period(ipi,ipc,freq_c)
346 IF (ncycle==0) THEN
347 betate = betate_0
348 IF (debug(11)==1.AND.ispmd==0) THEN
349 idf = 251
350 filname='P_ES.TMP'
351 OPEN(unit=idf,file=filname,status='UNKNOWN',form='FORMATTED')
352 write(idf,*)
353 . '# NCYCLE OMEGA(Retenu)IFIRST OMEGA_N(new) P_Eint P_Kin'
354 END IF
355C--to avoid division by zero /DT1
356 RETURN
357 END IF
358C------- /DEBUG/ADYREL for periods print out
359 IF (debug(11)==1.AND.ispmd==0) THEN
360 idf = 251
361 ELSE
362 idf = 0
363 END IF
364 ifirst = nint(nfirst)
365 ei_tol=em12
366 IF (iresp==1) ei_tol=em07
367C-----0.5%w_max 5*0.5 for AMS
368 IF(idtmins==0) THEN
369 f_max=em02/dt1
370 ELSE
371 f_max=zep05/dt1
372 END IF
373 IF(freq_c<zero) freq_c = -freq_c*f_max
374 f_0=em02*f_max
375 betate_n = f_max
376 IF (ipi==1) THEN
377 fi = one/pimax
378 betate_n = min(f_max,fi)
379 IF (betate==zero.AND.eint_0>ei_tol) THEN
380 IF (ncycle>=nc_act) THEN
381 betate = min(f_0,betate_n)
382 ifirst = 1
383 END IF
384C------exception for increasing----
385 ELSEIF (ifirst==1) THEN
386 betate = betate_n
387 ifirst = ifirst+1
388 ELSE
389 betate = min(betate,betate_n)
390 END IF
391 IF (idf>0) write(idf,1000)-ncycle,betate,ifirst,fi,pimax ,pcmax
392 END IF
393 IF (ipc==1) THEN
394 fc = one/pcmax
395 betate_n = min(f_max,fc)
396 IF (betate==zero.AND.eint_0>ei_tol) THEN
397 IF (ncycle>=nc_act) THEN
398 betate = min(f_0,betate_n)
399 ifirst = 1
400 END IF
401 ELSEIF (ifirst==1) THEN
402 betate = betate_n
403 ifirst = ifirst+1
404 ELSE
405C----limiting the using of PCMAX not take it if Ke/Eint is very small
406 IF (encin_0/max(em20,eint_0)>em03.AND.betate<betate_n*three_half) THEN
407 betate_m = half*betate
408 betate = min(betate,betate_n)
409 betate = max(betate,betate_m)
410 END IF
411 END IF
412 IF (idf>0) write(idf,1000)ncycle,betate,ifirst,fc,pimax ,pcmax
413 END IF
414C--Case high periods for both Eint and Ke
415 IF (betate==zero.AND.eint_0>ei_tol.AND.ncycle>=nc_act) THEN
416C------initial frequency is 1% of F_MAX otherwise damping too high especially w/ gravity
417 betate=f_0
418 fi = one/pimax
419 ifirst = 1
420 IF (idf>0) write(idf,1000)ncycle,betate,ifirst,fi,pimax ,pcmax
421 END IF
422C--Case BETATE too high w/ large periods for both Eint and Ke
423 IF (ifirst>=1.AND.(ipc+ipi)==0) THEN
424 fi = one/max(pimax,pcmax)
425 IF (betate>1.1*fi) THEN
426 betate=fi
427 IF (idf>0) write(idf,1000)ncycle,betate,-ifirst,fi,pimax ,pcmax
428 END IF
429 END IF
430 IF (ifirst==1.AND.irun>1.AND.mcheck==0) ifirst = 2
431 nfirst = ifirst
432C
433 RETURN
434 1000 FORMAT(i8,5(g14.7)/)
435 END
subroutine butterworth(dt, freq, x2, x1, x, fx2, fx1, fx)
Definition butterworth.F:31
#define my_real
Definition cppsort.cpp:32
#define min(a, b)
Definition macros.h:20
#define max(a, b)
Definition macros.h:21
subroutine e_period(ipi, ipc, f_fil)
Definition static.F:212
subroutine static(v, vr, a, ar, ms, in, igrnod, weight_md, wfext)
Definition static.F:33
subroutine ener_w0
Definition static.F:313