From 9400838859f325637a097e1aadd7e134c6565987 Mon Sep 17 00:00:00 2001 From: Cameron Yang Date: Fri, 8 Nov 2019 09:33:05 -0800 Subject: Added Maxwell-Juttner getMomentum function to IM.H and PI.cpp. --- Source/Initialization/PlasmaInjector.cpp | 22 ++++++++++++++++++++++ 1 file changed, 22 insertions(+) (limited to 'Source/Initialization/PlasmaInjector.cpp') diff --git a/Source/Initialization/PlasmaInjector.cpp b/Source/Initialization/PlasmaInjector.cpp index af04c6b31..417959ab3 100644 --- a/Source/Initialization/PlasmaInjector.cpp +++ b/Source/Initialization/PlasmaInjector.cpp @@ -269,6 +269,28 @@ void PlasmaInjector::parseMomentum (ParmParse& pp) // Construct InjectorMomentum with InjectorMomentumGaussian. inj_mom.reset(new InjectorMomentum((InjectorMomentumGaussian*)nullptr, ux_m, uy_m, uz_m, ux_th, uy_th, uz_th)); + } else if (mom_dist_s == "maxwell_juttner"){ + Real beta = 0.; + Real theta = 10.; + int dir = 0; + std::string direction = "x"; + pp.query("beta", beta); + pp.query("theta", theta); + pp.query("direction", direction); + if(direction == "x"){ + dir = 0; + }else if(direction == "y"){ + dir = 1; + }else if(direction == "z"){ + dir = 2; + }else{ + std::stringstream stringstream; + stringstream<<"Direction "<< direction <<" is not recognzied. Please enter x, y, or z."; + direction = stringstream.str(); + amrex::Abort(direction.c_str()); + } + // Construct InjectorMomentum with InjectorMomentumJuttner. + inj_mom.reset(new InjectorMomentum((InjectorMomentumJuttner*)nullptr, theta, beta, dir)); } else if (mom_dist_s == "radial_expansion") { Real u_over_r = 0.; pp.query("u_over_r", u_over_r); -- cgit v1.2.3 From 29fbc3c68ae74607696c044978e2966f81e5b7ad Mon Sep 17 00:00:00 2001 From: Cameron Yang Date: Fri, 8 Nov 2019 10:45:17 -0800 Subject: Corrected mistakes and made formatting changes. --- Source/Initialization/InjectorMomentum.H | 8 ++++---- Source/Initialization/PlasmaInjector.cpp | 10 +++++----- 2 files changed, 9 insertions(+), 9 deletions(-) (limited to 'Source/Initialization/PlasmaInjector.cpp') diff --git a/Source/Initialization/InjectorMomentum.H b/Source/Initialization/InjectorMomentum.H index 9b92dea74..5784559f6 100644 --- a/Source/Initialization/InjectorMomentum.H +++ b/Source/Initialization/InjectorMomentum.H @@ -54,10 +54,10 @@ struct InjectorMomentumJuttner InjectorMomentumJuttner(amrex::Real theta, amrex::Real beta, int dir) noexcept : t(theta), b(beta), G(1/std::sqrt(1-pow(beta,2))), d(dir) {} - + AMREX_GPU_HOST_DEVICE amrex::XDim3 - getMomentum (amrex::Real, amrex::Real, amrex::Real) const noexcept + getMomentum (amrex::Real x, amrex::Real y, amrex::Real z) const noexcept { amrex::Real x1, x2; amrex::Real g, ux; @@ -65,7 +65,7 @@ struct InjectorMomentumJuttner x1 = 0.; g = 0.; u[d] = 0.; - while( u[d]-g <= x1) + while(u[d]-g <= x1) { u[d] = -t*std::log(amrex::Random()*amrex::Random()*amrex::Random()); g = std::sqrt(1+std::pow(u[d],2)); @@ -80,7 +80,7 @@ struct InjectorMomentumJuttner } ux = G*(ux+g*b); x2 = amrex::Random(); - u[(d+1)%3] = 2*u[d]*std::sqrt(x1*(1-x1))*std::cos(2*M_PI*x2); + u[(d+1)%3] = 2*u[d]*std::sqrt(x1*(1-x1))*std::sin(2*M_PI*x2); u[(d+2)%3] = 2*u[d]*std::sqrt(x1*(1-x1))*std::cos(2*M_PI*x2); u[d] = ux; return amrex::XDim3 {u[0],u[1],u[2]}; diff --git a/Source/Initialization/PlasmaInjector.cpp b/Source/Initialization/PlasmaInjector.cpp index 417959ab3..de9b6279c 100644 --- a/Source/Initialization/PlasmaInjector.cpp +++ b/Source/Initialization/PlasmaInjector.cpp @@ -277,15 +277,15 @@ void PlasmaInjector::parseMomentum (ParmParse& pp) pp.query("beta", beta); pp.query("theta", theta); pp.query("direction", direction); - if(direction == "x"){ + if(direction == "x" || direction == "X"){ dir = 0; - }else if(direction == "y"){ + } else if (direction == "y" || direction == "Y"){ dir = 1; - }else if(direction == "z"){ + } else if (direction == "z" || direction == "Z"){ dir = 2; - }else{ + } else{ std::stringstream stringstream; - stringstream<<"Direction "<< direction <<" is not recognzied. Please enter x, y, or z."; + stringstream << "Direction " << direction << " is not recognzied. Please enter x, y, or z."; direction = stringstream.str(); amrex::Abort(direction.c_str()); } -- cgit v1.2.3 From 0c592033a2dab4c465da853052b1140faa4930ce Mon Sep 17 00:00:00 2001 From: Cameron Yang Date: Fri, 8 Nov 2019 15:35:04 -0800 Subject: Untabify --- Source/Initialization/InjectorMomentum.H | 64 ++++++++++++++++---------------- Source/Initialization/PlasmaInjector.cpp | 42 ++++++++++----------- 2 files changed, 53 insertions(+), 53 deletions(-) (limited to 'Source/Initialization/PlasmaInjector.cpp') diff --git a/Source/Initialization/InjectorMomentum.H b/Source/Initialization/InjectorMomentum.H index 404aacaa8..108ec6913 100644 --- a/Source/Initialization/InjectorMomentum.H +++ b/Source/Initialization/InjectorMomentum.H @@ -52,36 +52,36 @@ private: struct InjectorMomentumJuttner { InjectorMomentumJuttner(amrex::Real theta, amrex::Real beta, int dir) noexcept - : t(theta), b(beta), d(dir) - {} + : t(theta), b(beta), d(dir) + {} AMREX_GPU_HOST_DEVICE amrex::XDim3 getMomentum (amrex::Real x, amrex::Real y, amrex::Real z) const noexcept { - amrex::Real x1, x2, g; - amrex::Real u [3]; - x1 = 0.; - g = 0.; - u[d] = 0.; - while(u[d]-g <= x1) - { - u[d] = -t*std::log(amrex::Random()*amrex::Random()*amrex::Random()); - g = std::sqrt(1+std::pow(u[d],2)); - x1 = t*std::log(amrex::Random()); - } - x1 = amrex::Random(); - x2 = amrex::Random(); - u[(d+1)%3] = 2*u[d]*std::sqrt(x1*(1-x1))*std::sin(2*M_PI*x2); - u[(d+2)%3] = 2*u[d]*std::sqrt(x1*(1-x1))*std::cos(2*M_PI*x2); - u[d] = u[d]*(2*x1-1); - x1 = amrex::Random(); - if(-b*u[d]/g>x1) - { - u[d] = -u[d]; - } - u[d] = 1/std::sqrt(1-pow(b,2))*(u[d]+g*b); - return amrex::XDim3 {u[0],u[1],u[2]}; + amrex::Real x1, x2, g; + amrex::Real u [3]; + x1 = 0.; + g = 0.; + u[d] = 0.; + while(u[d]-g <= x1) + { + u[d] = -t*std::log(amrex::Random()*amrex::Random()*amrex::Random()); + g = std::sqrt(1+std::pow(u[d],2)); + x1 = t*std::log(amrex::Random()); + } + x1 = amrex::Random(); + x2 = amrex::Random(); + u[(d+1)%3] = 2*u[d]*std::sqrt(x1*(1-x1))*std::sin(2*M_PI*x2); + u[(d+2)%3] = 2*u[d]*std::sqrt(x1*(1-x1))*std::cos(2*M_PI*x2); + u[d] = u[d]*(2*x1-1); + x1 = amrex::Random(); + if(-b*u[d]/g>x1) + { + u[d] = -u[d]; + } + u[d] = 1/std::sqrt(1-pow(b,2))*(u[d]+g*b); + return amrex::XDim3 {u[0],u[1],u[2]}; } private: @@ -169,9 +169,9 @@ struct InjectorMomentum // This constructor stores a InjectorMomentumJuttner in union object. InjectorMomentum (InjectorMomentumJuttner* t, - amrex::Real theta, amrex::Real beta, int dir) - : type(Type::juttner), - object(t, theta, beta, dir) + amrex::Real theta, amrex::Real beta, int dir) + : type(Type::juttner), + object(t, theta, beta, dir) { } // This constructor stores a InjectorMomentumCustom in union object. @@ -215,7 +215,7 @@ struct InjectorMomentum { return object.gaussian.getMomentum(x,y,z); } - case Type::juttner: + case Type::juttner: { return object.juttner.getMomentum(x,y,z); } @@ -258,8 +258,8 @@ private: amrex::Real a_uz_m, amrex::Real a_ux_th, amrex::Real a_uy_th, amrex::Real a_uz_th) noexcept : gaussian(a_ux_m,a_uy_m,a_uz_m,a_ux_th,a_uy_th,a_uz_th) {} - Object (InjectorMomentumJuttner*, - amrex::Real t, amrex::Real b, int dir) noexcept + Object (InjectorMomentumJuttner*, + amrex::Real t, amrex::Real b, int dir) noexcept : juttner(t,b,dir) {} Object (InjectorMomentumRadialExpansion*, amrex::Real u_over_r) noexcept @@ -272,7 +272,7 @@ private: InjectorMomentumConstant constant; InjectorMomentumCustom custom; InjectorMomentumGaussian gaussian; - InjectorMomentumJuttner juttner; + InjectorMomentumJuttner juttner; InjectorMomentumRadialExpansion radial_expansion; InjectorMomentumParser parser; }; diff --git a/Source/Initialization/PlasmaInjector.cpp b/Source/Initialization/PlasmaInjector.cpp index de9b6279c..32711cd31 100644 --- a/Source/Initialization/PlasmaInjector.cpp +++ b/Source/Initialization/PlasmaInjector.cpp @@ -270,27 +270,27 @@ void PlasmaInjector::parseMomentum (ParmParse& pp) inj_mom.reset(new InjectorMomentum((InjectorMomentumGaussian*)nullptr, ux_m, uy_m, uz_m, ux_th, uy_th, uz_th)); } else if (mom_dist_s == "maxwell_juttner"){ - Real beta = 0.; - Real theta = 10.; - int dir = 0; - std::string direction = "x"; - pp.query("beta", beta); - pp.query("theta", theta); - pp.query("direction", direction); - if(direction == "x" || direction == "X"){ - dir = 0; - } else if (direction == "y" || direction == "Y"){ - dir = 1; - } else if (direction == "z" || direction == "Z"){ - dir = 2; - } else{ - std::stringstream stringstream; - stringstream << "Direction " << direction << " is not recognzied. Please enter x, y, or z."; - direction = stringstream.str(); - amrex::Abort(direction.c_str()); - } - // Construct InjectorMomentum with InjectorMomentumJuttner. - inj_mom.reset(new InjectorMomentum((InjectorMomentumJuttner*)nullptr, theta, beta, dir)); + Real beta = 0.; + Real theta = 10.; + int dir = 0; + std::string direction = "x"; + pp.query("beta", beta); + pp.query("theta", theta); + pp.query("direction", direction); + if(direction == "x" || direction == "X"){ + dir = 0; + } else if (direction == "y" || direction == "Y"){ + dir = 1; + } else if (direction == "z" || direction == "Z"){ + dir = 2; + } else{ + std::stringstream stringstream; + stringstream << "Direction " << direction << " is not recognzied. Please enter x, y, or z."; + direction = stringstream.str(); + amrex::Abort(direction.c_str()); + } + // Construct InjectorMomentum with InjectorMomentumJuttner. + inj_mom.reset(new InjectorMomentum((InjectorMomentumJuttner*)nullptr, theta, beta, dir)); } else if (mom_dist_s == "radial_expansion") { Real u_over_r = 0.; pp.query("u_over_r", u_over_r); -- cgit v1.2.3