OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
mdtsph.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!|| mdtsph ../engine/source/materials/mat_share/mdtsph.F
25!||--- called by ------------------------------------------------------
26!|| m1law ../engine/source/materials/mat/mat001/m1law.F
27!|| m1lawi ../engine/source/materials/mat/mat001/m1lawi.F
28!|| m1lawtot ../engine/source/materials/mat/mat001/m1lawtot.f
29!|| m22law ../engine/source/materials/mat/mat022/m22law.F
30!|| m24law ../engine/source/materials/mat/mat024/m24law.F
31!|| m2law ../engine/source/materials/mat/mat002/m2law.F
32!|| mmain ../engine/source/materials/mat_share/mmain.F90
33!|| mulaw ../engine/source/materials/mat_share/mulaw.F90
34!|| sboltlaw ../engine/source/elements/solid/solide/sboltlaw.F
35!|| usermat_solid ../engine/source/materials/mat_share/usermat_solid.F
36!||====================================================================
37 SUBROUTINE mdtsph(
38 1 PM, OFF, RHO, RK,
39 2 T, RE, STI, DT2T,
40 3 NELTST, ITYPTST, OFFG, GEO,
41 4 PID, MUMAX, SSP, VOL,
42 5 VD2, DELTAX, VIS, D1,
43 6 D2, D3, PNEW, PSH,
44 7 MAT, NGL, QVIS, SSP_EQ,
45 8 G_DT, DTSPH, NEL, ITY,
46 9 JTUR, JTHE)
47C-----------------------------------------------
48C I m p l i c i t T y p e s
49C-----------------------------------------------
50#include "implicit_f.inc"
51#include "comlock.inc"
52C-----------------------------------------------
53C G l o b a l P a r a m e t e r s
54C-----------------------------------------------
55#include "mvsiz_p.inc"
56C-----------------------------------------------
57C C o m m o n B l o c k s
58C-----------------------------------------------
59#include "com08_c.inc"
60#include "scr02_c.inc"
61#include "scr18_c.inc"
62#include "param_c.inc"
63#include "cong1_c.inc"
64#include "units_c.inc"
65#include "scr07_c.inc"
66C-----------------------------------------------
67C D u m m y A r g u m e n t s
68C-----------------------------------------------
69 INTEGER, INTENT(IN) :: NEL
70 INTEGER, INTENT(IN) :: ITY
71 INTEGER, INTENT(IN) :: JTUR
72 INTEGER, INTENT(IN) :: JTHE
73 INTEGER NELTST,ITYPTST,PID(*),MAT(*), NGL(*)
74 my_real :: DT2T
75
76 my_real
77 . PM(NPROPM,*), OFF(*), RHO(*), RK(*), T(*),
78 . re(*),sti(*),offg(*),geo(npropg,*),mumax(*),
79 . vol(*), vd2(*), deltax(*), ssp(*), vis(*),
80 . psh(*), pnew(*),qvis(*) ,ssp_eq(*), d1(*),
81 . d2(*), d3(*)
82 my_real, INTENT(INOUT) :: dtsph(1:nel)
83 INTEGER,INTENT(IN) :: G_DT
84C-----------------------------------------------
85C L o c a l V a r i a b l e s
86C-----------------------------------------------
87 INTEGER I, MT
88 my_real
89 . AL(MVSIZ),DTX(MVSIZ), QX(MVSIZ), CX(MVSIZ), QXMATER(MVSIZ),
90 . QA, QB, VISI, FACQ,
91 . cns1, cns2, sph, ak1, bk1, ak2, bk2, tli, akk, xmu, tmu, rpr,
92 . atu
93C-----------------------------------------------
94C S o u r c e L i n e s
95C-----------------------------------------------
96 IF(impl==zero)THEN
97 DO i=1,nel
98 cx(i)=ssp(i)+sqrt(vd2(i))
99 ENDDO
100 visi=one
101 facq=one
102 ELSE
103 DO i=1,nel
104 cx(i)=sqrt(vd2(i))
105 ENDDO
106 visi=zero
107 facq=zero
108 ENDIF
109
110 !not a bug : only law 24 uses CNS1 & CNS2
111 !(they are not yet available with SPH).
112 DO i=1,nel
113 al(i)=zero
114 IF(off(i)<1.) cycle
115 al(i)=vol(i)**third
116 ENDDO
117
118 mt = mat(1)
119 DO i=1,nel
120 qa =facq*geo(14,pid(i))
121 qb =facq*geo(15,pid(i))
122 cns1=geo(16,pid(i))
123 cns2=geo(17,pid(i))*ssp(i)*al(i)*rho(i)
124 psh(i)=pm(88,mt)
125 pnew(i)=zero
126 qxmater(i)=cns1*ssp(i) + visi*(two*vis(i)+cns2) / max(em20,rho(i)*deltax(i))
127 qx(i)=qb*ssp(i) + qa*mumax(i) + qxmater(i)
128 qvis(i)=zero
129 ENDDO
130
131 DO i=1,nel
132 dtx(i)=deltax(i)/max(em20,qx(i)+sqrt(qx(i)*qx(i)+cx(i)*cx(i)))
133 !preparing material sound speed for nodal time step computation:
134 ssp_eq(i) = max(em20,qxmater(i)+sqrt(qxmater(i)*qxmater(i)+cx(i)*cx(i)))
135 ENDDO
136
137 IF(jthe/=0)THEN
138 mt = mat(1)
139 sph = pm(69,mt)
140 ak1 = pm(75,mt)
141 bk1 = pm(76,mt)
142 ak2 = pm(77,mt)
143 bk2 = pm(78,mt)
144 tli = pm(80,mt)
145 DO i=1,nel
146 IF(t(i)<tli)THEN
147 akk=ak1+bk1*t(i)
148 ELSE
149 akk=ak2+bk2*t(i)
150 ENDIF
151 IF(jtur/=0)THEN
152 xmu = rho(i)*pm(24,mt)
153 tmu = pm(81,mt)
154 rpr = pm(95,mt)
155 atu=rpr*tmu*rk(i)*rk(i)/(max(em15,re(i)*vol(i))*xmu)
156 akk=akk*(one + atu)
157 ENDIF
158 dtx(i) = min(dtx(i),half*deltax(i)*deltax(i)*sph/akk)
159 ENDDO
160 ENDIF
161
162 DO i=1,nel
163 sti(i) = zero
164 IF(off(i)==zero) cycle
165 sti(i) = two*rho(i) * vol(i) / (dtx(i)*dtx(i))
166 dtx(i)= dtfac1(ity)*dtx(i)
167 !dt2 replaced by dt2t
168 IF(nodadt==0)dt2t= min(dtx(i),dt2t)
169 ENDDO
170
171 IF(g_dt/=zero)THEN
172 DO i=1,nel
173 dtsph(i) = dtx(i)
174 ENDDO
175 ENDIF
176
177
178 IF(nodadt==0)THEN
179 IF(idtmin(ity)==1)THEN
180 DO 170 i=1,nel
181 IF(dtx(i)>dtmin1(ity).OR.off(i)==zero) GO TO 170
182 tstop = tt
183#include "lockon.inc"
184 WRITE(iout,*) ' **ERROR : TIME STEP LESS OR EQUAL DTMIN FOR SPH PARTICLE'
185 WRITE(istdo,*)' **ERROR : TIME STEP LESS OR EQUAL DTMIN FOR SPH PARTICLE'
186#include "lockoff.inc"
187 170 CONTINUE
188 ELSEIF(idtmin(ity)==2)THEN
189 DO 270 i=1,nel
190 IF(dtx(i)>dtmin1(ity).OR.off(i)==zero) GO TO 270
191 off(i) = zero
192#include "lockon.inc"
193 WRITE(iout,*) ' -- DELETE SPH PARTICLE',ngl(i)
194 WRITE(istdo,*)' -- DELETE SPH PARTICLE',ngl(i)
195#include "lockoff.inc"
196 270 CONTINUE
197 ELSEIF(idtmin(ity)==5)THEN
198 DO 570 i=1,nel
199 IF(dtx(i)>dtmin1(ity).OR.off(i)==0.) GO TO 570
200 mstop = 2
201#include "lockon.inc"
202 WRITE(iout,*)
203 . ' **ERROR : TIME STEP LESS OR EQUAL DTMIN FOR SPH PARTICLE'
204 WRITE(istdo,*)
205 . ' **ERROR : TIME STEP LESS OR EQUAL DTMIN FOR SPH PARTICLE'
206#include "lockoff.inc"
207 570 CONTINUE
208 ENDIF
209
210 DO i=1,nel
211 IF(dtx(i)>dt2t.OR.off(i)==zero) cycle
212 !nelts and itypts replaced by neltst and itypst
213 neltst =ngl(i)
214 ityptst=ity
215 ENDDO
216
217 ENDIF
218
219 RETURN
220 END
subroutine dtsph(ssp, pm, geo, pid, mat, rho0, vis, deltax, vol, dtx)
Definition dtsph.F:44
subroutine m1lawtot(pm, off, sig, eint, rho, qold, vol, stifn, dt2t, neltst, ityptst, offg, geo, pid, amu, mumax, mat, ngl, ssp, dvol, aire, vnew, vd2, deltax, vis, d1, d2, d3, d4, d5, d6, pnew, psh, qnew, ssp_eq, sold1, sold2, sold3, sold4, sold5, sold6, mssa, dmels, conde, mfxx, mfxy, mfxz, mfyx, mfyy, mfyz, mfzx, mfzy, mfzz, offg0, vol_avg, epsth, dtel, g_dt, nel, etotsh, iselect, ipm, rhoref, rhosp, sigl, ity, ismstr, jtur, jthe, jcvt, jsph, jsms, npg, glob_therm, numgeo, igeo)
Definition m1lawtot.F:63
#define min(a, b)
Definition macros.h:20
#define max(a, b)
Definition macros.h:21
subroutine mdtsph(pm, off, rho, rk, t, re, sti, dt2t, neltst, ityptst, offg, geo, pid, mumax, ssp, vol, vd2, deltax, vis, d1, d2, d3, pnew, psh, mat, ngl, qvis, ssp_eq, g_dt, dtsph, nel, ity, jtur, jthe)
Definition mdtsph.F:47