Warwick-Plasma / epoch

Particle-in-cell code for plasma physics simulations
https://epochpic.github.io
GNU General Public License v3.0
184 stars 58 forks source link

How does absorption work? #522

Open ktangtar opened 1 year ago

ktangtar commented 1 year ago

Hi,

I have a simulation with no particles. I define two lasers at the x boundary. After they are both in, epoch gives me these two values: Total_Field_Energy_in_Simulation_J = 26.5354 Absorption/Total_Laser_Energy_Injected_J = 16.7829

I expect these to be the same, or the first one less than the second one. I would think "Total_Field_Energy_in_Simulation_J" may be straightforward to calculate, proportional to the square of the fields.

For "Total_Laser_Energy_Injected_J ", is it only looking at the energy of one laser? Or maybe calculating the energy of each laser separately and then adding them up? I think the correct way would be to add up both the fields of the laser at the boundary, and calculate the energy from that?

Regards, Kavin Tangtartharakul


begin:constant 
####-----laser beam parameters
  laser_lambda    = 0.8 * micron
  laser_intensity = 3.5e21 # a0~40
  laser_omega     = 2.0 * pi * c / laser_lambda
  T_wave          = laser_lambda / c
  xf              = 5.0*micron

  w0_1              = 2.9 * micron # 
  rayl_1            = pi*w0_1^2/laser_lambda
  sG_1              = xf/rayl_1
  wb_1              = w0_1*sqrt(1+sG_1^2) 
  n_1               = 2.13

  w0_2              = 8.7 * micron # 
  rayl_2            = pi*w0_2^2/laser_lambda
  sG_2              = xf/rayl_2
  wb_2              = w0_2*sqrt(1+sG_2^2) 
  n_2               = 10.05
####-----plasma parameters
  n_crit          = critical(laser_omega)
  den_max         = 0.0*n_crit
####-----output parameters
  T_ave           = 2.0*T_wave
  T_snap          = 50.0*femto
  T_snap_freq     = 1/2.73*T_wave
  T_snap_infreq   = 100.0*femto
  T_snap_dist     = 40.0*femto
end:constant

begin:control
  nx = 105*15
  ny = 30*15
  nz = 30*15

  x_min = -5.0*micron
  x_max =  100.0*micron
  y_min = -15.0*micron
  y_max =  15.0*micron
  z_min = -15.0*micron
  z_max =  15.0*micron

####----------------------duration of simulation
  t_end = 400.1*femto

  dlb_threshold = 0.4
#  restart_snapshot = roll0001.sdf
  use_random_seed = T
end:control

begin:boundaries
  bc_x_min = simple_laser
  bc_x_max = open
  bc_y_min = open 
  bc_y_max = open
  bc_z_min = open
  bc_z_max = open
end:boundaries

begin:laser
  boundary = x_min 
  lambda = laser_lambda
  intensity_w_cm2 = laser_intensity / (1+sG_1^2) / n_1

  t_profile = gauss(time,50.0e-15,25.0e-15/sqrt(2*loge(2))) #intensity FWHM 25.0fs, intensity w_0 as input
  profile = exp(-(sqrt(y*y+z*z)/wb_1)^2)

  phase = sG_1*(sqrt(y*y+z*z)/w0_1)^2/(1+sG_1^2)
  pol = 0
end:laser

begin:laser
  boundary = x_min 
  lambda = laser_lambda
  intensity_w_cm2 = laser_intensity / (1+sG_2^2) / n_2

  t_profile = gauss(time,50.0e-15,25.0e-15/sqrt(2*loge(2))) #intensity FWHM 25.0fs, intensity w_0 as input
  profile = exp(-(sqrt(y*y+z*z)/wb_2)^2)

  phase = sG_2*(sqrt(y*y+z*z)/w0_2)^2/(1+sG_2^2)
  pol = 0
end:laser

begin:output
  name = ey_fields_freq
  file_prefix = ey_fields_freq
  dump_first = T
  time_start = 0.0*femto
  dt_snapshot = T_snap/4
  grid = always + single 
  ey = always + single + snapshot
end:output

begin:output
name = absorp
  file_prefix = absorp
  dt_snapshot= T_snap_freq
# dt_average = T_ave
  dump_first = T
  absorption = always + single 
  total_energy_sum = always + single
end:output
Status-Mirror commented 1 year ago

For Absorption/Total_Laser_Energy_Injected_J, each MPI rank keeps a variable called laser_inject_local. This takes the laser intensity at each boundary cell (as specified in the laser block), multiplies it by cell area and simulation timestep, and estimates the injected energy. The laser_inject_local sums energies from all cells and time-steps in the local rank. This is done in laser.F90, and an MPI_ALLREDUCE command in io/diagnostics.F90 produces the final global output.

For Total_Field_Energy_in_Simulation_J, the calculation is performed in the calc_total_energy_sum subroutine of io/calc_df.F90. This calculates the sums of $E^2$ and $B^2$, and combines them to get the total EM energy density.

I was playing around with a more simple input deck (written below), and I found the following. If we inject two uniform lasers into a 1D simulation as in the below example, we get:

Here, the Absorption/Total_Laser_Energy_Injected_J sums each laser contribution independently. The Total_Field_Energy_in_Simulation_J considers the super-posed fields. As another example, consider if we lauched the two lasers in anti-phase. Each laser may be firing into the simulation window, but they would interfere destructively and there would be no field energy in the simulation window. This is why the two numbers may differ.

The Absorption/Fraction_of_Laser_Energy_Absorbed would appear to be quite useless here. In addition to being energetically impossible, it's also misleading as the simulation domain is too large for laser energy to escape.

Perhaps it would be a good idea to switch the Absorption/Total_Laser_Energy_Injected_J to calculate the injected Poynting flux, or maybe the documentation should be altered to make the distinction between the diagnostics more clear. I'll leave the issue open while we figure out what we want to do.

Cheers, Stuart

begin:control
  nx = 800
  x_min = 0
  x_max = 40.0e-6
  t_end = 100.0e-15
end:control

begin:boundaries
  bc_x_min = simple_laser
  bc_x_max = open
end:boundaries

begin:laser
  boundary = x_min 
  lambda = 1.0e-6
  intensity_w_cm2 = 1.0e19
end:laser

begin:laser
  boundary = x_min 
  lambda = 1.0e-6
  intensity_w_cm2 = 1.0e19
end:laser

begin:output
  dt_snapshot= 1.0e-15
  absorption = always + single 
  total_energy_sum = always + single
  ey = always
end:output
ktangtar commented 1 year ago

Thanks for the very thorough explanation!