Added function for truelove-density, and took out the unecessary use of switch on-off module for accretion.

This commit is contained in:
JackHsu
2019-08-06 16:17:11 +08:00
parent 632819764e
commit aebf4f11ae
2 changed files with 21 additions and 15 deletions

View File

@@ -33,7 +33,6 @@
#define LFORCING (0) #define LFORCING (0)
#define LUPWD (0) #define LUPWD (0)
#define LSINK (1) #define LSINK (1)
#define LACCRETION(1)
#define AC_THERMAL_CONDUCTIVITY (AcReal(0.001)) // TODO: make an actual config parameter #define AC_THERMAL_CONDUCTIVITY (AcReal(0.001)) // TODO: make an actual config parameter

View File

@@ -76,23 +76,30 @@ sink_gravity(int3 globalVertexIdx){
return force_gravity; return force_gravity;
} }
#endif #endif
// Now this is more of a suedo-code, just write out some intuitions. // Now this is more of a suedo-code, just write out some intuitions.
#if LACCRETION #if LSINK
Scalar Scalar
mass_accrection(int3 globalVertexIdx){ truelove_density(in Scalar lnrho){
Scalar mass = (Scalar){() const Scalar rho = exp(lnrho);
Scalar Jeans_length = pow(((3.14 * cv_sound * cv_sound) / (AC_G_const * lnrho), 0.5)); const Scalar Jeans_length_squared = (M_PI * cs2_sound) / (AC_G_const * rho);
Scalar Truelove_Jeans_lnrho = ((3.14) * Jeans_length * Jeans_length * cv2_sound) / (AC_G_const * dsx * dsy * dsz); const Scalar TJ_rho = ((M_PI) * ((dsx * dsx) / Jeans_length_squared) * cs2_sound) / (AC_G_const * dsx * dsx);
const Vector sink_pos = (Vector){DCONST_REAL(AC_sink_pos_x),
DCONST_REAL(AC_sink_pos_y),
DCONST_REAL(AC_sink_pos_z)}; Scalar accrection_rho = rho - TJ_rho;
const Scalar sink_mass = DCONST_REAL(AC_M_sink); if (accrection_rho < 0){
if (lnrho >= Truelove_Jeans_lnrho) accrection_rho = Scalar(0);
{
sink_mass = sink_mass + (Truelove_Jeans_lnrho * (dsx * dsy * dsz))
} }
return sink_mass;
return accretion_rho;
} }
#endif
//TODO: basic structure of this part is as follows
// update_accretion_buffer() <--> accretion_profile() <--> truelove_density()
#if LENTROPY #if LENTROPY
Vector Vector
momentum(int3 globalVertexIdx, in Vector uu, in Scalar lnrho, in Scalar ss, in Vector aa) { momentum(int3 globalVertexIdx, in Vector uu, in Scalar lnrho, in Scalar ss, in Vector aa) {