Laser group velocity is greater than light speed with field_order=4

[XiangyanAn]

Hello,
I just found that with field_order=4 in the control block, the group velocity of a laser pulse in the simulation could be greater than the speed of light.
I'm using the release epoch-4.19.0, which is unzipped as epoch-4.18.0-2022-12-14, and taking 1d simulation. The input.deck is :

begin:constant
lambda0 = 800e-9
k0 = 2 * pi / lambda0
omega0 = c * k0
tau0 = 2 * pi / omega0

a0 = 1.4
tau_laser = 8.0 * lambda0 / c
delay = tau_laser * 3.0

end:constant

begin:control
x_min = -30.0 * lambda0
x_max = 30.0 * lambda0

nx = 60 * 30

t_end = delay + (0 - x_min) / c + 4000.0 * tau0

field_order = 4
end:control

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

begin:window
move_window = T
window_v_x = c
window_start_time = (x_max - x_min) / c
bc_x_min_after_move = simple_outflow
bc_x_max_after_move = simple_outflow
end:window

begin:laser
boundary = x_min
amp = a0 * me * omega0 * c / qe
lambda = lambda0
polarisation_angle = 0.0
phase = 0.0
t_profile = gauss(time, delay, tau_laser)
end:laser

begin:output
dt_snapshot = t_end / 50
full_dump_every = 0

grid = full
ex = full
ey = full
ez = full
bx = full
by = full
bz = full
end:output

After simulation, the laser pulse electric fields of the files 0001.sdf and 0050.sdf are, respectively,

Since the window is moving at the light speed, we can see that the pulse is faster than light in vacuum.
Moreover, I also calculated the group velocity according to the centroid position of the laser pulse: $x_c = \int x E_y^2\mathrm d x/ \int E_y^2\mathrm d x$. The difference of the centroid position and time of 0001.sdf and 0050.sdf also gives, $(x_{c,50}-x_{c,1})/[c(t_{50}-t_{1})]=1.0036>1$.

I'm not sure whether it is a feature of the fourth field order, or a bug.
Thans for any help.

[Status-Mirror]

Hey @XiangyanAn,

Are you sure there's a problem here? When I calculate your t_end, I get around 10.8ps. In 10.8ps, the start of the laser can travel around 3.24mm - which is beyond your window. This is what we observe, isn't it?

Cheers,
Stuart

[XiangyanAn]

Hi @Status-Mirror ,

We can notice that there is a delay of $t_d=24\tau_0$ for the laser pulse peak and it is injected at the $x_\text{min}$ boundary. Therefore, the pulse peak will arrive at the position $x=0$ at the time $t=t_d+\left|x_\text{min}\right|/c$. Then the laser pulse will travel for $4000\tau_0=10.674 \text{ps}$ and it will arrive at 3.2 mm at the simulation end, which should be still in the window.

However, as we observed, it is not. I also calculated the group velocity quantitatively in the previous comment, which is also greater than the light speed.

Thanks,
Xiangyan An

[Status-Mirror]

Hey Xiangyan An,

Sorry, I misread your initial input deck - I was assuming your x_min was 0. I can confirm you do have a laser which has travelled too far. I was playing around with some values, and I think the issue is due to resolution. When I double the cell-count, I find the laser travels less distance. The second-order field solver seems less sensitive to the cell-size, so this may be better for vacuum simulations.

The 4th-order field solver calculates the gradients for Maxwell's equation using more cells. In a noisy plasma environment, this can smooth out the spatial field gradients and reduce noice. It seems this method is doing more harm than good in your vacuum example though.

Hope this helps,
Stuart

[TomGoffrey]

Numerical dispersion can result in group (and phase) velocities greater than c (or indeed less than). I would suggest reading "A systematic approach to numerical dispersion in Maxwell solvers" by Blinne et al (2018) as a starting point to understand more.

[XiangyanAn]

Hi @TomGoffrey,

Thanks for your help.
I used to think that numerical dispersion would only result in group velocities less than $c$. Otherwise, there might be causality problems. I have also observed group velocity less than $c$ with field_order = 2. That's why I'm surprised by the group velocity greater than $c$ with field_order = 4.
I now have a deeper understanding of numerical dispersion. Thanks again for your reminder.

Also thank @Status-Mirror for the further detailed check. This problem can indeed be alleviated by using more cells.

[amarsarari]

Hey Xiangyan An,

Sorry, I misread your initial input deck - I was assuming your x_min was 0. I can confirm you do have a laser which has travelled too far. I was playing around with some values, and I think the issue is due to resolution. When I double the cell-count, I find the laser travels less distance. The second-order field solver seems less sensitive to the cell-size, so this may be better for vacuum simulations.

The 4th-order field solver calculates the gradients for Maxwell's equation using more cells. In a noisy plasma environment, this can smooth out the spatial field gradients and reduce noice. It seems this method is doing more harm than good in your vacuum example though.

Hope this helps, Stuart

[amarsarari]

Hi sir can you please provide me matlab code to see this electric field Ey vs x

[Status-Mirror]

Hey @amarsarari,

First I open MATLAB, find the EPOCH directory in the internal file explorer, and right click to select "Add to path -> Selected folders and subfolders". This adds all the SDF files from the EPOCH directory to the path, so we can open them.

Once I've done this, I run this script:

% plot_Ey.m
% Plots the Ey field from a specific SDF file of a 1D EPOCH simulation

%% User input

file_name = "0050.sdf";

%% Extract data

% Pull variables from SDF file
data = GetDataSDF(file_name);
Ey = data.Electric_Field.Ey.data;
x_edges = data.Grid.Grid.x;

% Format variables
x_centres = 0.5*(x_edges(2:end) + x_edges(1:end-1));

%% Create plot

plot(x_centres * 1.0e6, Ey, 'LineWidth', 2);
box on;
grid on;
xlabel('x [\mum]');
ylabel('Ey [V/m]');
ax = gca;
ax.FontSize = 16;

For the multiple lines shown in my figure, you can repeat this process and stack the lines by running hold on, and re-running this script with a new filename. A legend can be added using the legend("Line 1 name", "Line 2 name", "Line 3 name") command.

Hope this helps,
Stuart