Télécharger racsti.eso

Retour à la liste

Numérotation des lignes :

  1. C RACSTI SOURCE BECC 11/05/26 21:16:16 6981
  2. SUBROUTINE RACSTI(
  3. & PC_L, GAM_L,
  4. & PC_R, GAM_R,
  5. & PMIN,
  6. & RHO_L, P_L, U_L,
  7. & RHO_R, P_R, U_R,
  8. & D_L, D_B, D_A, D_R,
  9. & RHO_B, P_B, U_B,
  10. & RHO_A, P_A, U_A,
  11. & LOGDEB, LOGAN, LOGNC, LOGVAC)
  12. *
  13. *************************************************************************
  14. *
  15. * project :
  16. *
  17. * name : racsti
  18. *
  19. * description : euler equations for a mixture of stiffened gases
  20. * flux in the non-reactive case.
  21. *
  22. * field-by-field decomposition (entropy-respecting).
  23. *
  24. * language : fortran 77
  25. *
  26. * author : a. beccantini den/dm2s/sfme/ltmf
  27. *
  28. *************************************************************************
  29. *
  30. * called by :
  31. *
  32. *
  33. *************************************************************************
  34. *
  35. ***** input
  36. *
  37. * pc_l, gam_l = properties of the gas in the left
  38. *
  39. * pc_r, gam_r = properties of the gas in the right
  40. *
  41. * rho_lr, p_lr, w_lr
  42. * = density, pressure, velocity
  43. *
  44. * pmin = if p < pmin => vacuum
  45. *
  46. * logdeb = debugging ?
  47. *
  48. ***** input / output
  49. *
  50. * logvac = vacuum ?
  51. *
  52. ***** output
  53. *
  54. * rho_ba, p_ba, w_ba
  55. * = density, pressure, velocity
  56. *
  57. * d_l, d_b, d_a, d_r
  58. * = wave speeds d_l <= d_b <= d_a <= d_r
  59. *
  60. * logan = anomaly
  61. *
  62. * lognc = if true, the Newton for the Riemann problem
  63. * solution did not converge
  64. *
  65. *************************************************************************
  66. *
  67. * 26/11/2009 created
  68. * 25/05/2011 evolution in CAST3M
  69. *
  70. *************************************************************************
  71. *
  72. * n.b.: all variables are declared
  73. *
  74. C IMPLICIT NONE
  75. IMPLICIT INTEGER(I-N)
  76. INTEGER ITER, NITER
  77. PARAMETER (NITER=50)
  78. REAL*8
  79. & PMIN,
  80. & PC_L, GAM_L,
  81. & PC_R, GAM_R,
  82. & D_L, D_B, D_A, D_R,
  83. & RHO_L, P_L, U_L,
  84. & RHO_B, P_B, U_B,
  85. & RHO_A, P_A, U_A,
  86. & RHO_R, P_R, U_R
  87. & , EPSERR
  88. PARAMETER (EPSERR = 1.0D-6)
  89. *
  90. REAL*8 P_VAC, A_L, A_R
  91. & , P_ES, ROES_B, UES_B, DES_B, DER_B
  92. & , ROES_A, UES_A, DES_A, DER_A
  93. & , P_INT
  94. *
  95. * debugging ?
  96. *
  97. LOGICAL LOGDEB, LOGAN, LOGNC, LOGVAC
  98. *
  99. LOGNC = .TRUE.
  100. *
  101. ***** initialisation of the speed of the wave on the left
  102. * and on the right
  103. *
  104. CALL ASTIFF(RHO_L, P_L, GAM_L, PC_L, A_L, LOGDEB, LOGAN)
  105. CALL ASTIFF(RHO_R, P_R, GAM_R, PC_R, A_R, LOGDEB, LOGAN)
  106. IF ( LOGAN ) THEN
  107. WRITE(*,*) 'SUBROUTINE RACSTI.F'
  108. WRITE(*,*) 'INITIALISATION OF THE WAVE SPEEDS'
  109. WRITE(*,*) 'ANOMALY DETECTED'
  110. GOTO 9999
  111. ENDIF
  112. D_L = U_L - A_L
  113. D_R = U_R + A_R
  114. *
  115. ******** states
  116. *
  117. * l
  118. * ***
  119. * |* r
  120. * | * b *********
  121. * | ******** *|
  122. * | | | * a * |
  123. * | | | ************* |
  124. * | | | | | |
  125. * | | | | | |
  126. * d_l d_b un_b un_a d_a d_r
  127. *
  128. *
  129. * if no vacuum, un_b = un_a
  130. *
  131. P_VAC = MIN( P_L, P_R )
  132. P_VAC = MAX( P_VAC, PMIN)
  133. P_VAC = 1.0D-8 * P_VAC
  134. CALL FSTRAL( PC_L, GAM_L, P_VAC, RHO_L, P_L, U_L,
  135. & U_B, RHO_B, D_B, LOGDEB, LOGAN)
  136. CALL FSTRAR(PC_R, GAM_R, P_VAC, RHO_R, P_R, U_R,
  137. & U_A, RHO_A, D_A, LOGDEB, LOGAN)
  138. IF ( LOGAN ) THEN
  139. WRITE(*,*) 'SUBROUTINE RACSTI.F'
  140. WRITE(*,*) 'COMPUTATION OF THE VACUUM STATE'
  141. WRITE(*,*) 'ANOMALY DETECTED'
  142. GOTO 9999
  143. ENDIF
  144. *
  145. IF (LOGVAC) THEN
  146. *
  147. * From elsewhere, we know that we are in the vacuum...
  148. *
  149. P_B = P_VAC
  150. P_A = P_VAC
  151. LOGNC = .FALSE.
  152.  
  153. C
  154. ELSEIF ( U_B .LT. U_A) THEN
  155. *
  156. * vacuum
  157. *
  158. P_B = P_VAC
  159. P_A = P_VAC
  160. LOGVAC = .TRUE.
  161. LOGNC = .FALSE.
  162. ELSE
  163. LOGVAC = .FALSE.
  164. *
  165. *************************************************************************
  166. ***** intersection (states l, b, a, r) **********************************
  167. *************************************************************************
  168. *
  169. *
  170. ***** initialization of p_int and of the states b, a
  171. *
  172. *
  173. P_INT = 0.5D0 * (P_R + P_L)
  174. *
  175. CALL FSTERL( PC_L, GAM_L, P_INT, RHO_L, P_L, U_L,
  176. & U_B, RHO_B, D_B, LOGDEB, LOGAN)
  177. CALL FSTERR( PC_R, GAM_R, P_INT, RHO_R, P_R, U_R,
  178. & U_A, RHO_A, D_A, LOGDEB, LOGAN)
  179. *
  180. IF ( LOGAN ) THEN
  181. WRITE(*,*) 'SUBROUTINE RACSTI.F'
  182. WRITE(*,*) 'INITIALIZATION OF THE STATES'
  183. WRITE(*,*) 'ANOMALY DETECTED'
  184. GOTO 9999
  185. ENDIF
  186. *
  187. DO ITER = 1, NITER, 1
  188. *
  189. * pmin
  190. *
  191. P_VAC = 1D-8 * MIN (P_L, P_R, P_INT)
  192. P_VAC = MAX( P_VAC, PMIN)
  193. *
  194. IF (.TRUE.) THEN
  195. *
  196. * exact evaluation of the derivative
  197. *
  198. CALL DERL(PC_L, GAM_L, P_INT, RHO_L, P_L, U_L, DER_B,
  199. $ LOGDEB, LOGAN)
  200. IF ( LOGAN ) THEN
  201. WRITE(*,*) 'SUBROUTINE RACSTI.F'
  202. WRITE(*,*) 'EVALUATION OF THE DERIVATES'
  203. WRITE(*,*) 'ANOMALY DETECTED'
  204. GOTO 9999
  205. ENDIF
  206. CALL DERR(PC_R, GAM_R, P_INT, RHO_R, P_R, U_R, DER_A,
  207. $ LOGDEB, LOGAN)
  208. *
  209. C WRITE(*,*) 'DER_B, DER_A ex = ', DER_B, DER_A
  210. C
  211. IF ( LOGAN ) THEN
  212. WRITE(*,*) 'SUBROUTINE RACSTI.F'
  213. WRITE(*,*) 'EVALUATION OF THE DERIVATES'
  214. WRITE(*,*) 'ANOMALY DETECTED'
  215. GOTO 9999
  216. ENDIF
  217. ELSE
  218. *
  219. * numerical evaluation of the derivatives
  220. *
  221. P_ES = 1.001D0 * P_INT
  222. CALL FSTERL( PC_L, GAM_L, P_ES, RHO_L, P_L, U_L,
  223. & UES_B, ROES_B, DES_B, LOGDEB, LOGAN)
  224. IF ( LOGAN ) THEN
  225. WRITE(*,*) 'SUBROUTINE RACSTI.F'
  226. WRITE(*,*) 'EVALUATION OF THE DERIVATES'
  227. WRITE(*,*) 'ANOMALY DETECTED'
  228. GOTO 9999
  229. ENDIF
  230. DER_B = (UES_B - U_B) / (P_ES - P_INT)
  231. CALL FSTERR( PC_R, GAM_R, P_ES, RHO_R, P_R, U_R,
  232. & UES_A, ROES_A, DES_A, LOGDEB, LOGAN)
  233. IF ( LOGAN ) THEN
  234. WRITE(*,*) 'SUBROUTINE RACSTI.F'
  235. WRITE(*,*) 'EVALUATION OF THE DERIVATES'
  236. WRITE(*,*) 'ANOMALY DETECTED'
  237. GOTO 9999
  238. ENDIF
  239. DER_A = (UES_A - U_A) / (P_ES - P_INT)
  240. C WRITE(*,*) 'DER_B, DER_A num = ', DER_B, DER_A
  241. C STOP
  242. ENDIF
  243. IF ((DER_B .GT. 0.0D0)) THEN
  244. WRITE(*,*) 'SUBROUTINE RACSTI.F'
  245. WRITE (*,*) 'DER_B =', DER_B
  246. LOGAN = .TRUE.
  247. GOTO 9999
  248. ENDIF
  249. IF ((DER_A .LT. 0.0D0)) THEN
  250. WRITE(*,*) 'SUBROUTINE RACSTI.F'
  251. WRITE (*,*) 'DER_A =', DER_A
  252. LOGAN = .TRUE.
  253. GOTO 9999
  254. ENDIF
  255. IF ((DER_B - DER_A) .GE. 0.0D0) THEN
  256. WRITE(*,*) 'SUBROUTINE RACSTI.F'
  257. WRITE (*,*) 'DER_B - DER_A =', (DER_B - DER_A)
  258. LOGAN = .TRUE.
  259. GOTO 9999
  260. ENDIF
  261. C
  262. P_ES = P_INT - (U_B - U_A) / (DER_B - DER_A)
  263. P_ES = MAX (P_ES, P_VAC)
  264. CALL FSTERL( PC_L, GAM_L, P_ES, RHO_L, P_L, U_L,
  265. & U_B, RHO_B, D_B, LOGDEB, LOGAN)
  266. CALL FSTERR( PC_R, GAM_R, P_ES, RHO_R, P_R, U_R,
  267. & U_A, RHO_A, D_A, LOGDEB, LOGAN)
  268. IF ( LOGAN ) THEN
  269. WRITE(*,*) 'SUBROUTINE RACSTI.F'
  270. WRITE(*,*) 'COMPUTATION OF THE INTERSECTION'
  271. WRITE(*,*) 'ANOMALY DETECTED'
  272. GOTO 9999
  273. ENDIF
  274. IF ((ABS (P_ES - P_INT) / P_INT) .LT. EPSERR) THEN
  275. LOGNC = .FALSE.
  276. P_A = P_ES
  277. P_B = P_ES
  278. IF (P_B .GT. P_L) THEN
  279. D_L = D_B
  280. ENDIF
  281. IF (P_A .GT. P_R) THEN
  282. D_R = D_A
  283. ENDIF
  284. IF (LOGDEB) THEN
  285. WRITE(*,*) 'SUBROUTINE RACSTI.F'
  286. WRITE(*,*) 'CONVERGENCE ACHIEVED'
  287. WRITE(*,*) 'ITER =', ITER
  288. ENDIF
  289. GOTO 9999
  290. ENDIF
  291. P_INT = P_ES
  292. IF (LOGDEB) THEN
  293. WRITE(*,*) 'SUBROUTINE RACSTI'
  294. WRITE(*,*) 'P_INT = ', P_INT
  295. WRITE(*,*) 'RHO_B, U_B', RHO_B, U_B
  296. WRITE(*,*) 'RHO_A, U_A', RHO_A, U_A
  297. ENDIF
  298. ENDDO
  299. P_A = P_ES
  300. P_B = P_ES
  301. IF (LOGNC) THEN
  302. WRITE(*,*) 'SUBROUTINE RACSTI.F'
  303. WRITE(*,*) 'CONVERGENCE NOT ACHIEVED'
  304. WRITE(*,*) 'ITER =', ITER
  305. WRITE(*,*) 'ERROR =', (ABS (P_ES - P_INT) / P_INT)
  306. ENDIF
  307. ENDIF
  308. *
  309. 9999 RETURN
  310. END
  311.  
  312.  

© Cast3M 2003 - Tous droits réservés.
Mentions légales