What boundary conditions did PICNIC use when solving the Poisson equation

[qixinimp]

Dear Didier,

we are studying the effect of boundary on beam. According to the Instructions of Tracewin, the inner potential of the beam are calculated by solving Poisson equation. A particle-in-cell approach based on numerical-calculation of the interaction between cubes (PICNIC) is used. The PICNIC mesh extends to ±3.5*X, ±3.5*Y and ±3.5*Z. X, Y and Z are the x, y and z beam rms-sizes, respectively. Electric field at positions outside the mesh is calculated as that of a gaussian beam with the same rms-sizes.

Here is my question:

1) When solving Poisson equation, what is the boundary conditions used in PICNIC? The instructions of Tracewin don’t mention that. I looked up the references and still couldn't find it. Could you give me some instruction or references?

2) The instructions says electric field at positions outside the mesh is calculated as that of a gaussian beam with the same rms-sizes. In a real beam line, the x, y and z beam rms-sizes may be different. In this case, how does the software choose the appropriate rms? And how much does it affect the results?

Thanks
Xin Qi

[npichoff]

Dear Xin Qi,

The PICNIC subroutine is considering open conditions (no (image) charge outside the mesh).

When the particle is outside the mesh, the field is calculated as if the beam had a gaussian distribution with the same charge and the same 3 RMS sizes of the beam (but all three RMS sizes could be different). This field can be calculated with a simple numerical integration. This calculation being longer than for particles into mesh, it is better if it concerns only a few particles.

We have shown that, when far enough (let say about 4 sigmas), the field shape is not depending too much from the inner particle detailled distribution.
Our recommandation is to adjust the size of the mesh according to the number of particles and the number of lattices:
- For very fast calculation (fast design...), 10k parcticles with a mesh size of 3.5 sigmas and a number of lattice of 2*(7*7*7) is generally Ok.
- For fast calculation (optim design...), 100k particles with a mesh size of 4 sigmas and a number of lattice of 2*(10*10*10) is generally OK.
- For long calculation (finalisation...), 1M particles with a mesh size of 5 sigmas and a number of lattice of 2*(15*15*15) is generally OK.
- For very long calculation (halo...), 10M or particles with a mesh size of 6 sigmas and a number of lattice of 2*(20*20*20) is generally OK.

Of course, if your beam is much bigger(in average) in one direction than in other (for example longitudinal), it is better to use a higher number of lattices in this direction (the more square the lattice, the better the model).

As each configuration can be specific, my last recommandation would be to test different parameters to check the sensitivity/stability (convergence) of your calculation with respect to these parameters.

Best regards.

Nicolas.

[qixinimp]

Dear Nicolas,
Thank you very much!
I still have another personal question. When using Tracewin, we all marvel at the speed of the software. As far as I known, the Poisson equation with open boundary condition should be integrated to infinity by, Green function or some other methods. These methods will consume a lot of computing time, So what method does PICNIC used? Since the software runs so fast, Is that FFT ? Can you kindly give me some tips or literature？
We look forward to your reply.

Best wishes!

Xin

[npichoff]

Dear Xin Qi,

Here is the way PICNIC3D works:
- At each time step, the beam RMS sizes (Xrms, Yrms and Zrms) are calculated (particles too far from the beam c.o.g. are excluded from the calculation).
- The mesh sizes is adjusted to +/-xW*Wrms, (W for X, Y or Z; xW from 3.5 to 20 as suggested before).
- The mesh is cut in 2*Nx*2*Ny*2*Nz lattices.
- Considering a generic node of this mesh, the contributions, to the field, of a lattice as if it was uniformly charged by one particle, are calculated and tabulated in a interaction table, according to its relative distance to the node. The way it is done depends on the distance between the node and the mesh :
o if the node is attached to the mesh, an numerical integration is done considering a varying charge distribution in the 8 closest lattices (in order to avoid a charge discontinuity on the node),
o if the lattices are relativily close, an simple numerical integration is done considering a uniform charge distribution in the lattices,
o if the lattices are far, the charge is simply considered at the center of the lattice.
- The number of particles are counted in the lattice (either considering the full charge in the lattice where the particle stands or spread in the 8 closest lattices.
- the field at each node is then calculated using the charge distribution in each lattices and its contribution given by the interaction table.
- the field on each particle is interpolated from the field on each node.
- the field outside the mesh is calculated by a simple numerical integration as if the beam was gaussian in all direction.

The calculation time is then :
T = k1*Np + k2*(2Nx*2Ny*2Nz*(2*Nx+1)*(2*Ny+1)*(2*Nz+1)
k1 : Proportionnality coefficient given by the particle spread in the mesh and field calculation for each particle,
k2 : Proportionnality coefficient given by the field calculation on each node.

Then, in therory, the calculation time should be longer than FFT or a multigrid algorithm (TOUTATIS) where :
T = k1*Np + k3*(2Nx*2Ny*2Nz)*log((2Nx*2Ny*2Nz))

Nevertheless, the secret in PICNIC is:
- The tabulation, which allows to reduce significantly the coefficient k2 with respect to k3,
- The reduced mesh size, which allows to reduce significantly the number of meshs for the same resolution.

A static image charge could be easily added in cylindrical pipes (maybe we already did that?) without increasing significantly the calculation time. We already implemented the image from neighboor bunchs when the beam is flat (Zrms << X/Y rms) or is close to continuous (de-bunching).

You know now almost all of PICNIC secrets .

Best regards.

Nicolas.

[qixinimp]

Dear Nicolas,

I got it！Thank you very much! It is so kind of you!

Your reply has solved my long-standing question. Before that, I did not known the boundary conditions used in Tracewin, so I tried to solve the open boundary problem with Professor JI Qiang's algorium（Computer Physics Communications 219 (2017) 255–260）. However, the algorium speed can not be improved. According to your explanation, I realize that the algorithm adopts cylindrical coordinates. The number of meshes cannot be reduced due to axisymmetry. It seems that variable meshes in Cartesian coordinate system used in PICNIC is the best way to greatly reduce the mesh number. That's really a genius idea！

I still have a problem. In Professor JI Qiang's algorithm, the free boundary condition in cylindrical coordinates is solved by Green's function integration, so how do you realize the free boundary condition in Cartesian coordinates？

Best regards！

Xin Qi

[npichoff]

Dear Xin Qi,

PICNIC considers only the contribution of charges standing in the mesh (Gauss theorem : E = 1/(4*pi*eps0 d^2)).
No other charge (outside the mesh) is considered, except the charges from neighbooring bunches in the conditions described in my preceeding post.
Then PICNIC don't need to set boundary conditions.

I imagine that, if I had to do so (for example in FTT or multigrid algorithm), I would have used the boundary conditions given by a gaussian beam with the same RMS size . It is probably better than to put the same voltage or the same field all along the boundary.
Maybe you could try in your algorithm ?

Best regards.

Nicolas.

[qixinimp]

Dear Nicolas,

Thank you! I get it now! I've been wrong. I have always tried to understand PICNIC as a numerical algorithm that solve Poisson equation with some kind of boundary conditions. Now I realize Picnic did not use a similar scheme.

Now I almost fully understand Picnic. But I still have some questions, so I drew a picture according to my understanding.
1 In my understanding, are the mesh, lattice and node right, as shown in the figure?
2 For a attached node as shown in the figure, an numerical integration is done considering a varying charge distribution in the 8 closest lattices to avoid the Integral singularity. But how does PICNIC chose the varying charge distribution in the 8 closest lattices?
3 For a close node, there is no Integral singularity, so a simple integration is enough，Is my understanding right？
4 For a node far alway, lattice can be viewed as a point charge, but how does PICNIC distinguish a node that is near or far?

In my point of view, PICNIC is really a genius scheme. Under this scheme, the amount of calculation has been greatly reduced.

I was planning to solve Poisson's equation by spectral method, which is also very fast. According to your suggestion, the boundary conditions may be an analytical solution of a gaussian beam with the same RMS size, or semi-analytical solution of a gaussian beam with different RMS sizes. When I finish it, the speed as well as accuracy will be compared with Tracewin. I'll contact you then.

Sincerely yours,

Xin

[npichoff]

Dear Xin,

1- Yes the mesh is rectangular along/orthogonal to the beam average direction,

2- In the closest lattices, PICNIC considers that :
- the total charge in each lattice is given by the number of particles in it,
- the charge density in each lattice is linearly varying to the same density (the sum of the charges in the 8 closest lattices divided by the volume of these lattices) on the central node.
Normally, you should not need to do that if you plane to use spectral algorithm.

3-Yes it is right.

4-If I remember well, PICNIC considers that a lattice is far enough if its average distance is 10 greater than the size of the lattice.

I wish you good pleasure and success in your development. I'm impatient to here from your results and benchmarks.

Have a nice day.

Best regards.

Nicolas.

[warson]

Dear Nicolas,
I'm also interested in the detail of PICNIC algorithm and thanks for your answer. However, I still have some questions to ask you.

1. What's the expression of the linearly varying charge density? Is it in the form of ax+by+cz, where x,y,z is the coordinate axes, and if so, how did you determine the coefficients?

2. If the field calculation for the closest lattices is dependent on the charge distribution, how did you make the interaction table for the closest lattices since the table should be made before calculating the charge distribution in my opinion?

I'm looking forward to your reply.

Best regards,
Warson

[npichoff]

Dear Warson,

1- if (dx, dy, dz) is the signed distance with respect to the closest cell (in fact the cell where the particle stands). In fraction of mesh size (Xc, Yc and Zc), it means that -0.5 < dx/Xc, dy/Yc, dz/Zc < 0.5). The particle charge is deposited in the eight closest cells, with weights which is the mulitplication of:
- |dx/Xc| or 1-|dx/Xc|,
- by |dy/Yc| or 1-|dy/Yc|,
- by |dz/Zc| or 1-|dz/Zc|,
For example, for the closest cell (i,j,k), the weight is : (1-|dx/Xc|)*(1-|dy/Yc|)*(1-|dz/Zc|) and for the furthest cell (i+sign(dx), j+sign(dy), k+sign(dz)), the weight is : |dx/Xc|*|dy/Yc|*|dz/Zc|.
It is as if the particle where represented by a uniformily charged cubic rectangle with lengths Xc, Yc, Zc and we use as weight the faction of the particle volume in each cell.

2- the influence of the charge distribution in each lattice (considered as uniform) is tabulated, according to its normalised distance with respect to a node. Only one eigth of the space is tabulated (the other parts can be deduced by a change of sign). The closer the mesh from the node, the more precise the calculation is made. For neighboor lattices, the calculation considers the gradient between the 8 lattices in order to avoid infinite field because of density discontinuity from lattice to lattice.

Now, you have all the secrets of PICNIC...
I hope it is clear enough.

Best regards.

Nicolas.