diff --git a/doc/wiki/Runtime-Parameters-related/Runtime-Parameters:-All.md b/doc/wiki/Runtime-Parameters-related/Runtime-Parameters:-All.md index 9129bf3665..f06172bb2f 100644 --- a/doc/wiki/Runtime-Parameters-related/Runtime-Parameters:-All.md +++ b/doc/wiki/Runtime-Parameters-related/Runtime-Parameters:-All.md @@ -327,7 +327,7 @@ For variables with `Default/Min/Max` labeled as `Depend`, click the parameter na | OPT__RES_PHASE | 0 | None | None | restriction on phase [0] ##ELBDM ONLY## | | [[ OPT__REUSE_MEMORY \| Runtime-Parameters:-Refinement#OPT__REUSE_MEMORY ]] | 2 | 0 | 2 | reuse patch memory to reduce memory fragmentation: (0=off, 1=on, 2=aggressive) [2] | | [[ OPT__RHO_INT_SCHEME \| Runtime-Parameters:-Interpolation#OPT__RHO_INT_SCHEME ]] | INT_CQUAD | 1 | 7 | ghost-zone mass density for the Poisson solver [4] | -| [[ OPT__SAME_INTERFACE_B \| Runtime-Parameters:-Hydro#OPT__SAME_INTERFACE_B ]] | 0 | None | None | ensure B field consistency on the shared interfaces between sibling patches (for debugging) [0] ##MHD ONLY## | +| [[ OPT__SAME_INTERFACE_B \| Runtime-Parameters:-Hydro#OPT__SAME_INTERFACE_B ]] | -1 | None | 1 | ensure B field consistency on the shared interfaces between sibling patches (for debugging) (<0=auto, 0=off, 1=on) [-1] ##MHD ONLY## | | [[ OPT__SELF_GRAVITY \| Runtime-Parameters:-Gravity#OPT__SELF_GRAVITY ]] | 1 | None | None | add self-gravity [1] | | [[ OPT__SORT_PATCH_BY_LBIDX \| Runtime-Parameters:-Miscellaneous#OPT__SORT_PATCH_BY_LBIDX ]] | Depend | Depend | Depend | sort patches to improve bitwise reproducibility [SERIAL:0, LOAD_BALACNE:1] | | [[ OPT__TIMING_BALANCE \| Runtime-Parameters:-Miscellaneous#OPT__TIMING_BALANCE ]] | 0 | None | None | record the max/min elapsed time in various code sections for checking load balance [0] | diff --git a/doc/wiki/Runtime-Parameters-related/Runtime-Parameters:-Hydro.md b/doc/wiki/Runtime-Parameters-related/Runtime-Parameters:-Hydro.md index 9b5c78250a..197b6aca58 100644 --- a/doc/wiki/Runtime-Parameters-related/Runtime-Parameters:-Hydro.md +++ b/doc/wiki/Runtime-Parameters-related/Runtime-Parameters:-Hydro.md @@ -155,7 +155,8 @@ solutions. * **Description:** Maximum number of iterations to reduce [MINMOD_COEFF](#MINMOD_COEFF) when data reconstruction fails. It improves code stability but -may break strict conservation. +may break strict conservation and magnetic field consistency at +patch interfaces. * **Restriction:** @@ -196,7 +197,7 @@ Only applicable when enabling the compilation option [[--dual | Installation:-Option-List#--dual]]. -* #### `OPT__SAME_INTERFACE_B`   (0=off, 1=on)   [0] +* #### `OPT__SAME_INTERFACE_B`   (-1=set to default, 0=off, 1=on)   [depend] * **Description:** Ensure adjacent patches on the same level have exactly the same magnetic field on their shared interfaces (including round-off errors). diff --git a/example/test_problem/Template/Input__Parameter b/example/test_problem/Template/Input__Parameter index d7742ab30f..07f3651b45 100644 --- a/example/test_problem/Template/Input__Parameter +++ b/example/test_problem/Template/Input__Parameter @@ -257,7 +257,7 @@ OPT__1ST_FLUX_CORR -1 # correct unphysical results (defined # (<0=auto, 0=off, 1=3D, 2=3D+1D) [-1] ##MHM/MHM_RP/CTU ONLY## OPT__1ST_FLUX_CORR_SCHEME -1 # Riemann solver for OPT__1ST_FLUX_CORR (-1=auto, 0=none, 1=Roe, 2=HLLC, 3=HLLE, 4=HLLD) [-1] DUAL_ENERGY_SWITCH 2.0e-2 # apply dual-energy if E_int/E_kin < DUAL_ENERGY_SWITCH [2.0e-2] ##DUAL_ENERGY ONLY## -OPT__SAME_INTERFACE_B 0 # ensure B field consistency on the shared interfaces between sibling patches (for debugging) [0] ##MHD ONLY## +OPT__SAME_INTERFACE_B -1 # ensure B field consistency on the shared interfaces between sibling patches (for debugging) (<0=auto, 0=off, 1=on) [-1] ##MHD ONLY## # fluid solver in ELBDM (MODEL==ELBDM only) diff --git a/include/Global.h b/include/Global.h index ac1d5eb66e..cf42027e92 100644 --- a/include/Global.h +++ b/include/Global.h @@ -133,7 +133,7 @@ extern bool OPT__FIXUP_ELECTRIC, OPT__CK_INTERFACE_B, OPT__OUTPUT_CC extern bool OPT__OUTPUT_DIVMAG; extern int OPT__CK_DIVERGENCE_B; extern double UNIT_B; -extern bool OPT__SAME_INTERFACE_B; +extern int OPT__SAME_INTERFACE_B; extern OptInitMagByVecPot_t OPT__INIT_BFIELD_BYVECPOT; #endif diff --git a/include/Typedef.h b/include/Typedef.h index fe1670b400..9748c0a8d0 100644 --- a/include/Typedef.h +++ b/include/Typedef.h @@ -561,6 +561,13 @@ const LoadParaMode_t LOAD_HDF5_OUTPUT = 2; +// whether to enforce consistency of magnetic field at the patch interface +typedef int SameInterfaceB_t; +const SameInterfaceB_t + SAME_INTERFACE_B_NO = 0, + SAME_INTERFACE_B_YES = 1; + + // function pointers typedef real (*EoS_GUESS_t) ( const real Con[], real* const Constant, const double AuxArray_Flt[], const int AuxArray_Int[], const real *const Table[EOS_NTABLE_MAX] ); diff --git a/src/Init/Init_GAMER.cpp b/src/Init/Init_GAMER.cpp index d6d48c9313..46a95c3262 100644 --- a/src/Init/Init_GAMER.cpp +++ b/src/Init/Init_GAMER.cpp @@ -246,7 +246,7 @@ void Init_GAMER( int *argc, char ***argv ) // ensure B field consistency on the shared interfaces between sibling patches # if ( MODEL == HYDRO && defined MHD ) - if ( OPT__SAME_INTERFACE_B ) + if ( OPT__SAME_INTERFACE_B == SAME_INTERFACE_B_YES ) for (int lv=0; lvAdd( "DUAL_ENERGY_SWITCH", &DUAL_ENERGY_SWITCH, 2.0e-2, 0.0, NoMax_double ); # endif # ifdef MHD - ReadPara->Add( "OPT__SAME_INTERFACE_B", &OPT__SAME_INTERFACE_B, false, Useless_bool, Useless_bool ); + ReadPara->Add( "OPT__SAME_INTERFACE_B", &OPT__SAME_INTERFACE_B, -1, NoMin_int, 1 ); # endif # endif // #if ( MODEL == HYDRO ) diff --git a/src/Init/Init_ResetParameter.cpp b/src/Init/Init_ResetParameter.cpp index 4ecd59889f..46ad214160 100644 --- a/src/Init/Init_ResetParameter.cpp +++ b/src/Init/Init_ResetParameter.cpp @@ -324,6 +324,15 @@ void Init_ResetParameter() # endif +# ifdef MHD + if ( OPT__SAME_INTERFACE_B < 0 ) + { + OPT__SAME_INTERFACE_B = ( MINMOD_MAX_ITER > 0 ) ? SAME_INTERFACE_B_YES : SAME_INTERFACE_B_NO; + + PRINT_RESET_PARA( OPT__SAME_INTERFACE_B, FORMAT_INT, "" ); + } +# endif + // text format parameters // --> The current strategy is to read the integer in between % and . to determine the string length. // For example, the format %20.16e will give a length of 20. However, if only checking the string after %, diff --git a/src/Main/EvolveLevel.cpp b/src/Main/EvolveLevel.cpp index 3305379464..fcaadfae0c 100644 --- a/src/Main/EvolveLevel.cpp +++ b/src/Main/EvolveLevel.cpp @@ -303,6 +303,22 @@ void EvolveLevel( const int lv, const double dTime_FaLv ) } // if ( OPT__OVERLAP_MPI ) ... else ... +# if ( MODEL == HYDRO && defined MHD ) + if ( OPT__SAME_INTERFACE_B == SAME_INTERFACE_B_YES ) + { + TIMING_FUNC( Buf_GetBufferData( lv, SaveSg_Flu, SaveSg_Mag, NULL_INT, DATA_GENERAL, + _ENGY, _MAG, Flu_ParaBuf, USELB_YES ), + Timer_GetBuf[lv][0], TIMER_ON ); + + TIMING_FUNC( MHD_SameInterfaceB( lv ), + Timer_Flu_Advance[lv], TIMER_ON ); + + TIMING_FUNC( Buf_GetBufferData( lv, SaveSg_Flu, SaveSg_Mag, NULL_INT, DATA_GENERAL, + _ENGY, _MAG, Flu_ParaBuf, USELB_YES ), + Timer_GetBuf[lv][0], TIMER_ON ); + } // if ( OPT__SAME_INTERFACE_B == SAME_INTERFACE_B_YES ) +# endif // #if ( MODEL == HYDRO && defined MHD ) + amr->FluSg [lv] = SaveSg_Flu; amr->FluSgTime[lv][SaveSg_Flu] = TimeNew; # ifdef MHD diff --git a/src/Main/Main.cpp b/src/Main/Main.cpp index a495c688a5..96176dfbdb 100644 --- a/src/Main/Main.cpp +++ b/src/Main/Main.cpp @@ -123,7 +123,7 @@ bool OPT__FIXUP_ELECTRIC, OPT__CK_INTERFACE_B, OPT__OUTPUT_CC_MA bool OPT__OUTPUT_DIVMAG; int OPT__CK_DIVERGENCE_B; double UNIT_B; -bool OPT__SAME_INTERFACE_B; +int OPT__SAME_INTERFACE_B; OptInitMagByVecPot_t OPT__INIT_BFIELD_BYVECPOT; #endif @@ -698,7 +698,7 @@ int main( int argc, char *argv[] ) TIMING_FUNC( Flu_CorrAfterAllSync(), Timer_Main[6], TIMER_ON ); # if ( MODEL == HYDRO && defined MHD ) - if ( OPT__SAME_INTERFACE_B ) + if ( OPT__SAME_INTERFACE_B == SAME_INTERFACE_B_YES ) { if ( OPT__VERBOSE && MPI_Rank == 0 ) Aux_Message( stdout, " MHD_SameInterfaceB ... " ); diff --git a/src/Model_Hydro/MHD_SameInterfaceB.cpp b/src/Model_Hydro/MHD_SameInterfaceB.cpp index 657cf0fb08..e3f86503c1 100644 --- a/src/Model_Hydro/MHD_SameInterfaceB.cpp +++ b/src/Model_Hydro/MHD_SameInterfaceB.cpp @@ -8,7 +8,7 @@ //------------------------------------------------------------------------------------------------------- // Function : MHD_SameInterfaceB // Description : Ensure that adjacent patches on the same level (i.e., sibling patches) have *exactly* -// the same B field on their shared interfaces +// the same B field and B energy on their shared interfaces // --> Even round-off errors are the same // // Note : 1. Applied to both real and buffer patches @@ -32,24 +32,64 @@ void MHD_SameInterfaceB( const int lv ) # endif + const int FluSg = amr->FluSg[lv]; const int MagSg = amr->MagSg[lv]; // iterate over all real and buffer patches # pragma omp parallel for schedule( runtime ) for (int PID=0; PIDnum[lv]; PID++) { + real ***Emag_old = NULL; + + Aux_AllocateArray3D( Emag_old, 3, PS1, PS1 ); + // use the B field on the +x/+y/+z sides to overwrite that on the -x/-y/-z sides // --> skip s=1/3/5 for (int s=0; s<6; s+=2) { + const int d = s/2; const int SibPID = amr->patch[0][lv][PID]->sibling[s]; + if ( SibPID < 0 ) continue; + // some buffer patches may have magnetic == NULL --> skip them - if ( SibPID >= 0 && - amr->patch[MagSg][lv][ PID]->magnetic != NULL && - amr->patch[MagSg][lv][SibPID]->magnetic != NULL ) - MHD_CopyPatchInterfaceBField( lv, PID, s, MagSg ); - } + if ( amr->patch[MagSg][lv][ PID]->magnetic == NULL || + amr->patch[MagSg][lv][SibPID]->magnetic == NULL ) continue; + + int ii, jj, kk; + + for (int j=0; jpatch[FluSg][lv][PID]->fluid[ENGY][kk][jj][ii] += (Emag - Emag_old[d][j][i]); + } + } // for (int s=0; s<6; s+=2) + + Aux_DeallocateArray3D( Emag_old ); } // for (int PID=0; PIDnum[lv]; PID++) } // FUNCTION : MHD_SameInterfaceB