ComputationalRadiationPhysics / picongpu

Performance-Portable Particle-in-Cell Simulations for the Exascale Era :sparkles:
https://picongpu.readthedocs.io
Other
710 stars 218 forks source link

Absorber and symmetry #3148

Closed cbontoiu closed 4 years ago

cbontoiu commented 4 years ago

1) I try to understand what the absorber means namely ABSORBER_CELLS and ABSORBER_STRENGTH. I couldn't find information in the manual, apart from guidance on the usage. My assumption is that absorber makes sense only when a laser or background field is present and enables interaction between this and the matter as declared in density.param. Since in my case matter and interactions are present for the whole simulation domain, a declaration such as:

constexpr uint32_t ABSORBER_CELLS[3][2] = {
    {0, 256},  /*x direction [negative,positive]*/
    {0, 256},  /*y direction [negative,positive]*/
    {0, 1280} /*z direction [negative,positive]*/
};

In the snapshots attached I used TBG_gridSize="384 384 32" and TBG_periodic="--periodic 1 0 1" so there is a mismatch. The simulation domain along y (laser direction) is 384 cells long while the interaction region (absrber cells) is only 256. However, I see that the interaction happens really within the first half of the region. For the second half I see only the laser pulse passing unaffected.

Screenshot

Also, how is the ABSORBER_STRENGTH defined since I always see values like 1e-03? Is it fine to work with:

constexpr float_X ABSORBER_STRENGTH[3][2] = {
    {1.0e-3, 1.0e-3}, /*x direction [negative,positive]*/
    {1.0e-3, 1.0e-3}, /*y direction [negative,positive]*/
    {1.0e-3, 1.0e-3}  /*z direction [negative,positive]*/
};

though absorber cells may not be defined in both sides for all three directions?

2) I once tried to align the material density with respect to the x-vertical axis. With respect to the image below, I have used negative coordinates for x (upper part) trying to get a circular shape but when I ran with something like TBG_gridSize="N N M" TBG_periodic="--periodic 1 0 1" only the part hashed with red appeared.

circle

Is it necessary to work with a centred target? If yes, how can I do it?

Many thanks

sbastrakov commented 4 years ago

About the absorber. The parameters in question are for exponential damping absorber, the only absorber available in version 0.4.3. Please note that the dev branch has a more advanced convolutional PML absorber which does not require nearly as thick absorbing layer. However, this is, unfortunately, not available in 0.4.3, so the rest of my message concerns the exponential damping absorber.

I believe the TBG_periodic="--periodic 1 0 1" is processed first, and so the code will only apply absorber along the Y axis in this example, even if you set non-zero thickness for X or Z. (However, it is still required to specify the parameters for all 6 boundaries to keep the structure intact, as you correctly did).

ABSORBER_CELLS defines thickness, in number of cells, of the absorbing layer along each of the 6 boundaries of your simulation domain, the order is given by the comments. Please note that in case it exceeds the domain size, the domain will be extended to fit the absorber. Geometrically, the absorber starts at the outer boundary of the domain and goes inside the domain. For your example, 256 cells were requested for the Y positive absorber, and your domain has 384 cells in Y, so it is mostly covered with the absorber. This is probably suboptimal. Generally, it is desirable to have no particles in the absorber area (since their behavior would be not physical in the absorber area). For this one would virtually split the domain into absorber area (near the outer border) and no-absorber area (internal part of the domain) and only put particles in the internal part.

Having ABSORBER_CELLS[...][...] == 0 at some axis and direction means there is no absorber there, so the corresponding value of ABSORBER_STRENGTH is not used, but you still have to specify it to keep the structure intact. The values of ABSORBER_STRENGTH define the decay rate inside the absorbing layer. The absorber introduces two kinds of reflections: from the boundary between the absorber and the internal area, and from the outer boundary of the absorber. Increasing the absorber strength will increase the first kind of reflections, and decrease the second, and vice versa. Since one would like to minimize the combined effect of such reflections, there is a tradeoff and the strength value should be not too small and not too large. I don't think there is a good recipe of what value to use, 1.0e-3 is reasonable in our experience. You might try to adjust it in case there are visible reflections.

sbastrakov commented 4 years ago

Regarding the second part about symmetry. Could you attach your set of .param files which results in this asymmetry? I understand you might have done it in another issue, and sorry for asking again for potentially the exact same data, however, this way makes this issue more consistent and guarantees we look at the most up-to-date version.

sbastrakov commented 4 years ago

Just based on your comments about negative x I am confused. The PIConGPU global coordinate system is set up so that the domain always starts at (0, 0, 0) and goes to the positive direction, with the size being a product of cell size and number of cells. So in terms of this coordinate system there can be no negative coordinates. Of course, inside your density, and other setups you can use any coordinate system and in the end transform it to ours. However, operating in our coordinates, a target in the center of the domain has coordinates equal to cell size * half number of cells in the domain

cbontoiu commented 4 years ago

It is clear now. I was confused by one of the earlier replies to my "Periodicity" inquiry saying that "the asymmetry is coming from the non centered target." So I thought that the code can handle negative coordinates but I didn't see how to handle this with the only positive TBG_gridSize="m m n", so I turned my attention to the absorber which allows for the negative sides. Here is a snapshot of what I got. For me the y-axis is horizontal here and the laser comes from the left to the right. e_png_yx_0 5_000000 In my density file I had: // vertical axis ............ constexpr double X_O = 0; // middle constexpr double X_A = X_O - 1(tube_out_diam + tubes_gap)cossin45; // upper-left constexpr double X_B = X_O - 1(tube_out_diam + tubes_gap)cossin45; // upper-right constexpr double X_C = X_O + 1(tube_out_diam + tubes_gap)cossin45; // lower-right constexpr double X_D = X_O + 1(tube_out_diam + tubes_gap)cossin45; // lower-left // horizontal axis .......... constexpr double Y_O = 1(tube_out_diam + tubes_gap)cossin45; // middle constexpr double Y_A = 0; // upper-left constexpr double Y_B = Y_O + (tube_out_diam + tubes_gap)cossin45; // upper-right constexpr double Y_C = Y_O + (tube_out_diam + tubes_gap)cossin45; // lower-right constexpr double Y_D = 0; // lower-left and for the cfg file: TBG_gridSize="384 384 32" TBG_periodic="--periodic 1 0 1" and as you can see the upper left and upper right tubes are not rendered which is just expected, but I couldn't be sure without your clarifications.

Thanks. I don't think the param files are necessary in this case.

sbastrakov commented 4 years ago

Sorry for confusing names. The absorber also does not go into negative, we just call the border that's near min value of x / y /z "negative". So this "negative border" absorber starts at 0 and goes towards the positive direction. I believe your middle in X (X_O) just needs to be in the middle of the domain.

cbontoiu commented 4 years ago

I read all I could about the absorber and thanks for the explanation above. Still I don't know what is its meaning and if I should use YeePML< CurrentInterpolation >: standard Yee solver with PML absorber. And, in this case how? Do I still need the ABSORBER_CELLS[3][2], or maybe something else? I am working with the @develop version.

My worry comes from the fact that a plane wave laser could create a quadrupole pattern as that in the picture attached.

index

The quadrupole pattern develops as the laser interacts with the annular region of carbon atoms and spreads until it occupies the whole simulation domain. There were 32 absorber cells along each of the 6 cube faces. It seems that the laser field is constrained in the simulation box and resonates?

cfg.txt grid.txt laser.txt

sbastrakov commented 4 years ago

Yes, so in dev you could set Solver to YeePML< ... > in fieldSolver.param. In this case the PML absorber will be used, it has its own set of parameters in pml.param file. The parameters for the exponential damping absorber is then ignored. The default values in pml.param should be reasonable, maybe you only need to change the thickness set by NUM_CELLS variable there.

sbastrakov commented 4 years ago

Please note that in case your setup directory does not have a certain .param file, the default version (in include/picongpu/param) will be used. To override this, copy a .param file there and modify or use pic-edit to create it first.

cbontoiu commented 4 years ago

I knew about the param files and everytime I need to modify something I bring a copy of the corresponding .param file from the repository inside spack. Still with the YeePML solver (8 cells and then 16 cells) I get the situation below. Only the BSIEffective and the ADKLinPol are enabled as ionization methods. The laser is polarizaed as LINEAR_X and introduced at the 33rd cell because initPlaneY = 33u but similar things appeared with initPlaneY = 0 Puzzling.

Screenshot_2020-01-29_17-50-06

cbontoiu commented 4 years ago

Also I wonder why the laser field amplitude is wiggling? Is there a precision problem with my data?

laser

sbastrakov commented 4 years ago

Thanks for reporting @cbontoiu , I will take a look at the PML weird results later today. About wiggling: i probably do not follow, since on the picture from your last post the amplitude seems the same along each y = const (horizontal) line?

cbontoiu commented 4 years ago

Thanks for replying. The quadrupole pattern appeared also with the classical absorber. About the wiggling if you choose any value of x and look upwards along the y axis you can notice variations of the colours. I understand the there is ramping of the pulse until it reaches a plateau but this ramping has variations of ~4e10 V/m if you compare the first to colour stripes.

sbastrakov commented 4 years ago

Concerning the latter picture. As far as I can see (and if i did not make a mistake in calculations), your spatial step is 6.19e-11 meters and so the domain size is roughly 13.8 nm. However, the pulse wavelength is 75 nm. So what we see is just the very start of the first wave of the pulse (and even this will not fit your simulation volume).

PrometheusPi commented 4 years ago

@cbontoiu as @sbastrakov correctly pointed out at 33 iterations you laser propagated delta_t * 33 * c = 1.2 * 10^-19 s * 33 * 2.9979e8= 1.2nm << lambda_0 = 75nm. The initialization plane is at 2nm (33 cells), so your laser field should go from ~0.8nm till 3.2nm - which you see in red in the plot above.

Since we usually are not looking at the laser pulse rising edge at such a fine temporal resolution the lines you see might originate from how we initialize the laser field. But my first guess is that this is coming from how you interpolate your data with imshow(?) (pcolormesh) when plotting. Between 0nm and 2nm, there should only be 33 cells but your gray lines suggest a much higher resolution - so might this be a plotting issue (Could you select nearest interpolation)?

cbontoiu commented 4 years ago

Thank you! I will try with the interpolation. Meanwhile, is there a way to change to background of the pictures produced with png.param from black to white please?

sbastrakov commented 4 years ago

Not sure if we support it, @psychocoderHPC do you have any idea? I guess most graphics programs (even MS Paint) support such a color replacement manually, although that would be not automatic.

psychocoderHPC commented 4 years ago

Thank you! I will try with the interpolation. Meanwhile, is there a way to change to background of the pictures produced with png.param from black to white please?

You can try to change in png.param

namespace preParticleDensCol = colorScales::grayInv;

to

namespace preParticleDensCol = colorScales::gray;

Please keep in mind the pngs are only a preview.

cbontoiu commented 4 years ago

One remark here: Although I set PARAM_FIELDSOLVER YeePML in fieldSolver.param and this means that number of cells for the absorber is controlled now by NUM_CELLS in the pml.param file, there is still an error at compilation time if the laser is introduced let's say at the 9th cell while the corresponding ABSORBER_CELLS number is 16 in grid.param. My understanding is that when using YeePML the ABSORBER_CELLS in grid.param is ignored as the pml.param file takes over. Is it true? If yes, the error should not appear.

sbastrakov commented 4 years ago

@cbontoiu you are right, there is a bug with that check for laser init plane always using ABSORBER_CELLS in grid.param regardless of the field solver being used. I will check if this is easy to fix. For now, as a workaround please set ABSORBER_CELLS so that it passes that check. And PML would use NUM_CELLS from pml.param anyways.

Thanks for reporting!

sbastrakov commented 4 years ago

@cbontoiu FYI a fix for the bug you pointed out has just been merged to dev, thanks!