diff options
Diffstat (limited to 'Example')
| -rw-r--r-- | Example/Langmuir/langmuir.py | 78 | ||||
| -rw-r--r-- | Example/Langmuir/langmuir_PICMI.py | 36 | ||||
| -rw-r--r-- | Example/laser_acceleration/laser_acceleration_PICMI.py | 199 | ||||
| -rw-r--r-- | Example/plasma_acceleration/plasma_acceleration_PICMI.py | 54 |
4 files changed, 217 insertions, 150 deletions
diff --git a/Example/Langmuir/langmuir.py b/Example/Langmuir/langmuir.py new file mode 100644 index 000000000..da9e3d8ad --- /dev/null +++ b/Example/Langmuir/langmuir.py @@ -0,0 +1,78 @@ +import numpy as np +from pywarpx import * + +nx = 64 +ny = 64 +nz = 64 + +xmin = -20.e-6 +ymin = -20.e-6 +zmin = -20.e-6 +xmax = +20.e-6 +ymax = +20.e-6 +zmax = +20.e-6 + +dx = (xmax - xmin)/nx +dy = (ymax - ymin)/ny +dz = (zmax - zmin)/nz + +# Maximum number of time steps +max_step = 40 + +# number of grid points +amr.n_cell = [nx, ny, nz] + +# Maximum allowable size of each subdomain in the problem domain; +# this is used to decompose the domain for parallel calculations. +amr.max_grid_size = 32 + +# Maximum level in hierarchy (for now must be 0, i.e., one level in total) +amr.max_level = 0 + +amr.plot_int = 1 # How often to write plotfiles. "<= 0" means no plotfiles. + +# Geometry +geometry.coord_sys = 0 # 0: Cartesian +geometry.is_periodic = [1, 1, 1] # Is periodic? +geometry.prob_lo = [xmin, ymin, zmin] # physical domain +geometry.prob_hi = [xmax, ymax, zmax] + +# Verbosity +warpx.verbose = 1 +warpx.do_moving_window = 0 +warpx.moving_window_dir = 'x' +warpx.moving_window_v = 0.0 # in units of the speed of light + +# Algorithms +algo.current_deposition = 3 +algo.charge_deposition = 0 +algo.field_gathering = 0 +algo.particle_pusher = 0 + +# CFL +warpx.cfl = 1.0 + +particles.nspecies = 1 +particles.species_names = 'electrons' + +electrons.charge = '-q_e' +electrons.mass = 'm_e' +electrons.injection_style = "NUniformPerCell" +electrons.num_particles_per_cell_each_dim = [2, 2, 2] + +electrons.xmin = xmin +electrons.xmax = 0.e-6 +electrons.ymin = ymin +electrons.ymax = ymax +electrons.zmin = zmin +electrons.zmax = zmax + +electrons.profile = 'constant' +electrons.density = 1.e25 # number of electrons per m^3 + +electrons.momentum_distribution_type = "constant" +electrons.ux = 0.01 # ux = gamma*beta_x + +# --- Initialize the simulation +warpx.write_inputs('inputs_from_python', max_step=max_step) + diff --git a/Example/Langmuir/langmuir_PICMI.py b/Example/Langmuir/langmuir_PICMI.py new file mode 100644 index 000000000..dd2156157 --- /dev/null +++ b/Example/Langmuir/langmuir_PICMI.py @@ -0,0 +1,36 @@ +import numpy as np +from pywarpx import PICMI + +nx = 64 +ny = 64 +nz = 64 + +xmin = -20.e-6 +ymin = -20.e-6 +zmin = -20.e-6 +xmax = +20.e-6 +ymax = +20.e-6 +zmax = +20.e-6 + +grid = PICMI.Grid(nx=nx, ny=ny, nz=nz, + xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, zmin=zmin, zmax=zmax, + bcxmin='periodic', bcxmax='periodic', bcymin='periodic', bcymax='periodic', bczmin='periodic', bczmax='periodic', + moving_window_velocity=[0., 0., 0.], + max_grid_size=32, max_level=0, coord_sys=0) + +solver = PICMI.EM_solver(current_deposition_algo = 3, + charge_deposition_algo = 0, + field_gathering_algo = 0, + particle_pusher_algo = 0) + +electrons = PICMI.Species(type='electron', name='electrons') + +plasma = PICMI.Plasma(species=electrons, density=1.e25, xmax=0., vxmean=0.1, number_per_cell_each_dim=[2,2,2]) + +sim = PICMI.Simulation(verbose = 1, + timestep_over_cfl = 1.0, + max_step = 40, + plot_int = 1) + +sim.write_inputs(inputs_name='inputs_from_PICMI') + diff --git a/Example/laser_acceleration/laser_acceleration_PICMI.py b/Example/laser_acceleration/laser_acceleration_PICMI.py index 761ff7cfe..307d3e7c8 100644 --- a/Example/laser_acceleration/laser_acceleration_PICMI.py +++ b/Example/laser_acceleration/laser_acceleration_PICMI.py @@ -2,167 +2,66 @@ import numpy as np from pywarpx import PICMI #from warp import PICMI -from mpi4py import MPI -comm = MPI.COMM_WORLD +nx = 64 +ny = 64 +nz = 480 -def get_parallel_indices(Np, rank, size): - ''' +xmin = -30.e-6 +xmax = +30.e-6 +ymin = -30.e-6 +ymax = +30.e-6 +zmin = -56.e-6 +zmax = +12.e-6 - This decomposes the arrays of particle attributes into subarrays - for parallel processing. - - ''' - Navg = Np / size - Nleft = Np - Navg * size - if (rank < Nleft): - lo = rank*(Navg + 1) - hi = lo + Navg + 1 - else: - lo = rank * Navg + Nleft - hi = lo + Navg - return lo, hi - - -def load_plasma(ncells, domain_min, domain_max, injected_plasma_density, injected_plasma_ppc, plasma_min, plasma_max): - ''' - - Sets initial conditions for the plasma acceleration setup. - - ''' - - comm.Barrier() - - dx = (domain_max - domain_min) / ncells - - nplasma_cells = ((plasma_max - plasma_min)/dx + 0.5).astype('l') + 1 - iplasma_min = ((plasma_min - domain_min)/dx).astype('l')*dx + domain_min - - # Species 0 - the plasma - plasma_weight = injected_plasma_density * dx[0] * dx[1] * dx[2] / injected_plasma_ppc - - # Fill the entire domain with particles. Only a subset of these points - # will be selected for each species - Z, Y, X, P = np.mgrid[0:nplasma_cells[2], 0:nplasma_cells[1], 0:nplasma_cells[0], 0:injected_plasma_ppc] - Z = Z.flatten() - Y = Y.flatten() - X = X.flatten() - P = P.flatten() - - particle_shift = (0.5 + P) / injected_plasma_ppc - - xp = (X + particle_shift)*dx[0] + iplasma_min[0] - yp = (Y + particle_shift)*dx[1] + iplasma_min[1] - zp = (Z + particle_shift)*dx[2] + iplasma_min[2] - - # now do the plasma species - plasma_locs = np.logical_and(xp >= plasma_min[0], xp < plasma_max[0]) - plasma_locs = np.logical_and(plasma_locs, yp >= plasma_min[1]) - plasma_locs = np.logical_and(plasma_locs, zp >= plasma_min[2]) - plasma_locs = np.logical_and(plasma_locs, yp < plasma_max[1]) - plasma_locs = np.logical_and(plasma_locs, zp < plasma_max[2]) - - plasma_xp = xp[plasma_locs] - plasma_yp = yp[plasma_locs] - plasma_zp = zp[plasma_locs] - - plasma_uxp = np.zeros_like(plasma_xp) - plasma_uyp = np.zeros_like(plasma_xp) - plasma_uzp = np.zeros_like(plasma_xp) - - # --- and weights - Np = plasma_xp.shape[0] - plasma_wp = np.full(Np, plasma_weight) - - lo, hi = get_parallel_indices(Np, comm.rank, comm.size) - - electrons.add_particles(x=plasma_xp[lo:hi], y=plasma_yp[lo:hi], z=plasma_zp[lo:hi], - ux=plasma_uxp[lo:hi], uy=plasma_uyp[lo:hi], uz=plasma_uzp[lo:hi], - w=plasma_wp[lo:hi], unique_particles=True) - - comm.Barrier() - -def inject_new_plasma(ncells, domain_min, domain_max, num_shift, direction, injected_plasma_density, injected_plasma_ppc, plasma_min, plasma_max): - ''' - - This injects fresh plasma into the domain after the moving window has been upated. - - ''' - - comm.Barrier() - - dx = (domain_max - domain_min) / ncells - - plasma_min_inject = plasma_min.copy() - plasma_max_inject = plasma_max.copy() - if (num_shift > 0): - plasma_min_inject[direction] = domain_max[direction] - plasma_max_inject[direction] = domain_max[direction] + num_shift*dx[direction] - else: - plasma_min_inject[direction] = domain_min[direction] - num_shift*dx[direction] - plasma_max_inject[direction] = domain_min[direction] - - load_plasma(ncells, domain_min, domain_max, injected_plasma_density, injected_plasma_ppc, plasma_min_inject, plasma_max_inject) - -ncells = np.array([64, 64, 480])/2 -domain_min = np.array([-30.e-6, -30.e-6, -56.e-6]) -domain_max = np.array([ 30.e-6, 30.e-6, 12.e-6]) -dx = (domain_max - domain_min) / ncells - -moving_window_velocity = np.array([0., 0., PICMI.clight]) +moving_window_velocity = [0., 0., PICMI.c] injected_plasma_density = 1.e23 -injected_plasma_ppc = 4 -plasma_min = np.array([-20.e-6, -20.e-6, 0.0e-6]) -plasma_max = np.array([ 20.e-6, 20.e-6, domain_max[2]]) +number_per_cell_each_dim = [2, 2, 2] +plasma_min = [-20.e-6, -20.e-6, 0.0e-6] +plasma_max = [ 20.e-6, 20.e-6, zmax] -grid = PICMI.Grid(nx=ncells[0], ny=ncells[1], nz=ncells[2], - xmin=domain_min[0], xmax=domain_max[0], ymin=domain_min[1], ymax=domain_max[1], zmin=domain_min[2], zmax=domain_max[2], +grid = PICMI.Grid(nx=nx, ny=ny, nz=nz, xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, zmin=zmin, zmax=zmax, bcxmin='periodic', bcxmax='periodic', bcymin='periodic', bcymax='periodic', bczmin='open', bczmax='open', moving_window_velocity = moving_window_velocity, max_grid_size=32, max_level=0, coord_sys=0) -solver = PICMI.EM_solver(current_deposition = 3, - charge_deposition = 0, - field_gathering = 0, - particle_pusher = 0) +solver = PICMI.EM_solver(current_deposition_algo = 3, + charge_deposition_algo = 0, + field_gathering_algo = 0, + particle_pusher_algo = 0) -laser = PICMI.Gaussian_laser(antenna_z0 = 9.e-6, # This point is on the laser plane - antenna_zvec = 1., # The plane normal direction - pol_angle = PICMI.pi/2., # The main polarization vector - E0 = 16.e12, # Maximum amplitude of the laser field (in V/m) +t_peak = 30.e-15 # The time at which the laser reaches its peak at the antenna (in seconds) +focal_distance = 100.e-6 # Focal distance from the antenna (in meters) +antenna_z0 = 9.e-6 # This point is on the laser plane +laser = PICMI.Gaussian_laser(wavelength = 0.8e-6, # The wavelength of the laser (in meters) waist = 5.e-6, # The waist of the laser (in meters) duration = 15.e-15, # The duration of the laser (in seconds) - t_peak = 30.e-15, # The time at which the laser reaches its peak (in seconds) - focal_position = 109.e-6, # Focal position (m) - wavelength = 0.8e-6, # The wavelength of the laser (in meters) - em_solver = solver) - -electrons = PICMI.Species(type=PICMI.Electron, name='electrons') - -# Maximum number of time steps -max_step = 1000 - -sim = PICMI.Simulation(plot_int = 200, + pol_angle = np.pi/2., # The main polarization vector + focal_position = focal_distance + antenna_z0, # Focal position (m) + E0 = 16.e12, # Maximum amplitude of the laser field (in V/m) + z0 = antenna_z0 - PICMI.c*t_peak, # Position of the laser centroid in Z at time 0 + ) + +laser = PICMI.Laser_antenna(laser = laser, + antenna_z0 = antenna_z0, # This point is on the laser plane + antenna_zvec = 1., # The plane normal direction + ) + +electrons = PICMI.Species(type='electron', name='electrons') + +plasma = PICMI.Plasma(species = electrons, + density = injected_plasma_density, + xmin = plasma_min[0], + xmax = plasma_max[0], + ymin = plasma_min[1], + ymax = plasma_max[1], + zmin = plasma_min[2], + zmax = plasma_max[2], + number_per_cell_each_dim = number_per_cell_each_dim) + +sim = PICMI.Simulation(plot_int = 100, verbose = 1, - cfl = 1.0) - -load_plasma(ncells, domain_min, domain_max, injected_plasma_density, injected_plasma_ppc, plasma_min, plasma_max) - -direction = np.argmax(abs(moving_window_velocity)) -old_x = grid.getmins()[direction] -new_x = old_x -for i in range(1, max_step + 1): - - # check whether the moving window has updated - num_shift = int( 0.5 + (new_x - old_x) / dx[direction]) - if (num_shift != 0): - inject_new_plasma(ncells, domain_min, domain_max, num_shift, direction, injected_plasma_density, injected_plasma_ppc, plasma_min, plasma_max) - domain_min[direction] += num_shift*dx[direction] - domain_max[direction] += num_shift*dx[direction] - - sim.step(1) - - old_x = new_x - new_x = grid.getmins()[direction] + cfl = 1.0, + max_step = 1000) -sim.finalize() +sim.write_inputs(inputs_name = 'inputs_from_PICMI') diff --git a/Example/plasma_acceleration/plasma_acceleration_PICMI.py b/Example/plasma_acceleration/plasma_acceleration_PICMI.py new file mode 100644 index 000000000..f3bf6a574 --- /dev/null +++ b/Example/plasma_acceleration/plasma_acceleration_PICMI.py @@ -0,0 +1,54 @@ +import numpy as np +from pywarpx import PICMI +#from warp import PICMI + +nx = 64 +ny = 64 +nz = 64 + +xmin = -200.e-6 +xmax = +200.e-6 +ymin = -200.e-6 +ymax = +200.e-6 +zmin = -200.e-6 +zmax = +200.e-6 + +moving_window_velocity = [0., 0., PICMI.c] + +number_per_cell_each_dim = [2, 2, 1] + +grid = PICMI.Grid(nx=nx, ny=ny, nz=nz, xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax, zmin=zmin, zmax=zmax, + bcxmin='periodic', bcxmax='periodic', bcymin='periodic', bcymax='periodic', bczmin='open', bczmax='open', + moving_window_velocity = moving_window_velocity, + max_grid_size=32, max_level=0, coord_sys=0) + +solver = PICMI.EM_solver(current_deposition_algo = 3, + charge_deposition_algo = 0, + field_gathering_algo = 0, + particle_pusher_algo = 0) + +beam = PICMI.Species(type='electron', name='beam') +plasma = PICMI.Species(type='electron', name='plasma') + +beam_distribution = PICMI.Plasma(species = beam, + density = 1.e23, + xmin = -20.e-6, xmax = +20.e-6, + ymin = -20.e-6, ymax = +20.e-6, + zmin = -150.e-6, zmax = -100.e-6, + vzmean = 1.e9, + number_per_cell_each_dim = number_per_cell_each_dim) + +plasma_distribution = PICMI.Plasma(species = plasma, + density = 1.e22, + xmin = -200.e-6, xmax = +200.e-6, + ymin = -200.e-6, ymax = +200.e-6, + zmin = 0., + number_per_cell_each_dim = number_per_cell_each_dim, + fill_in = True) + +sim = PICMI.Simulation(plot_int = 2, + verbose = 1, + cfl = 1.0, + max_step = 1000) + +sim.write_inputs(inputs_name = 'inputs_from_PICMI') |