OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
tabulated_mod Module Reference

Functions/Subroutines

subroutine tabulated (iflag, nel, pmin, off, eint, mu, espe, dvol, df, vnew, psh, pnew, dpdm, dpde, npf, tf, snpc, stf, eos_struct)

Function/Subroutine Documentation

◆ tabulated()

subroutine tabulated_mod::tabulated ( integer, intent(in) iflag,
integer, intent(in) nel,
intent(in) pmin,
dimension(nel), intent(in) off,
dimension(nel), intent(inout) eint,
dimension(nel), intent(in) mu,
dimension(nel), intent(in) espe,
dimension(nel), intent(in) dvol,
dimension(nel), intent(in) df,
dimension(nel), intent(in) vnew,
dimension(nel), intent(inout) psh,
dimension(nel), intent(inout) pnew,
dimension(nel), intent(inout) dpdm,
dimension(nel), intent(inout) dpde,
integer, dimension(snpc), intent(in) npf,
dimension(stf), intent(in) tf,
integer, intent(in) snpc,
integer, intent(in) stf,
type(eos_param_ ), intent(in) eos_struct )

Definition at line 41 of file tabulated.F.

46C-----------------------------------------------
47C M o d u l e s
48C-----------------------------------------------
49 USE constant_mod , ONLY : zero, em15, half, one, two
50 USE eos_param_mod , ONLY : eos_param_
51C-----------------------------------------------
52C D e s c r i p t i o n
53C-----------------------------------------------
54C This subroutine contains numerical solving of TABULATED Equation Of State (/EOS/TABULATED)
55C P(mu,E) = A(mu) + B(mu)*E
56C where A and B are two user functions (/FUNCT)
57C
58!----------------------------------------------------------------------------
59!! \details STAGGERED SCHEME IS EXECUTED IN TWO PASSES IN EOSMAIN : IFLG=0 THEN IFLG=1
60!! \details COLLOCATED SCHEME IS DOING A SINGLE PASS : IFLG=2
61!! \details
62!! \details STAGGERED SCHEME
63!! \details EOSMAIN / IFLG = 0 : DERIVATIVE CALCULATION FOR SOUND SPEED ESTIMATION c[n+1] REQUIRED FOR PSEUDO-VISCOSITY (DPDE:partial derivative, DPDM:total derivative)
64!! \details MQVISCB : PSEUDO-VISCOSITY Q[n+1]
65!! \details MEINT : INTERNAL ENERGY INTEGRATION FOR E[n+1] : FIRST PART USING P[n], Q[n], and Q[n+1] CONTRIBUTIONS
66!! \details EOSMAIN / IFLG = 1 : UPDATE P[n+1], T[N+1]
67!! \details INTERNAL ENERGY INTEGRATION FOR E[n+1] : LAST PART USING P[n+1] CONTRIBUTION
68!! \details (second order integration dE = -P.dV where P = 0.5(P[n+1] + P[n]) )
69!! \details COLLOCATED SCHEME
70!! \details EOSMAIN / IFLG = 2 : SINGLE PASS FOR P[n+1] AND DERIVATIVES
71!----------------------------------------------------------------------------
72C-----------------------------------------------
73C I m p l i c i t T y p e s
74C-----------------------------------------------
75 IMPLICIT NONE
76#include "my_real.inc"
77C-----------------------------------------------
78C D u m m y A r g u m e n t s
79C-----------------------------------------------
80 INTEGER, INTENT(IN) :: SNPC,STF
81 INTEGER, INTENT(IN) :: IFLAG, NEL,NPF(SNPC)
82 my_real, INTENT(IN) :: off(nel), mu(nel), espe(nel), dvol(nel), df(nel), vnew(nel), tf(stf)
83 my_real,INTENT(INOUT) :: psh(nel), eint(nel), dpdm(nel), dpde(nel), pnew(nel)
84 my_real,INTENT(IN) :: pmin
85 TYPE(EOS_PARAM_ ),INTENT(IN) :: EOS_STRUCT
86C-----------------------------------------------
87C L o c a l V a r i a b l e s
88C-----------------------------------------------
89 INTEGER I
90 my_real :: e0,aa,bb,dvv,pp
91 my_real :: xscale_a,xscale_b,fscale_a,fscale_b
92 INTEGER :: A_fun_id, B_fun_id
93 my_real :: res_a(nel),res_b(nel),deri_a(nel),deri_b(nel)
94 my_real,EXTERNAL :: finter
95C-----------------------------------------------
96C S o u r c e L i n e s
97C-----------------------------------------------
98 e0 = eos_struct%E0
99 psh(1:nel) = eos_struct%PSH
100 xscale_a = eos_struct%UPARAM(1)
101 xscale_b = eos_struct%UPARAM(2)
102 fscale_a = eos_struct%UPARAM(3)
103 fscale_b = eos_struct%UPARAM(4)
104 a_fun_id = eos_struct%IPARAM(1)
105 b_fun_id = eos_struct%IPARAM(2)
106
107 IF(iflag == 0) THEN
108 ! both A_fun_id & B_fun_id cannot be 0. This is ensured by the reader during the Starter process
109 IF(a_fun_id == 0)THEN
110 DO i=1,nel
111 res_a(i) = zero
112 deri_a(i) = zero
113 res_b(i) = fscale_b*finter(b_fun_id,mu(i),npf,tf,deri_b(i))
114 ENDDO
115 ELSEIF(b_fun_id == 0)THEN
116 DO i=1,nel
117 res_a(i) = fscale_a*finter(a_fun_id,mu(i),npf,tf,deri_a(i))
118 res_b(i) = zero
119 deri_b(i) = zero
120 ENDDO
121 ELSE
122 DO i=1,nel
123 res_a(i) = fscale_a*finter(a_fun_id,mu(i),npf,tf,deri_a(i))
124 res_b(i) = fscale_b*finter(b_fun_id,mu(i),npf,tf,deri_b(i))
125 ENDDO
126 ENDIF
127 DO i=1,nel
128 pp = res_a(i) + res_b(i) * espe(i) - psh(i) ! A(MU(I))+B(MU(I))*ESPE(I)
129 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
130 dpde(i) = res_b(i) ! B(MU(I)) !partial derivative
131 pnew(i) = max(pp,pmin)*off(i) ! P(mu[n+1],E[n])
132 ENDDO
133
134C-----------------------------------------------
135 ELSEIF(iflag == 1) THEN
136 ! both A_fun_id & B_fun_id cannot be 0. This is ensured by the reader during the Starter process
137 IF(a_fun_id == 0)THEN
138 DO i=1,nel
139 res_a(i) = zero
140 res_b(i) = fscale_b*finter(b_fun_id,mu(i),npf,tf,deri_b(i))
141 ENDDO
142 ELSEIF(b_fun_id == 0)THEN
143 DO i=1,nel
144 res_a(i) = fscale_a*finter(a_fun_id,mu(i),npf,tf,deri_a(i))
145 res_b(i) = zero
146 ENDDO
147 ELSE
148 DO i=1,nel
149 res_a(i) = fscale_a*finter(a_fun_id,mu(i),npf,tf,deri_a(i))
150 res_b(i) = fscale_b*finter(b_fun_id,mu(i),npf,tf,deri_b(i))
151 ENDDO
152 ENDIF
153 DO i=1,nel
154 aa = res_a(i)
155 bb = res_b(i)
156 dpde(i) = bb
157 dvv = half*dvol(i)*df(i) / max(em15,vnew(i))
158 pp = aa + bb * espe(i)
159 pnew(i) = (aa+bb*(espe(i)-psh(i)*dvv))/(one+bb*dvv)
160 pnew(i) = max(pnew(i),pmin )*off(i) ! P(mu[n+1],E[n+1])
161 eint(i) = eint(i) - half*dvol(i)*(pnew(i)+psh(i) )
162 ENDDO
163
164C-----------------------------------------------
165 ELSEIF (iflag == 2) THEN
166 ! both A_fun_id & B_fun_id cannot be 0. This is ensured by the reader during the Starter process
167 IF(a_fun_id == 0)THEN
168 DO i=1,nel
169 res_a(i) = zero
170 res_b(i) = fscale_b*finter(b_fun_id,mu(i),npf,tf,deri_b(i))
171 deri_a(i) = zero
172 ENDDO
173 ELSEIF(b_fun_id == 0)THEN
174 DO i=1,nel
175 res_a(i) = fscale_a*finter(a_fun_id,mu(i),npf,tf,deri_a(i))
176 res_b(i) = zero
177 deri_b(i) = zero
178 ENDDO
179 ELSE
180 DO i=1,nel
181 res_a(i) = fscale_a*finter(a_fun_id,mu(i),npf,tf,deri_a(i))
182 res_b(i) = fscale_b*finter(b_fun_id,mu(i),npf,tf,deri_b(i))
183 ENDDO
184 ENDIF
185 DO i=1, nel
186 IF (vnew(i) > zero) THEN
187 pp = res_a(i) + res_b(i)*espe(i) - psh(i)
188 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
189 dpde(i) = res_b(i) ! B(MU(I)) !partial derivative
190 pnew(i) = pp
191 ENDIF
192 ENDDO
193
194 ENDIF
195C-----------------------------------------------
196 RETURN
my_real
#define my_real
Definition cppsort.cpp:32
max
#define max(a, b)
Definition macros.h:21
constant_mod
Definition constant_mod.F:355