OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
tabulated.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!|| tabulated ../common_source/eos/tabulated.F
25!||--- called by ------------------------------------------------------
26!|| eosmain ../common_source/eos/eosmain.F
27!||====================================================================
28 SUBROUTINE tabulated(IFLAG,NEL ,PM ,OFF ,EINT ,MU ,
29 2 ESPE ,DVOL ,DF ,VNEW ,MAT ,PSH ,
30 3 PNEW ,DPDM ,DPDE ,
31 4 NPF ,TF)
32C-----------------------------------------------
33C M o d u l e s
34C-----------------------------------------------
35C
36C-----------------------------------------------
37C D e s c r i p t i o n
38C-----------------------------------------------
39C This subroutine contains numerical solving of TABULATED Equation Of State (/EOS/TABULATED)
40C P(mu,E) = A(mu) + B(mu)*E
41C where A and B are two user functions (/FUNCT)
42C
43!----------------------------------------------------------------------------
44!! \details STAGGERED SCHEME IS EXECUTED IN TWO PASSES IN EOSMAIN : IFLG=0 THEN IFLG=1
45!! \details COLLOCATED SCHEME IS DOING A SINGLE PASS : IFLG=2
46!! \details
47!! \details STAGGERED SCHEME
48!! \details EOSMAIN / IFLG = 0 : DERIVATIVE CALCULATION FOR SOUND SPEED ESTIMATION c[n+1] REQUIRED FOR PSEUDO-VISCOSITY (DPDE:partial derivative, DPDM:total derivative)
49!! \details MQVISCB : PSEUDO-VISCOSITY Q[n+1]
50!! \details MEINT : INTERNAL ENERGY INTEGRATION FOR E[n+1] : FIRST PART USING P[n], Q[n], and Q[n+1] CONTRIBUTIONS
51!! \details EOSMAIN / IFLG = 1 : UPDATE P[n+1], T[N+1]
52!! \details INTERNAL ENERGY INTEGRATION FOR E[n+1] : LAST PART USING P[n+1] CONTRIBUTION
53!! \details (second order integration dE = -P.dV where P = 0.5(P[n+1] + P[n]) )
54!! \details COLLOCATED SCHEME
55!! \details EOSMAIN / IFLG = 2 : SINGLE PASS FOR P[n+1] AND DERIVATIVES
56!----------------------------------------------------------------------------
57C-----------------------------------------------
58C I m p l i c i t T y p e s
59C-----------------------------------------------
60#include "implicit_f.inc"
61#include "comlock.inc"
62C-----------------------------------------------
63C C o m m o n B l o c k s
64C-----------------------------------------------
65#include "param_c.inc"
66#include "tabsiz_c.inc"
67#include "com04_c.inc"
68#include "com06_c.inc"
69#include "com08_c.inc"
70#include "vect01_c.inc"
71#include "scr06_c.inc"
72C-----------------------------------------------
73C D u m m y A r g u m e n t s
74C-----------------------------------------------
75 INTEGER, INTENT(IN) :: MAT(NEL), IFLAG, NEL,NPF(SNPC)
76 my_real, INTENT(INOUT) :: PM(NPROPM,NUMMAT),
77 . OFF(NEL) ,EINT(NEL) ,MU(NEL) ,
78 . espe(nel) ,dvol(nel) ,df(nel) ,
79 . vnew(nel) ,pnew(nel) ,dpdm(nel),
80 . dpde(nel) ,
81 . psh(nel), tf(stf)
82C-----------------------------------------------
83C L o c a l V a r i a b l e s
84C-----------------------------------------------
85 INTEGER I, MX
86 my_real :: E0,AA,BB,DVV,PP
87 my_real :: XSCALE_A,XSCALE_B,FSCALE_A,FSCALE_B
88 INTEGER :: A_fun_id, B_fun_id
89 my_real :: res_a(nel),res_b(nel),deri_a(nel),deri_b(nel),pc
90 my_real,EXTERNAL :: finter
91C-----------------------------------------------
92C S o u r c e L i n e s
93C-----------------------------------------------
94 IF(iflag == 0) THEN
95 mx = mat(1)
96 e0 = pm(23,mx)
97 psh(1:nel) = pm(88,mx)
98 xscale_a = pm(33,mx)
99 xscale_b = pm(34,mx)
100 a_fun_id = pm(35,mx)
101 b_fun_id = pm(36,mx)
102 pc = pm( 37,mx)
103 fscale_a = pm(160,mx)
104 fscale_b = pm(161,mx)
105
106 ! both A_fun_id & B_fun_id cannot be 0. This is ensured by the reader during the Starter process
107 IF(a_fun_id == 0)THEN
108 DO i=1,nel
109 res_a(i) = zero
110 deri_a(i) = zero
111 res_b(i) = fscale_b*finter(b_fun_id,mu(i),npf,tf,deri_b(i))
112 ENDDO
113 ELSEIF(b_fun_id == 0)THEN
114 DO i=1,nel
115 res_a(i) = fscale_a*finter(a_fun_id,mu(i),npf,tf,deri_a(i))
116 res_b(i) = zero
117 deri_b(i) = zero
118 ENDDO
119 ELSE
120 DO i=1,nel
121 res_a(i) = fscale_a*finter(a_fun_id,mu(i),npf,tf,deri_a(i))
122 res_b(i) = fscale_b*finter(b_fun_id,mu(i),npf,tf,deri_b(i))
123 ENDDO
124 ENDIF
125
126 DO i=1,nel
127 pp = res_a(i) + res_b(i) * espe(i) - psh(i) ! A(MU(I))+B(MU(I))*ESPE(I)
128 dpdm(i) = deri_a(i)+deri_b(i)*espe(i) + res_b(i)*(pp+psh(i))/( (one+mu(i))*(one+mu(i)) ) ! A'(MU0) + B'(MU0)*E0+B(MU0)/(ONE+MU0)/(ONE+MU0)*P0 !total derivative
129 dpde(i) = res_b(i) ! B(MU(I)) !partial derivative
130 pnew(i) = max(pp,pc)*off(i) ! P(mu[n+1],E[n])
131 ENDDO
132
133C-----------------------------------------------
134 ELSEIF(iflag == 1) THEN
135 mx = mat(1)
136 e0 = pm(23,mx)
137 psh(1:nel) = pm(88,mx)
138 xscale_a = pm(33,mx)
139 xscale_b = pm(34,mx)
140 a_fun_id = pm(35,mx)
141 b_fun_id = pm(36,mx)
142 pc = pm( 37,mx)
143 fscale_a = pm(160,mx)
144 fscale_b = pm(161,mx)
145 ! both A_fun_id & B_fun_id cannot be 0. This is ensured by the reader during the Starter process
146 IF(a_fun_id == 0)THEN
147 DO i=1,nel
148 res_a(i) = zero
149 res_b(i) = fscale_b*finter(b_fun_id,mu(i),npf,tf,deri_b(i))
150 ENDDO
151 ELSEIF(b_fun_id == 0)THEN
152 DO i=1,nel
153 res_a(i) = fscale_a*finter(a_fun_id,mu(i),npf,tf,deri_a(i))
154 res_b(i) = zero
155 ENDDO
156 ELSE
157 DO i=1,nel
158 res_a(i) = fscale_a*finter(a_fun_id,mu(i),npf,tf,deri_a(i))
159 res_b(i) = fscale_b*finter(b_fun_id,mu(i),npf,tf,deri_b(i))
160 ENDDO
161 ENDIF
162 DO i=1,nel
163 aa = res_a(i)
164 bb = res_b(i)
165 dpde(i) = bb
166 dvv = half*dvol(i)*df(i) / max(em15,vnew(i))
167 pp = aa + bb * espe(i)
168 pnew(i) = (aa+bb*(espe(i)-psh(i)*dvv))/(one+bb*dvv)
169 pnew(i) = max(pnew(i),pc )*off(i) ! P(mu[n+1],E[n+1])
170 eint(i) = eint(i) - half*dvol(i)*(pnew(i)+psh(i) )
171 ENDDO
172
173C-----------------------------------------------
174 ELSEIF (iflag == 2) THEN
175 mx = mat(1)
176 e0 = pm(23,mx)
177 psh(1:nel) = pm(88,mx)
178 xscale_a = pm(33,mx)
179 xscale_b = pm(34,mx)
180 a_fun_id = pm(35,mx)
181 b_fun_id = pm(36,mx)
182 pc = pm( 37,mx)
183 fscale_a = pm(160,mx)
184 fscale_b = pm(161,mx)
185 ! both A_fun_id & B_fun_id cannot be 0. This is ensured by the reader during the Starter process
186 IF(a_fun_id == 0)THEN
187 DO i=1,nel
188 res_a(i) = zero
189 res_b(i) = fscale_b*finter(b_fun_id,mu(i),npf,tf,deri_b(i))
190 deri_a(i) = zero
191 ENDDO
192 ELSEIF(b_fun_id == 0)THEN
193 DO i=1,nel
194 res_a(i) = fscale_a*finter(a_fun_id,mu(i),npf,tf,deri_a(i))
195 res_b(i) = zero
196 deri_b(i) = zero
197 ENDDO
198 ELSE
199 DO i=1,nel
200 res_a(i) = fscale_a*finter(a_fun_id,mu(i),npf,tf,deri_a(i))
201 res_b(i) = fscale_b*finter(b_fun_id,mu(i),npf,tf,deri_b(i))
202 ENDDO
203 ENDIF
204 DO i=1, nel
205 IF (vnew(i) > zero) THEN
206 pp = res_a(i) + res_b(i)*espe(i) - psh(i)
207 dpdm(i) = deri_a(i)+deri_b(i)*espe(i) + res_b(i)*(pp+psh(i))/( (one+mu(i))*(one+mu(i)) ) ! A'(MU0) + B'(MU0)*E0+B(MU0)/(ONE+MU0)/(ONE+MU0)*P0 !total derivative
208 dpde(i) = res_b(i) ! B(MU(I)) !partial derivative
209 pnew(i) = pp
210 ENDIF
211 ENDDO
212
213 ENDIF
214C-----------------------------------------------
215 RETURN
216 END
max
#define max(a, b)
Definition macros.h:21
tabulated
subroutine tabulated(iflag, nel, pm, off, eint, mu, espe, dvol, df, vnew, mat, psh, pnew, dpdm, dpde, npf, tf)
Definition tabulated.F:32