diff options
Diffstat (limited to 'Examples/Tests/initial_distribution')
| -rwxr-xr-x | Examples/Tests/initial_distribution/analysis_distribution.py | 77 | ||||
| -rw-r--r-- | Examples/Tests/initial_distribution/inputs | 88 |
2 files changed, 163 insertions, 2 deletions
diff --git a/Examples/Tests/initial_distribution/analysis_distribution.py b/Examples/Tests/initial_distribution/analysis_distribution.py index 3da4cf423..d43ced2bb 100755 --- a/Examples/Tests/initial_distribution/analysis_distribution.py +++ b/Examples/Tests/initial_distribution/analysis_distribution.py @@ -12,6 +12,8 @@ # 3 denotes maxwell-juttner distribution. # 4 denotes gaussian position distribution. # 5 denotes maxwell-juttner distribution w/ spatially varying temperature +# 6 denotes maxwell-boltzmann distribution w/ constant velocity +# 7 denotes maxwell-boltzmann distribution w/ spatially-varying velocity # The distribution is obtained through reduced diagnostic ParticleHistogram. import numpy as np @@ -180,5 +182,80 @@ print('Maxwell-Juttner parser temperature difference:', f5_error) assert(f5_error < tolerance) +#============================================== +# maxwell-boltzmann with constant bulk velocity +#============================================== + +# load data +bin_value_g, bin_data_g = read_reduced_diags_histogram("h6.txt")[2:] +bin_value_uy, bin_data_uy = read_reduced_diags_histogram("h6uy.txt")[2:] + +# Expected values for beta and u = beta*gamma +beta_const = 0.2 +g_const = 1. / np.sqrt(1. - beta_const * beta_const) +uy_const = beta_const * g_const +g_bin_size = 0.004 +g_bin_min = 1. +uy_bin_size = 0.04 +uy_bin_min = -1. +V = 8.0 # volume in m^3 +n = 1.0e21 # number density in 1/m^3 + +f_g = np.zeros_like(bin_value_g) +i_g = int(np.floor((g_const - g_bin_min) / g_bin_size)) +f_g[i_g] = n * V +f_peak = np.amax(f_g) + +f_uy = np.zeros_like(bin_value_uy) +i_uy = int(np.floor((-uy_const - uy_bin_min) / uy_bin_size)) +f_uy[i_uy] = n * V + +f6_error = np.sum(np.abs(f_g - bin_data_g) + np.abs(f_uy - bin_data_uy)) \ + / bin_value_g.size / f_peak + +print('Maxwell-Boltzmann constant velocity difference:', f6_error) + +assert(f6_error < tolerance) + +#============================================ +# maxwell-boltzmann with parser bulk velocity +#============================================ + +# load data +bin_value_g, bin_data_g = read_reduced_diags_histogram("h7.txt")[2:] +bin_value_uy, bin_data_uy_neg = read_reduced_diags_histogram("h7uy_neg.txt")[2:] +bin_data_uy_pos = read_reduced_diags_histogram("h7uy_pos.txt")[3] + +# Expected values for beta and u = beta*gamma +beta_const = 0.2 +g_const = 1. / np.sqrt(1. - beta_const * beta_const) +uy_const = beta_const * g_const +g_bin_size = 0.004 +g_bin_min = 1. +uy_bin_size = 0.04 +uy_bin_min = -1. +V = 8.0 # volume in m^3 +n = 1.0e21 # number density in 1/m^3 + +f_g = np.zeros_like(bin_value_g) +i_g = int(np.floor((g_const - g_bin_min) / g_bin_size)) +f_g[i_g] = n * V +f_peak = np.amax(f_g) + +f_uy_neg = np.zeros_like(bin_value_uy) +i_uy_neg = int(np.floor((uy_const - uy_bin_min) / uy_bin_size)) +f_uy_neg[i_uy_neg] = n * V / 2. + +f_uy_pos = np.zeros_like(bin_value_uy) +i_uy_pos = int(np.floor((-uy_const - uy_bin_min) / uy_bin_size)) +f_uy_pos[i_uy_pos] = n * V / 2. + +f7_error = np.sum(np.abs(f_g - bin_data_g) + np.abs(f_uy_pos - bin_data_uy_pos) \ + + np.abs(f_uy_neg - bin_data_uy_neg)) / bin_value_g.size / f_peak + +print('Maxwell-Boltzmann parser velocity difference:', f7_error) + +assert(f7_error < tolerance) + test_name = filename[:-9] # Could also be os.path.split(os.getcwd())[1] checksumAPI.evaluate_checksum(test_name, filename) diff --git a/Examples/Tests/initial_distribution/inputs b/Examples/Tests/initial_distribution/inputs index 0b47249df..d301c5751 100644 --- a/Examples/Tests/initial_distribution/inputs +++ b/Examples/Tests/initial_distribution/inputs @@ -29,7 +29,7 @@ algo.particle_shape = 1 ################################# ############ PLASMA ############# ################################# -particles.species_names = gaussian maxwell_boltzmann maxwell_juttner beam maxwell_juttner_parser +particles.species_names = gaussian maxwell_boltzmann maxwell_juttner beam maxwell_juttner_parser velocity_constant velocity_parser particles.rigid_injected_species = beam gaussian.charge = -q_e @@ -93,6 +93,32 @@ maxwell_juttner_parser.momentum_distribution_type = "maxwell_juttner" maxwell_juttner_parser.theta_distribution_type = "parser" maxwell_juttner_parser.theta_function(x,y,z) = "1.0 + heaviside(x,0)" +velocity_constant.charge = -q_e +velocity_constant.mass = m_e +velocity_constant.injection_style = "NRandomPerCell" +velocity_constant.num_particles_per_cell = 1000 +velocity_constant.profile = constant +velocity_constant.density = 1.0e21 +velocity_constant.momentum_distribution_type = "maxwell_boltzmann" +velocity_constant.theta_distribution_type = "constant" +velocity_constant.theta = 1e-9 +velocity_constant.beta_distribution_type = "constant" +velocity_constant.beta = 0.2 +velocity_constant.bulk_vel_dir = -y + +velocity_parser.charge = -q_e +velocity_parser.mass = m_e +velocity_parser.injection_style = "NRandomPerCell" +velocity_parser.num_particles_per_cell = 1000 +velocity_parser.profile = constant +velocity_parser.density = 1.0e21 +velocity_parser.momentum_distribution_type = "maxwell_boltzmann" +velocity_parser.theta_distribution_type = "constant" +velocity_parser.theta = 1e-9 +velocity_parser.beta_distribution_type = "parser" +velocity_parser.beta_function(x,y,z) = "-0.2 + 0.4 * heaviside(z,0)" +velocity_parser.bulk_vel_dir = -y + ################################# ########## DIAGNOSTIC ########### ################################# @@ -101,7 +127,9 @@ maxwell_juttner_parser.theta_function(x,y,z) = "1.0 + heaviside(x,0)" # 3 for maxwell-juttner with constant temperature # 4 for beam # 5 for maxwell-juttner with parser function temperature -warpx.reduced_diags_names = h1x h1y h1z h2x h2y h2z h3 h3_filtered h4x h4y h4z bmmntr h5_neg h5_pos +# 6 for maxwell-boltzmann with constant velocity +# 7 for maxwell-boltzmann with parser velocity +warpx.reduced_diags_names = h1x h1y h1z h2x h2y h2z h3 h3_filtered h4x h4y h4z bmmntr h5_neg h5_pos h6 h6uy h7 h7uy_pos h7uy_neg h1x.type = ParticleHistogram h1x.intervals = 1 @@ -223,6 +251,62 @@ h5_pos.bin_max = 12.0 h5_pos.histogram_function(t,x,y,z,ux,uy,uz) = "sqrt(1.0 + ux*ux + uy*uy + uz*uz)" h5_pos.filter_function(t,x,y,z,ux,uy,uz) = "x > 0.0" +h5_pos.type = ParticleHistogram +h5_pos.intervals = 1 +h5_pos.path = "./" +h5_pos.species = maxwell_juttner_parser +h5_pos.bin_number = 50 +h5_pos.bin_min = 1.0 +h5_pos.bin_max = 12.0 +h5_pos.histogram_function(t,x,y,z,ux,uy,uz) = "sqrt(1.0 + ux*ux + uy*uy + uz*uz)" +h5_pos.filter_function(t,x,y,z,ux,uy,uz) = "x > 0.0" + +h6.type = ParticleHistogram +h6.intervals = 1 +h6.path = "./" +h6.species = velocity_constant +h6.bin_number = 50 +h6.bin_min = 1.0 +h6.bin_max = 1.2 +h6.histogram_function(t,x,y,z,ux,uy,uz) = "sqrt(1.0 + ux*ux + uy*uy + uz*uz)" + +h6uy.type = ParticleHistogram +h6uy.intervals = 1 +h6uy.path = "./" +h6uy.species = velocity_constant +h6uy.bin_number = 50 +h6uy.bin_min = -1.0 +h6uy.bin_max = 1.0 +h6uy.histogram_function(t,x,y,z,ux,uy,uz) = "uy" + +h7.type = ParticleHistogram +h7.intervals = 1 +h7.path = "./" +h7.species = velocity_parser +h7.bin_number = 50 +h7.bin_min = 1.0 +h7.bin_max = 1.2 +h7.histogram_function(t,x,y,z,ux,uy,uz) = "sqrt(1.0 + ux*ux + uy*uy + uz*uz)" + +h7uy_pos.type = ParticleHistogram +h7uy_pos.intervals = 1 +h7uy_pos.path = "./" +h7uy_pos.species = velocity_parser +h7uy_pos.bin_number = 50 +h7uy_pos.bin_min = -1.0 +h7uy_pos.bin_max = 1.0 +h7uy_pos.histogram_function(t,x,y,z,ux,uy,uz) = "uy" +h7uy_pos.filter_function(t,x,y,z,ux,uy,uz) = "z > 0.0" + +h7uy_neg.type = ParticleHistogram +h7uy_neg.intervals = 1 +h7uy_neg.path = "./" +h7uy_neg.species = velocity_parser +h7uy_neg.bin_number = 50 +h7uy_neg.bin_min = -1.0 +h7uy_neg.bin_max = 1.0 +h7uy_neg.histogram_function(t,x,y,z,ux,uy,uz) = "uy" +h7uy_neg.filter_function(t,x,y,z,ux,uy,uz) = "z < 0.0" # our little beam monitor bmmntr.type = BeamRelevant |