Skip to content

Commit

Permalink
More code for outputting deltaP
Browse files Browse the repository at this point in the history
  • Loading branch information
jzuhone committed Jan 24, 2025
1 parent 563020d commit 637bad0
Show file tree
Hide file tree
Showing 9 changed files with 41 additions and 18 deletions.
1 change: 1 addition & 0 deletions include/Global.h
Original file line number Diff line number Diff line change
Expand Up @@ -403,6 +403,7 @@ extern double VISCOSITY_CONSTANT_COEFF;
extern double VISCOSITY_SPITZER_FRAC;
extern double VISCOSITY_COULOMB_LOG;
extern double VISCOSITY_MAX_DIFFUSIVITY;
extern bool OPT__OUTPUT_DELTAP;
#endif

#ifdef CONDUCTION
Expand Down
3 changes: 3 additions & 0 deletions include/HDF5_Typedef.h
Original file line number Diff line number Diff line change
Expand Up @@ -853,6 +853,9 @@ struct InputPara_t
# ifdef MHD
int Opt__Output_DivMag;
# endif
# ifdef VISCOSITY
int Opt__Output_DeltaP;
# endif
# ifdef SRHD
int Opt__Output_Lorentz;
int Opt__Output_3Velocity;
Expand Down
2 changes: 2 additions & 0 deletions include/Prototype.h
Original file line number Diff line number Diff line change
Expand Up @@ -159,6 +159,8 @@ real Hydro_DensDual2Pres( const real Dens, const real Dual, const real Gamma_m1,
#ifdef VISCOSITY
void Hydro_ComputeViscosity( real &visc_mu, real &visc_nu, const MicroPhy_t *MicroPhy,
const real Dens, const real Temp );
void Hydro_Compute_DeltaP( real Out[], const real FluIn[][CUBE(DER_NXT)], const real faceB[][CUBE(DER_NXT)],
const real Temp[], const int NGhost, const real dh, const MicroPhy_t *MicroPhy );
#endif
#ifdef CONDUCTION
void Hydro_ComputeConduction( real &cond_kappa, real &cond_chi, const MicroPhy_t *MicroPhy,
Expand Down
3 changes: 3 additions & 0 deletions src/Auxiliary/Aux_TakeNote.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -1658,6 +1658,9 @@ void Aux_TakeNote()
# endif
# ifdef MHD
fprintf( Note, "OPT__OUTPUT_DIVMAG % d\n", OPT__OUTPUT_DIVMAG );
# endif
# ifdef VISCOSITY
fprintf( Note, "OPT__OUTPUT_DELTAP % d\n", OPT__OUTPUT_DELTAP );
# endif
fprintf( Note, "OPT__OUTPUT_USER_FIELD % d\n", OPT__OUTPUT_USER_FIELD );
# ifdef SRHD
Expand Down
3 changes: 3 additions & 0 deletions src/Init/Init_ByRestart_HDF5.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -2397,6 +2397,9 @@ void Check_InputPara( const char *FileName, const int FormatVersion )
# ifdef MHD
LoadField( "Opt__Output_DivMag", &RS.Opt__Output_DivMag, SID, TID, NonFatal, &RT.Opt__Output_DivMag, 1, NonFatal );
# endif
# ifdef VISCOSITY
LoadField( "Opt__Output_DeltaP", &RS.Opt__Output_DeltaP, SID, TID, NonFatal, &RT.Opt__Output_DeltaP, 1, NonFatal );
# endif
# ifdef SRHD
LoadField( "Opt__Output_Lorentz", &RS.Opt__Output_Lorentz, SID, TID, NonFatal, &RT.Opt__Output_Lorentz, 1, NonFatal );
LoadField( "Opt__Output_3Velocity", &RS.Opt__Output_3Velocity, SID, TID, NonFatal, &RT.Opt__Output_3Velocity, 1, NonFatal );
Expand Down
3 changes: 3 additions & 0 deletions src/Init/Init_Load_Parameter.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -568,6 +568,9 @@ void Init_Load_Parameter()
# ifdef MHD
ReadPara->Add( "OPT__OUTPUT_DIVMAG", &OPT__OUTPUT_DIVMAG, false, Useless_bool, Useless_bool );
# endif
# ifdef VISCOSITY
ReadPara->Add( "OPT__OUTPUT_DELTAP", &OPT__OUTPUT_DELTAP, false, Useless_bool, Useless_bool );
# endif
# ifdef SRHD
ReadPara->Add( "OPT__OUTPUT_LORENTZ", &OPT__OUTPUT_LORENTZ, false, Useless_bool, Useless_bool );
ReadPara->Add( "OPT__OUTPUT_3VELOCITY", &OPT__OUTPUT_3VELOCITY, false, Useless_bool, Useless_bool );
Expand Down
1 change: 1 addition & 0 deletions src/Main/Main.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -376,6 +376,7 @@ double VISCOSITY_CONSTANT_COEFF;
double VISCOSITY_SPITZER_FRAC;
double VISCOSITY_COULOMB_LOG;
double VISCOSITY_MAX_DIFFUSIVITY;
bool OPT__OUTPUT_DELTAP;
#endif

// d. conduction
Expand Down
33 changes: 17 additions & 16 deletions src/Microphysics/Viscosity/CPU_AddViscousFlux.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -426,15 +426,18 @@ void Hydro_AddViscousFlux( const real g_ConVar[][ CUBE(FLU_NXT) ],


#ifdef MHD
void Hydro_Compute_DeltaP( real Out[], const real FluIn[], const real faceB[], const real Temp[],
const int NGhost, const real dh )
void Hydro_Compute_DeltaP( real Out[], const real FluIn[][CUBE(DER_NXT)], const real faceB[][CUBE(DER_NXT)],
const real Temp[], const int NGhost, const real dh, const MicroPhy_t *MicroPhy )
{
real BBdV, delta_p, divV;
real dens_L, dens_R, temp_L, temp_R, visc_nu;
real N_slope_N, T1_slope_N, T2_slope_N;
real N_slope_T1, T1_slope_T1, T2_slope_T1;
real N_slope_T2, T1_slope_T2, T2_slope_T2;
real B_amp, B_N_mean, B_T1_mean, B_T2_mean, B2;
real mu, mu_l, mu_r;
int idx_fN, idx_fT1, idx_fT2;
int sizeB_i, sizeB_j, sizeB_k;

const real _dh = (real)1.0 / dh;
const int didx[3] = { 1, DER_NXT, SQR(DER_NXT) };
Expand Down Expand Up @@ -465,9 +468,9 @@ void Hydro_Compute_DeltaP( real Out[], const real FluIn[], const real faceB[], c
{
const int TDir1 = (d+1)%3; // transverse direction 1
const int TDir2 = (d+2)%3; // transverse direction 2
const int dB = d + MAGX;
const int TD1B = TDir1 + MAGX;
const int TD2B = TDir2 + MAGX;
//const int dB = d + MAGX;
//const int TD1B = TDir1 + MAGX;
//const int TD2B = TDir2 + MAGX;

switch ( d )
{
Expand All @@ -487,8 +490,6 @@ void Hydro_Compute_DeltaP( real Out[], const real FluIn[], const real faceB[], c

}

const int idxol = IDX321 ( i, j, k, DER_NXT, );

const int stride_fc_BT1[3] = { 1, sizeB_k, sizeB_k*sizeB_i };
const int stride_fc_BT2[3] = { 1, sizeB_j, sizeB_j*sizeB_k };

Expand Down Expand Up @@ -562,15 +563,15 @@ void Hydro_Compute_DeltaP( real Out[], const real FluIn[], const real faceB[], c
) * _dh;

// compute the mean magnetic field at the face-centered flux location
B_N_mean = faceB[ dB][ idx_fN ];
B_T1_mean = (real)0.25*( faceB[TD1B][ idx_fT1 ] +
faceB[TD1B][ idx_fT1 + stride_fc_BT1[TDir1] ] +
faceB[TD1B][ idx_fT1 - stride_fc_BT1[d] ] +
faceB[TD1B][ idx_fT1 - stride_fc_BT1[d] + stride_fc_BT1[TDir1] ] );
B_T2_mean = (real)0.25*( faceB[TD2B][ idx_fT2 ] +
faceB[TD2B][ idx_fT2 + stride_fc_BT2[TDir2] ] +
faceB[TD2B][ idx_fT2 - stride_fc_BT2[d] ] +
faceB[TD2B][ idx_fT2 - stride_fc_BT2[d] + stride_fc_BT2[TDir2] ] );
B_N_mean = faceB[ d][ idx_fN ];
B_T1_mean = (real)0.25*( faceB[TDir1][ idx_fT1 ] +
faceB[TDir1][ idx_fT1 + stride_fc_BT1[TDir1] ] +
faceB[TDir1][ idx_fT1 - stride_fc_BT1[d] ] +
faceB[TDir1][ idx_fT1 - stride_fc_BT1[d] + stride_fc_BT1[TDir1] ] );
B_T2_mean = (real)0.25*( faceB[TDir2][ idx_fT2 ] +
faceB[TDir2][ idx_fT2 + stride_fc_BT2[TDir2] ] +
faceB[TDir2][ idx_fT2 - stride_fc_BT2[d] ] +
faceB[TDir2][ idx_fT2 - stride_fc_BT2[d] + stride_fc_BT2[TDir2] ] );

B2 = SQR(B_N_mean) + SQR(B_T1_mean) + SQR(B_T2_mean);
B_amp = SQRT( B2 );
Expand Down
10 changes: 8 additions & 2 deletions src/Output/Output_DumpData_Total_HDF5.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -1080,8 +1080,8 @@ void Output_DumpData_Total_HDF5( const char *FileName )
}
// compute and store the target derived field
const int PID = PID0 + LocalID;
Hydro_Compute_DeltaP( FieldData[PID][0][0], Der_FluIn[LocalID][0], Der_MagCC[0],
NDer, DER_NXT, DER_NXT, DER_NXT, DER_GHOST_SIZE, amr->dh[lv] );
Hydro_Compute_DeltaP( FieldData[PID][0][0], Der_FluIn[LocalID][0], Der_MagFC[LocalID][0],
Temp, DER_GHOST_SIZE, amr->dh[lv], MicroPhy );
} // for (int LocalID=0; LocalID<8; LocalID++)
} // for (int PID0=0; PID0<amr->NPatchComma[lv][1]; PID0+=8)
} // if ( v == DeltaPDumpIdx || v == KappaDumpIdx )
Expand Down Expand Up @@ -2864,6 +2864,9 @@ void FillIn_InputPara( InputPara_t &InputPara, const int NFieldStored, char Fiel
# ifdef MHD
InputPara.Opt__Output_DivMag = OPT__OUTPUT_DIVMAG;
# endif
# ifdef VISCOSITY
InputPara.Opt__Output_DeltaP = OPT__OUTPUT_DELTAP;
# endif
# ifdef SRHD
InputPara.Opt__Output_3Velocity = OPT__OUTPUT_3VELOCITY;
InputPara.Opt__Output_Lorentz = OPT__OUTPUT_LORENTZ;
Expand Down Expand Up @@ -3955,6 +3958,9 @@ void GetCompound_InputPara( hid_t &H5_TypeID, const int NFieldStored )
# ifdef MHD
H5Tinsert( H5_TypeID, "Opt__Output_DivMag", HOFFSET(InputPara_t,Opt__Output_DivMag ), H5T_NATIVE_INT );
# endif
# ifdef VISCOSITY
H5Tinsert( H5_TypeID, "Opt__Output_DeltaP", HOFFSET(InputPara_t,Opt__Output_DeltaP ), H5T_NATIVE_INT );
# endif
# ifdef SRHD
H5Tinsert( H5_TypeID, "Opt__Output_3Velocity", HOFFSET(InputPara_t,Opt__Output_3Velocity ), H5T_NATIVE_INT );
H5Tinsert( H5_TypeID, "Opt__Output_Lorentz", HOFFSET(InputPara_t,Opt__Output_Lorentz ), H5T_NATIVE_INT );
Expand Down

0 comments on commit 637bad0

Please sign in to comment.