cwpearson/astaroth
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Enabled upwinding by default and updated the model helical forcing with the hotfixed changes from earlier commits. Autotests kinda pass (we get 1 failure but this is likely due to inaccuracies of the trigonometric functions used in helical forcing. The error is very close to the acceptable error bound).

Browse Source
This commit is contained in:
jpekkila
2019-08-07 15:53:32 +03:00
parent 0bb568642f
commit b61617ee0f
2 changed files with 27 additions and 27 deletions
Download Patch File Download Diff File Expand all files Collapse all files

3
acc/mhd_solver/stencil_defines.h
View File

@@ -31,8 +31,7 @@
#define LENTROPY (1) #define LENTROPY (1)
#define LTEMPERATURE (0) #define LTEMPERATURE (0)
#define LFORCING (1) #define LFORCING (1)
#define LUPWD (0) #define LUPWD (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

29
src/standalone/model/model_rk3.cc
View File

@@ -760,23 +760,24 @@ helical_forcing(ModelScalar magnitude, ModelVector k_force, ModelVector xx, Mode
ModelVector ff_im, ModelScalar phi) ModelVector ff_im, ModelScalar phi)
{ {
(void)magnitude; // WARNING: unused (void)magnitude; // WARNING: unused
xx.x = xx.x * (2.0l * M_PI / (get(AC_dsx) * (get(AC_ny_max) - get(AC_ny_min)))); xx.x = xx.x*(2.0*M_PI/(get(AC_dsx)*get(AC_nx)));
xx.y = xx.y * (2.0l * M_PI / (get(AC_dsy) * (get(AC_ny_max) - get(AC_ny_min)))); xx.y = xx.y*(2.0*M_PI/(get(AC_dsy)*get(AC_ny)));
xx.z = xx.z * (2.0l * M_PI / (get(AC_dsz) * (get(AC_ny_max) - get(AC_ny_min)))); xx.z = xx.z*(2.0*M_PI/(get(AC_dsz)*get(AC_nz)));
ModelScalar cosl_phi = cosl(phi); ModelScalar cos_phi = cosl(phi);
ModelScalar sinl_phi = sinl(phi); ModelScalar sin_phi = sinl(phi);
ModelScalar cosl_k_dox_x = cosl(dot(k_force, xx)); ModelScalar cos_k_dox_x = cosl(dot(k_force, xx));
ModelScalar sinl_k_dox_x = sinl(dot(k_force, xx)); ModelScalar sin_k_dox_x = sinl(dot(k_force, xx));
// Phase affect only the x-component // Phase affect only the x-component
// ModelScalar real_comp = cosl_k_dox_x; //Scalar real_comp = cos_k_dox_x;
// ModelScalar imag_comp = sinl_k_dox_x; //Scalar imag_comp = sin_k_dox_x;
ModelScalar real_comp_phase = cosl_k_dox_x * cosl_phi - sinl_k_dox_x * sinl_phi; ModelScalar real_comp_phase = cos_k_dox_x*cos_phi - sin_k_dox_x*sin_phi;
ModelScalar imag_comp_phase = cosl_k_dox_x * sinl_phi + sinl_k_dox_x * cosl_phi; ModelScalar imag_comp_phase = cos_k_dox_x*sin_phi + sin_k_dox_x*cos_phi;
ModelVector force = (ModelVector){ff_re.x * real_comp_phase - ff_im.x * imag_comp_phase,
ff_re.y * real_comp_phase - ff_im.y * imag_comp_phase, ModelVector force = (ModelVector){ ff_re.x*real_comp_phase - ff_im.x*imag_comp_phase,
ff_re.z * real_comp_phase - ff_im.z * imag_comp_phase}; ff_re.y*real_comp_phase - ff_im.y*imag_comp_phase,
ff_re.z*real_comp_phase - ff_im.z*imag_comp_phase};
return force; return force;
} }