From 89128af44b58b394854a3e3d8cd4868f886d882f Mon Sep 17 00:00:00 2001 From: JackHsu Date: Thu, 8 Aug 2019 15:53:41 +0800 Subject: [PATCH] added update_accretion_buffer. --- acc/mhd_solver/stencil_defines.h | 1 + acc/mhd_solver/stencil_process.sps | 43 +++++++++++++++--------------- config/astaroth.conf | 2 +- 3 files changed, 24 insertions(+), 22 deletions(-) diff --git a/acc/mhd_solver/stencil_defines.h b/acc/mhd_solver/stencil_defines.h index 0add153..6e22bc9 100644 --- a/acc/mhd_solver/stencil_defines.h +++ b/acc/mhd_solver/stencil_defines.h @@ -80,6 +80,7 @@ FUNC(AC_M_sink),\ FUNC(AC_M_sink_Msun),\ FUNC(AC_soft),\ + FUNC(AC_accretion_range),\ /* Run params */\ FUNC(AC_cdt), \ FUNC(AC_cdtv), \ diff --git a/acc/mhd_solver/stencil_process.sps b/acc/mhd_solver/stencil_process.sps index 65907c4..fce21bf 100644 --- a/acc/mhd_solver/stencil_process.sps +++ b/acc/mhd_solver/stencil_process.sps @@ -82,14 +82,14 @@ sink_gravity(int3 globalVertexIdx){ #if LSINK Scalar truelove_density(in Scalar lnrho){ - const Scalar rho = exp(lnrho); + const Scalar rho = exp(value(lnrho)); const Scalar Jeans_length_squared = (M_PI * cs2_sound) / (AC_G_const * rho); const Scalar TJ_rho = ((M_PI) * ((dsx * dsx) / Jeans_length_squared) * cs2_sound) / (AC_G_const * dsx * dsx); //TODO: dsx will cancel out, deal with it later for optimization. - Scalar accrection_rho = rho - TJ_rho; - if (accrection_rho < 0){ - accrection_rho = Scalar(0); + Scalar accretion_rho = rho - TJ_rho; + if (accretion_rho < 0){ + accretion_rho = Scalar(0); } return accretion_rho; @@ -97,23 +97,24 @@ truelove_density(in Scalar lnrho){ Scalar accretion_profile(int3 globalVertexIdx, in Scalar lnrho){ - // QUESTION: do I need to define grid_pos, sink_pos and distance again - // if the sink_gravity kernel will also be called once LSINK swtich is on? Seems redundant. - const Vector grid_pos = (Vector){(globalVertexIdx.x - nx_min) * dsx, - (globalVertexIdx.y - ny_min) * dsy, - (globalVertexIdx.z - nz_min) * dsz}; - const Vector sink_pos = (Vector){DCONST_REAL(AC_sink_pos_x), - DCONST_REAL(AC_sink_pos_y), - DCONST_REAL(AC_sink_pos_z)}; - const Scalar profile_range = DCONST_REAL(AC_accretion_range); - const Scalar accretion_distance = length(grid_pos - sink_pos); - if ((accretion_distance) <= profile_range){ + // QUESTION: do I need to define grid_pos, sink_pos and distance again + // if the sink_gravity kernel will also be called once LSINK swtich is on? Seems redundant. + const Vector grid_pos = (Vector){(globalVertexIdx.x - nx_min) * dsx, + (globalVertexIdx.y - ny_min) * dsy, + (globalVertexIdx.z - nz_min) * dsz}; + const Vector sink_pos = (Vector){DCONST_REAL(AC_sink_pos_x), + DCONST_REAL(AC_sink_pos_y), + DCONST_REAL(AC_sink_pos_z)}; + const Scalar profile_range = DCONST_REAL(AC_accretion_range); + const Scalar accretion_distance = length(grid_pos - sink_pos); + Scalar accretion_density; + if ((accretion_distance) <= profile_range){ // calculate accretion according to chosen criterion for the grid cell. - accretion_density = truelove_density(lnrho); - } else { - accretion_density = Scalar(0.0); - } - return accretion_density; + accretion_density = truelove_density(lnrho); + } else { + accretion_density = Scalar(0.0); + } + return accretion_density; } #endif @@ -417,7 +418,7 @@ solve(Scalar dt) { #if LSINK out_lnrho = log(exp(out_lnrho) - accretion_profile(globalVertexIdx, lnrho)); - out_accretion = accretion + (accretion_profile(globalVertexIdx,lnrho) * dsx * dsy * dsz); + out_accretion = value(accretion) + (accretion_profile(globalVertexIdx,lnrho) * dsx * dsy * dsz); #endif } diff --git a/config/astaroth.conf b/config/astaroth.conf index 4a127db..a003c8b 100644 --- a/config/astaroth.conf +++ b/config/astaroth.conf @@ -58,7 +58,7 @@ AC_soft = 0.12 // Accretion Parameters // profile_range is multiple of dsx -AC_profile_range = 2.0 +AC_accretion_range = 2.0 // Physical properties of the domain AC_unit_velocity = 1.0