From 65d69027be1bf10dafa9e0b74db7f0d2cb5c1096 Mon Sep 17 00:00:00 2001 From: Miikka Vaisala Date: Tue, 3 Sep 2019 17:48:56 +0800 Subject: [PATCH] Found an error in the gravitational constant. Now corrected! --- src/standalone/config_loader.cc | 7 +++---- src/standalone/model/host_memory.cc | 10 ++++++---- 2 files changed, 9 insertions(+), 8 deletions(-) diff --git a/src/standalone/config_loader.cc b/src/standalone/config_loader.cc index cc50883..3df7d12 100644 --- a/src/standalone/config_loader.cc +++ b/src/standalone/config_loader.cc @@ -125,7 +125,7 @@ update_config(AcMeshInfo* config) config->real_params[AC_gamma]; AcReal G_CONST_CGS = AcReal( - 6.674e-8); // g/cm3/s GGS definition //TODO define in a separate module + 6.674e-8); // cm^3/(g*s^2) GGS definition //TODO define in a separate module AcReal M_sun = AcReal(1.989e33); // g solar mass @@ -137,8 +137,6 @@ update_config(AcMeshInfo* config) - - config->real_params[AC_M_sink] = config->real_params[AC_M_sink_Msun] * M_sun / config->real_params[AC_unit_mass]; config->real_params[AC_M_sink_init] = config->real_params[AC_M_sink_Msun] * M_sun / @@ -147,7 +145,8 @@ update_config(AcMeshInfo* config) config->real_params[AC_G_const] = G_CONST_CGS / ((config->real_params[AC_unit_velocity] * config->real_params[AC_unit_velocity]) / (config->real_params[AC_unit_density] * - config->real_params[AC_unit_length])); + config->real_params[AC_unit_length] * + config->real_params[AC_unit_length] )); config->real_params[AC_sq2GM_star] = AcReal(sqrt(AcReal(2) * config->real_params[AC_GM_star])); diff --git a/src/standalone/model/host_memory.cc b/src/standalone/model/host_memory.cc index 21d64e3..ac47f28 100644 --- a/src/standalone/model/host_memory.cc +++ b/src/standalone/model/host_memory.cc @@ -227,16 +227,18 @@ simple_uniform_core(AcMesh* mesh) const double zorig = mesh->info.real_params[AC_zorig]; double xx, yy, zz, RR; double delx, dely, delz; - double tanhRR; + double core_profile; //TEMPORARY TEST INPUT PARAMETERS const double core_radius = DX*32.0; const double trans = DX*12.0; + const double epsilon = DX*2.0; const double vel_scale = mesh->info.real_params[AC_ampl_uu]; double abso_vel; RR = 1.0; printf("%e %e %e \n", RR, trans, core_radius); + for (int k = 0; k < mz; k++) { for (int j = 0; j < my; j++) { @@ -251,9 +253,9 @@ simple_uniform_core(AcMesh* mesh) delz = zz; RR = sqrt(delx*delx + dely*dely + delz*delz); - //tanhRR = double(-0.5)*tanh((core_radius-RR)/trans) + double(0.5) + //core_profile = double(-0.5)*tanh((core_radius-RR)/trans) + double(0.5) // + double(0.1); - tanhRR = double(1.0); + core_profile = pow(RR+epsilon, -2.0); //double(1.0); AcReal RR_inner_bound = mesh->info.real_params[AC_soft]/AcReal(2.0); if (RR >= RR_inner_bound) { @@ -268,7 +270,7 @@ simple_uniform_core(AcMesh* mesh) } - mesh->vertex_buffer[VTXBUF_LNRHO][idx] = log(exp(ampl_lnrho)*tanhRR); + mesh->vertex_buffer[VTXBUF_LNRHO][idx] = log(exp(ampl_lnrho)*core_profile); mesh->vertex_buffer[VTXBUF_UUX][idx] = -abso_vel * (yy / RR); mesh->vertex_buffer[VTXBUF_UUY][idx] = abso_vel * (xx / RR); mesh->vertex_buffer[VTXBUF_UUZ][idx] = double(0.0);