diff --git a/CMakeLists.txt b/CMakeLists.txt index 5077f46..5b50286 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -25,6 +25,10 @@ option(ARTEMIS_ENABLE_OPENMP "Enable OpenMP for artemis and parthenon" OFF) option(ARTEMIS_ENABLE_COMPILE_TIMING "Enable timing of compilation of artemis" ON) option(ARTEMIS_ENABLE_ASAN "Enable AddressSanitizer to detect memory errors" OFF) +# Force -03 optimization for RelWithDebInfo +set(CMAKE_CXX_FLAGS_RELWITHDEBINFO "-O3 -g -DNDEBUG" CACHE STRING "Force -O3 in RelWithDebInfo mode" FORCE) +set(CMAKE_C_FLAGS_RELWITHDEBINFO "-O3 -g -DNDEBUG" CACHE STRING "Force -O3 in RelWithDebInfo mode" FORCE) + # Timing of compilation if (ARTEMIS_ENABLE_COMPILE_TIMING) set_property(GLOBAL PROPERTY RULE_LAUNCH_COMPILE diff --git a/src/pgen/disk.hpp b/src/pgen/disk.hpp index 7c34ea3..65247de 100644 --- a/src/pgen/disk.hpp +++ b/src/pgen/disk.hpp @@ -39,6 +39,7 @@ #include "nbody/nbody.hpp" #include "utils/artemis_utils.hpp" #include "utils/eos/eos.hpp" +#include "utils/units.hpp" using ArtemisUtils::EOS; using ArtemisUtils::VI; @@ -60,6 +61,7 @@ struct DiskParams { Real alpha, nu0, nu_indx; Real mdot; Real temp_soft2; + Real kbmu; bool do_dust; bool nbody_temp; bool quiet_start; @@ -101,7 +103,8 @@ Real TempProfile(struct DiskParams pgen, const Real R, const Real z) { const Real H = R * pgen.h0 * std::pow(R / pgen.r0, pgen.flare); const Real ir1 = 1.0 / std::sqrt(R * R + pgen.temp_soft2); const Real omk2 = SQR(pgen.Omega0) * ir1 * ir1 * ir1; - const Real T0 = omk2 * H * H / pgen.Gamma; + // c_iso^2 = P/rho = kb/mu T = Omk^2 H^2 + const Real T0 = omk2 * H * H / (pgen.kbmu * pgen.Gamma); return T0 * std::pow(rho / rho0, pgen.Gamma - 1.0); } @@ -275,6 +278,9 @@ inline void InitDiskParams(MeshBlock *pmb, ParameterInput *pin) { disk_params.l0 = pin->GetOrAddReal("problem", "l0", 0.0); disk_params.dust_to_gas = pin->GetOrAddReal("problem", "dust_to_gas", 0.01); disk_params.temp_soft2 = pin->GetOrAddReal("problem", "temp_soft", 0.0); + const auto mu = gas_pkg->Param("mu"); + auto &constants = artemis_pkg->Param("constants"); + disk_params.kbmu = constants.GetKBCode() / (mu * constants.GetAMUCode()); disk_params.do_dust = params.Get("do_dust"); @@ -312,7 +318,7 @@ inline void InitDiskParams(MeshBlock *pmb, ParameterInput *pin) { disk_params.nu0 = disk_params.alpha * disk_params.gamma_gas * SQR(disk_params.h0 * disk_params.r0 * disk_params.Omega0); disk_params.nu_indx = 1.5 + disk_params.q; - } else if (vtype == "powerlaw") { + } else if ((vtype == "powerlaw") || (vtype == "constant")) { disk_params.nu0 = pin->GetReal("gas/viscosity", "nu"); disk_params.nu_indx = pin->GetOrAddReal("gas/viscosity", "r_exp", 0.0); } else {