hydrocode_Radial_Lag 0.3
This is an implementation of fully explict forward Euler scheme for multi-D radially symmetric compressible flows on Lagrangian coordinate
riemann_solver.h 文件参考

This file is the header file of several Riemann solvers and GRP solvers. 更多...

riemann_solver.h 的引用(Include)关系图:
此图展示该文件直接或间接的被哪些文件引用了:

浏览源代码.

宏定义

#define Riemann_solver_exact_single   Riemann_solver_exact_Ben
 Which solver is chosen as the exact Riemann solver for single-component flow. 更多...
 

函数

double Riemann_solver_exact (double *U_star, double *P_star, const double gammaL, const double gammaR, const double u_L, const double u_R, const double p_L, const double p_R, const double c_L, const double c_R, _Bool *CRW, const double eps, const double tol, int N)
 
double Riemann_solver_starPU (double *U_star, double *P_star, const double GammaL, const double GammaR, const double UL, const double UR, const double PL, const double PR, const double CL, const double CR, _Bool *CRW, const double eps, const double TOLPRE, const int NRITER)
 EXACT RIEMANN SOLVER FOR Two-Component γ-Law Gas 更多...
 
double Riemann_solver_exact_Ben (double *U_star, double *P_star, const double gamma, const double u_L, const double u_R, const double p_L, const double p_R, const double c_L, const double c_R, _Bool *CRW, const double eps, const double tol, const int N)
 
double Riemann_solver_exact_Toro (double *U_star, double *P_star, const double gamma, const double U_l, const double U_r, const double P_l, const double P_r, const double c_l, const double c_r, _Bool *CRW, const double eps, const double tol, const int N)
 
void linear_GRP_solver_LAG (double *D, double *U, const struct i_f_var *ifv_L, const struct i_f_var *ifv_R, const double eps, const double atc)
 
void linear_GRP_solver_Edir (double *D, double *U, const struct i_f_var *ifv_L, const struct i_f_var *ifv_R, const double eps, const double atc)
 
void linear_GRP_solver_Edir_Q1D (double *wave_speed, double *D, double *U, double *U_star, const struct i_f_var *ifv_L, const struct i_f_var *ifv_R, const double eps, const double atc)
 
void linear_GRP_solver_Edir_G2D (double *wave_speed, double *D, double *U, double *U_star, const struct i_f_var *ifv_L, const struct i_f_var *ifv_R, const double eps, const double atc)
 
void AcousticRLagTangent (double *dire, double *U_star, const struct i_f_var *ifv_L, const struct i_f_var *ifv_R, const double r, const double M, const double eps)
 A GRP solver for unsteady compressible inviscid two-component flow in tangential case. Lagrangian version (moving mesh) for cylindrical case. 更多...
 
void GRPsolverRLag (double *wave_speed, double *dire, double *U_star, const struct i_f_var *ifv_L, const struct i_f_var *ifv_R, const double r, const double M, const double eps, const double atc)
 A GRP solver for unsteady compressible inviscid two-component flow. Lagrangian version (moving mesh) cylindrical case. 更多...
 
void HLL_2D_solver (double *F, double *lambda_max, const struct i_f_var *ifv_L, const struct i_f_var *ifv_R)
 
void Roe_solver (double *F, double *lambda_max, const struct i_f_var *ifv_L, const struct i_f_var *ifv_R, const double delta)
 
void Roe_2D_solver (double *F, double *lambda_max, const struct i_f_var *ifv_L, const struct i_f_var *ifv_R, const double delta)
 
void Roe_HLL_solver (double *V_mk, double *F, double *lambda_max, const struct i_f_var *ifv_L, const struct i_f_var *ifv_R, const double delta)
 

详细描述

This file is the header file of several Riemann solvers and GRP solvers.

This header file declares functions in the folder 'Riemann_solver'.

在文件 riemann_solver.h 中定义.

宏定义说明

◆ Riemann_solver_exact_single

#define Riemann_solver_exact_single   Riemann_solver_exact_Ben

Which solver is chosen as the exact Riemann solver for single-component flow.

在文件 riemann_solver.h17 行定义.

函数说明

◆ AcousticRLagTangent()

void AcousticRLagTangent ( double *  dire,
double *  U_star,
const struct i_f_var ifv_L,
const struct i_f_var ifv_R,
const double  r,
const double  M,
const double  eps 
)

A GRP solver for unsteady compressible inviscid two-component flow in tangential case. Lagrangian version (moving mesh) for cylindrical case.

参数
[out]direthe temporal derivative of fluid variables in the Star Region.
[*, u, p]_t
[out]U_starthe Riemann solutions in the Star Region.
[*, u_star, p_star]
[in]ifv_LLeft States (rho_L, u_L, p_L, d_u_L, d_p_L, gammaL).
[in]ifv_RRight States (rho_R, u_R, p_R, d_u_R, d_p_R, gammaR).
  • d_u, d_p: spatial derivatives in coordinate r.
  • gamma: the constant of the perfect gas.
[in]rthe r-coordinate value.
[in]MSpatial dimension number for radially symmetric flow.
  • M=1: planar flow.
  • M=2: cylindrical flow. (√)
  • M=3: spherical flow.
[in]epsthe largest value could be seen as zero.

在文件 linear_grp_solver_radial_LAG.c30 行定义.

函数调用图:
这是这个函数的调用关系图:

◆ GRPsolverRLag()

void GRPsolverRLag ( double *  wave_speed,
double *  dire,
double *  U_star,
const struct i_f_var ifv_L,
const struct i_f_var ifv_R,
const double  r,
const double  M,
const double  eps,
const double  atc 
)

A GRP solver for unsteady compressible inviscid two-component flow. Lagrangian version (moving mesh) cylindrical case.

参数
[out]wave_speedthe velocity of left and right waves.
[out]direthe temporal derivative of fluid variables in the Star Region.
[rho_L, u, p, rho_R]_t
[out]U_starthe Riemann solutions in the Star Region.
[rho_star_L, u_star, p_star, rho_star_R]
[in]ifv_LLeft States (rho_L, u_L, p_L, d_rho_L, d_u_L, d_p_L, gammaL).
[in]ifv_RRight States (rho_R, u_R, p_R, d_rho_R, d_u_R, d_p_R, gammaR).
  • d_rho, d_u, d_p: r-spatial derivatives.
  • gamma: the constant of the perfect gas.
[in]rthe r-coordinate value.
[in]MSpatial dimension number for radially symmetric flow.
  • M=1: planar flow.
  • M=2: cylindrical flow. (√)
  • M=3: spherical flow.
[in]epsthe largest value could be seen as zero.
[in]atcParameter that determines the solver type.
  • INFINITY: acoustic approximation
  • eps: GRP solver(nonlinear + acoustic case)
  • -0.0: GRP solver(only nonlinear case)
参见
Theory is found in Reference [1].
[1] J. Li, T. Liu & Z. Sun, Implementation of the GRP scheme for computing radially symmetric compressible fluid flows. Journal of Computational Physics, 228.16: 5867–5887, 2009

在文件 linear_grp_solver_radial_LAG.c92 行定义.

函数调用图:
这是这个函数的调用关系图:

◆ HLL_2D_solver()

void HLL_2D_solver ( double *  F,
double *  lambda_max,
const struct i_f_var ifv_L,
const struct i_f_var ifv_R 
)

◆ linear_GRP_solver_Edir()

void linear_GRP_solver_Edir ( double *  D,
double *  U,
const struct i_f_var ifv_L,
const struct i_f_var ifv_R,
const double  eps,
const double  atc 
)

◆ linear_GRP_solver_Edir_G2D()

void linear_GRP_solver_Edir_G2D ( double *  wave_speed,
double *  D,
double *  U,
double *  U_star,
const struct i_f_var ifv_L,
const struct i_f_var ifv_R,
const double  eps,
const double  atc 
)

◆ linear_GRP_solver_Edir_Q1D()

void linear_GRP_solver_Edir_Q1D ( double *  wave_speed,
double *  D,
double *  U,
double *  U_star,
const struct i_f_var ifv_L,
const struct i_f_var ifv_R,
const double  eps,
const double  atc 
)

◆ linear_GRP_solver_LAG()

void linear_GRP_solver_LAG ( double *  D,
double *  U,
const struct i_f_var ifv_L,
const struct i_f_var ifv_R,
const double  eps,
const double  atc 
)

◆ Riemann_solver_exact()

double Riemann_solver_exact ( double *  U_star,
double *  P_star,
const double  gammaL,
const double  gammaR,
const double  u_L,
const double  u_R,
const double  p_L,
const double  p_R,
const double  c_L,
const double  c_R,
_Bool *  CRW,
const double  eps,
const double  tol,
int  N 
)

◆ Riemann_solver_exact_Ben()

double Riemann_solver_exact_Ben ( double *  U_star,
double *  P_star,
const double  gamma,
const double  u_L,
const double  u_R,
const double  p_L,
const double  p_R,
const double  c_L,
const double  c_R,
_Bool *  CRW,
const double  eps,
const double  tol,
const int  N 
)

◆ Riemann_solver_exact_Toro()

double Riemann_solver_exact_Toro ( double *  U_star,
double *  P_star,
const double  gamma,
const double  U_l,
const double  U_r,
const double  P_l,
const double  P_r,
const double  c_l,
const double  c_r,
_Bool *  CRW,
const double  eps,
const double  tol,
const int  N 
)

◆ Riemann_solver_starPU()

double Riemann_solver_starPU ( double *  U_star,
double *  P_star,
const double  GammaL,
const double  GammaR,
const double  UL,
const double  UR,
const double  PL,
const double  PR,
const double  CL,
const double  CR,
_Bool *  CRW,
const double  eps,
const double  TOLPRE,
const int  NRITER 
)

EXACT RIEMANN SOLVER FOR Two-Component γ-Law Gas

The purpose of this function is to compute the Riemann solution for pressure and velocity in the Star Region, for the time dependent one dimensional Euler equations for two-component γ-law gas.

参数
[out]U_star,P_starVelocity/Pressure in star region.
[in]UL,PL,CLInitial Velocity/Pressure/Sound_speed on left state.
[in]UR,PR,CRInitial Velocity/Pressure/Sound_speed on right state.
[in]GammaL,GammaRRatio of specific heats.
[out]CRWCentred Rarefaction Wave (CRW) Indicator of left and right waves.
  • true: CRW
  • false: Shock wave
[in]epsThe largest value can be seen as zero.
[in]TOLPRECondition value of 'gap' at the end of the iteration.
[in]NRITERMaximum iteration step (Recommended Value: 100).
返回
change: Relative pressure change after the last iteration.
作者
E. F. Toro
日期
February 1st 1999
参见
Theory is found in Chapter 4 of Reference [1].
[1] E. F. Toro, "Riemann Solvers and Numerical Methods for Fluid Dynamics". Springer-Verlag, Second Edition, 1999

在文件 riemann_solver_starPU.c115 行定义.

函数调用图:
这是这个函数的调用关系图:

◆ Roe_2D_solver()

void Roe_2D_solver ( double *  F,
double *  lambda_max,
const struct i_f_var ifv_L,
const struct i_f_var ifv_R,
const double  delta 
)

◆ Roe_HLL_solver()

void Roe_HLL_solver ( double *  V_mk,
double *  F,
double *  lambda_max,
const struct i_f_var ifv_L,
const struct i_f_var ifv_R,
const double  delta 
)

◆ Roe_solver()

void Roe_solver ( double *  F,
double *  lambda_max,
const struct i_f_var ifv_L,
const struct i_f_var ifv_R,
const double  delta 
)