Skip to content

Reintroduce Dexter rhoV with grtrans mitigation #88

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
c-prather wants to merge 1 commit into
base: master
Choose a base branch
from
Open

Reintroduce Dexter rhoV with grtrans mitigation #88

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
8 changes: 4 additions & 4 deletions src/model_radiation.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -37,9 +37,7 @@
#define E_LEUNG 5
#define E_CUSTOM 10
// Rotation
#define ROT_OLD 11
#define ROT_PIECEWISE 12
#define ROT_SHCHERBAKOV 13
#define ROT_DEXTER 11
// Debugging/internal
#define E_UNPOL 15

Expand Down Expand Up @@ -298,7 +296,9 @@ void jar_calc_dist(int dist, int pol, double X[NDIM], double Kcon[NDIM],
}

// ROTATIVITIES
paramsM.dexter_fit = 0; // Don't use the Dexter rhoV, as it's unstable at low temperature
if (dist != ROT_DEXTER) {
paramsM.dexter_fit = 0; // Don't use the Dexter rhoV, as it's unstable at low temperature
}
*rQ = rho_nu_fit(&paramsM, paramsM.STOKES_Q) * nu;
*rU = rho_nu_fit(&paramsM, paramsM.STOKES_U) * nu;
*rV = rho_nu_fit(&paramsM, paramsM.STOKES_V) * nu;
Expand Down
11 changes: 8 additions & 3 deletions src/symphony/maxwell_juettner_fits.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -226,7 +226,7 @@ double maxwell_juettner_rho_Q(struct parameters *params)

double k1 = gsl_sf_bessel_Kn(1, 1./params->theta_e);
double k2 = gsl_sf_bessel_Kn(2, 1./params->theta_e);
double k_ratio = (k2 > 0) ? k1/k2 : 1;
double k_ratio = ((params->dexter_fit && params->theta_e <= 0.01) || (k2 <= 0)) ? 1 : k1/k2;

double eps11m22 = jffunc * wp2 * pow(omega0, 2.)
/ pow(2.*params->pi * params->nu, 4.)
Expand Down Expand Up @@ -268,8 +268,13 @@ double maxwell_juettner_rho_V(struct parameters * params)
double fit_factor = 0;
if (params->dexter_fit && k2 > 0) { // TODO Further limit the usage here to match grtrans
// Jason Dexter (2016) fits using the modified difference factor g(X)
double shgmfunc = 0.43793091 * log(1. + 0.00185777 * pow(x, 1.50316886));
fit_factor = (k0 - shgmfunc) / k2; // TODO might be unstable way to phrase
if (params->theta_e > 0.01) {
double shgmfunc = 0.43793091 * log(1. + 0.00185777 * pow(x, 1.50316886));
double step = 0.5 + 0.5 * tanh((params->theta_e - 1.0) / 0.05);
fit_factor = (k0 - step*shgmfunc) / k2; // TODO might be unstable way to phrase
} else {
fit_factor = 1;
}
} else {
// Shcherbakov fits. Good to the smallest Thetae at high freq but questionable for low frequencies
double shgmfunc = 1 - 0.11*log(1 + 0.035*x);
Expand Down