gwm17/catima
1
0
Fork 0
mirror of https://github.com/gwm17/catima.git synced 2024-11-22 18:28:51 -05:00
Code Issues Projects Releases Wiki Activity

Fixed bug for units of dedx in gwm_integrators

Browse Source
This commit is contained in:
Gordon McCann 2022-06-09 09:36:48 -04:00
parent af4da593a8
commit 629690a6f9
1 changed files with 9 additions and 7 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

14
catima/gwm_integrators.cpp
View File

@ -12,6 +12,7 @@ namespace catima {
double x_step = 0.25*mat.thickness(); //g/cm^2 double x_step = 0.25*mat.thickness(); //g/cm^2
double x_traversed = 0.0; double x_traversed = 0.0;
double e_step = dedx(proj, mat, c)*x_step; double e_step = dedx(proj, mat, c)*x_step;
double A_recip = 1.0/proj.A;
while(true) while(true)
{ {
@ -19,12 +20,12 @@ namespace catima {
{ {
++depth; ++depth;
x_step *= 0.5; x_step *= 0.5;
e_step = dedx(proj, mat, c)*x_step; e_step = dedx(proj, mat, c)*x_step*A_recip;
} }
else if(x_step + x_traversed >= mat.thickness()) else if(x_step + x_traversed >= mat.thickness())
{ {
x_step = mat.thickness() - x_traversed; x_step = mat.thickness() - x_traversed;
e_step = dedx(proj, mat, c)*x_step; e_step = dedx(proj, mat, c)*x_step*A_recip;
e_final -= e_step; e_final -= e_step;
proj.T = e_final; proj.T = e_final;
return (e_in - e_final)*proj.A; return (e_in - e_final)*proj.A;
@ -35,7 +36,7 @@ namespace catima {
} }
else else
{ {
e_step = dedx(proj, mat, c)*x_step; e_step = dedx(proj, mat, c)*x_step*A_recip;
e_final -= e_step; e_final -= e_step;
proj.T = e_final; proj.T = e_final;
x_traversed += x_step; x_traversed += x_step;
@ -53,6 +54,7 @@ namespace catima {
double x_step = 0.25*mat.thickness(); //g/cm^2 double x_step = 0.25*mat.thickness(); //g/cm^2
double x_traversed = 0.0; double x_traversed = 0.0;
double e_step = dedx(proj, mat, c)*x_step; double e_step = dedx(proj, mat, c)*x_step;
double A_recip = 1.0/proj.A;
while(true) while(true)
{ {
@ -60,13 +62,13 @@ namespace catima {
{ {
++depth; ++depth;
x_step *= 0.5; x_step *= 0.5;
e_step = dedx(proj, mat, c)*x_step; e_step = dedx(proj, mat, c)*x_step*A_recip;
} }
else if(x_step + x_traversed >= mat.thickness()) else if(x_step + x_traversed >= mat.thickness())
{ {
x_step = mat.thickness() - x_traversed; x_step = mat.thickness() - x_traversed;
e_step = dedx(proj, mat, c)*x_step; e_step = dedx(proj, mat, c)*x_step;
e_initial += e_step; e_initial += e_step*A_recip;
proj.T = e_initial; proj.T = e_initial;
return (e_initial - e_out)*proj.A; return (e_initial - e_out)*proj.A;
} }
@ -77,7 +79,7 @@ namespace catima {
else else
{ {
e_step = dedx(proj, mat, c)*x_step; e_step = dedx(proj, mat, c)*x_step;
e_initial += e_step; e_initial += e_step*A_recip;
proj.T = e_initial; proj.T = e_initial;
x_traversed += x_step; x_traversed += x_step;
} }