Possible deviation when interpolating the fields felt by the particles
XiangyanAn opened this issue · 4 comments
I just found that the EPOCH codes are one-cell shifted when calculating the nearest boundary cell of the particles and interpolating the fields felt by the particles. Such a deviation would result in a numerical increase of a conserved quantity when the electrons interact with a circular-polarized laser. This is confirmed by reading the source code and checking the fields of the particles after setting the fields as E_{x,y,z}=B_{x,y,z}=x. The numerical increase is fixed by adding one to the variables cell_x2, cell_y2, and cell_z2.
I wrote a detailed description in the attached document.
Hey @XiangyanAn,
Thanks for the detailed write-up, this is definitely something worth checking up.
Your pdf file mentions a few different input decks - one with a circularly polarised beam, and another with the linear fields. Could you also send these if you still have them? I can have a go at recreating them if not.
Cheers,
Stuart
I'm still going through your report, but from a partial read I can see an EPOCH bug in the fields
block. I followed your method of setting
In this simulation, we have 10 cells ranging from 0m to 10m, such that we have:
- Cell edges (m) =
[0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10]
(11 edges) - Cell centres (m) =
[0.5, 1.5, 2.5, 3.5, 4.5, 5.5, 6.5, 7.5, 8.5, 9.5]
(10 cells)
Using the developer documentation stagger guide, or Fig 3b in your report, we expect the 10 cells of Ey
, Ez
and Bx
to match the cell centre numbers, and this is what we observe. You acknowledge this agreement in figures 2b, 2c, 2e.
For staggered field-points, the Yee grid suggests a stagger in the positive direction. If cell with index ix
has edges at 3m and 4m, then Ex(ix)
should be evaluated on the high
- Expected
$x$ stagger fields =[1, 2, 3, 4, 5, 6, 7, 8, 9, 10]
We observe:
- Ex =
[0, 1, 2, 3, 4, 5, 6, 7, 8, 0]
- By =
[0, 1, 2, 3, 4, 5, 6, 7, 8, 8.6780]
- Bz =
[0, 1, 2, 3, 4, 5, 6, 7, 8, 8.6780]
Even ignoring the strange value on the high fields
block.
With my understanding of the problem, this issue should only arise if you're using the fields
block. When setting up your circular polarisation in equations 1 and 2, did you use the fields
block, or the laser
block?
Cheers,
Stuart
begin:control
nx = 10
t_end = 1.0e-15
x_min = 0
x_max = 10
stdout_frequency = 100
end:control
begin:boundaries
bc_x_min = simple_laser
bc_x_max = open
end:boundaries
begin:fields
ex = x
ey = x
ez = x
bx = x
by = x
bz = x
end:fields
begin:output
dt_snapshot = t_end
ex = always
ey = always
ez = always
bx = always
by = always
bz = always
end:output
I'll also take a look at this tomorrow.
I'm still going through your report, but from a partial read I can see an EPOCH bug in the
fields
block. I followed your method of setting E and B to x, but in a 1D simulation. Input deck attached at the end.In this simulation, we have 10 cells ranging from 0m to 10m, such that we have:
- Cell edges (m) =
[0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10]
(11 edges)- Cell centres (m) =
[0.5, 1.5, 2.5, 3.5, 4.5, 5.5, 6.5, 7.5, 8.5, 9.5]
(10 cells)Using the developer documentation stagger guide, or Fig 3b in your report, we expect the 10 cells of
Ey
,Ez
andBx
to match the cell centre numbers, and this is what we observe. You acknowledge this agreement in figures 2b, 2c, 2e.For staggered field-points, the Yee grid suggests a stagger in the positive direction. If cell with index
ix
has edges at 3m and 4m, thenEx(ix)
should be evaluated on the high x cell-edge, or at 4m in this case. Hence, I expect the x-stagger fields to be at positions:
- Expected x stagger fields =
[1, 2, 3, 4, 5, 6, 7, 8, 9, 10]
We observe:
- Ex =
[0, 1, 2, 3, 4, 5, 6, 7, 8, 0]
- By =
[0, 1, 2, 3, 4, 5, 6, 7, 8, 8.6780]
- Bz =
[0, 1, 2, 3, 4, 5, 6, 7, 8, 8.6780]
Even ignoring the strange value on the high x simulation edge, we can see that these fields are being treated as if they evaluated on the low x cell-edge, which is inconsistent with how the rest of the code treats grid stagger. This is why you were able to "correct" the problem by adjusting the interpolation. I would instead argue that it is not the interpolation which is broken, but the
fields
block.With my understanding of the problem, this issue should only arise if you're using the
fields
block. When setting up your circular polarisation in equations 1 and 2, did you use thefields
block, or thelaser
block?Cheers, Stuart
begin:control nx = 10 t_end = 1.0e-15 x_min = 0 x_max = 10 stdout_frequency = 100 end:control begin:boundaries bc_x_min = simple_laser bc_x_max = open end:boundaries begin:fields ex = x ey = x ez = x bx = x by = x bz = x end:fields begin:output dt_snapshot = t_end ex = always ey = always ez = always bx = always by = always bz = always end:output
Thanks for your reminder. Since I wanted a simplest problem in physics to make a test, I did use the fields
block rather than the laser
block to make the plane wave.
I have read the codes of 'fields.F90' and the developer documentation again and made some more tests. I now understand that the Yee grid suggests a stagger in the positive direction for staggered field-points. This means that the index of Ex
, By
, and Bz
are 1 less than that of the corresponding xb
grid, while I took them as the same in the previously written document.
I think you are right. It is not the interpolation which is broken, but the fields
block.