ComputationalRadiationPhysics / picongpu

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

Dependence of simulation output from cell size #2227

Closed NastasiaM closed 7 years ago

NastasiaM commented 7 years ago

I simulated the plasma behavior and found that the result depends on the simulation cell size. Here are two simulations, one with cell size 2.205658e-9 * 2.205658e-9 * 2.546875e-9 and another - twice bigger in all dimensions. The whole picture looks quite same, but the maximum value varies a lot.

First picture with smaller cell size. Max value ~3.5e9 /bigdata/hplsim/external/mukhar40/colloidalMeltingHCPNEW_Int3.0_big_2ioniz/

ed1

And second with cell size twice bigger. /bigdata/hplsim/external/mukhar40/colloidalMeltingHCPNEW_Int3.0_2ioniz/

image

Is it some bug or I am doing something wrong? I understand that energy density is averaged over cell size but the result differs to much for me.

PrometheusPi commented 7 years ago

Hi @NastasiaM,

thank you for posting this issue report.

PIConGPU uses its own unit system in order to deal best with floating point limitations. Thus the physical quantities used and calculated in PIConGPU are scaled according to this unit system, which itself depends on quantities like the cell size, time step duration, expected density, etc. Thus a quantity in PIConGPU units has a different value depending on the resolution. However, after converting these quantities to SI units, there should be no dependency on resolution anymore. In hdf5, this conversion is done multiplying with the unitSI attribute.

Jugging from your plots, you retrieve your data using openPMD_viewer - is this correct? This should always do the conversion for you - but just to be sure, which version of openPMD_viewer are you using? There might be the chance, that openPMD_viewer or PIConGPU do something wrong with conversions here? Thus, could you additionally tell me what version of PIConGPU you are working with and whether your simulation is 2D or 3D?

Furthermore, such changes with resolution might also point toward a numerical issue like badly resolved plasma or laser evolution: Thus could you please give a brief summary of what you are simulating? What maximum plasma densities were there so far in the simulation? How good did you resolve time in both simulations? Why do you expect such high energy densities at this point in space?

Best, @PrometheusPi

PrometheusPi commented 7 years ago

Looking at your pictures, some more questions arose:

A first look at the energy density code did not reveal any misconceptions: the value is correctly divided by the cell volume and thus a density field, not a binning array.

NastasiaM commented 7 years ago

Morning!

1) Yes, it is openPMD_viewer 0.5.4. Actually, I wanted to ask you if there might be something wrong with energy conversion? 2)3D simulation 3) I simulate plasma creation in the colloidal crystal. Simulation parameters:

Absorbing in the top and bottom of the simulation box Other – periodic boundary conditions

Maximum possible electron density(H+ C4+) 2.31e23/cm3 Maximum possible plasma frequency 2.71e16 1/s Plasma wavelength 69 nm 30 points for one wavelength

GRID Δ Y=2.206 nm Δ X= Δ Z=2.547 nm Δ T=4.62 as 1800 electrons per cell (Interaction of electrons inside one cell are NOT calculated) 27 Macro particles per cell This is a smaller grid and the grid of the second simulation was twice bigger in all 4 (t,x,y,z) dimensions.

Simulation box looks like this: (also absorbing condition in the bottom of simulation box)

image

4)The spot of high energy is stationary. 5)What is a moving window simulation?)) If I don't know what is it probably I am not using it... 6) energy histogram - yes I run it with this flag

This is one for smaller simulation cell image And this is for bigger image

I hope this answers your questions and it will help you to answer mine :)

Best, Nastasia.

n01r commented 7 years ago

Hi Nastasia,

one problem might be your temporal resolution. In order to resolve the plasma frequency reasonably we usually aim for omega_pe * dt <= 0.1. But here according to your maximum plasma frequency and time step we arrive at 0.1252. If you now also double the length of your time step it will be 0.2502. That might already be too close to the lower resolution limit, especially if you get locally increased electron densities due to the plasma dynamics.

Best, Marco

PrometheusPi commented 7 years ago

Could you please give us read access to both simulations via

chmod -R  go=u-w /bigdata/hplsim/external/mukhar40/colloidalMeltingHCPNEW_Int3.0_big_2ioniz/
chmod -R  go=u-w /bigdata/hplsim/external/mukhar40/colloidalMeltingHCPNEW_Int3.0_2ioniz/

so that we can test our hypotheses directly?

NastasiaM commented 7 years ago

Hi! I gave you access yesterday morning. Hope you noticed it. If you have some progress in testing please write me! Best, Nastasia.

PrometheusPi commented 7 years ago

First analysis revealed that your temporal resolution is too large for the less well resolved simulation. Here I took the maximum electron density in the simulation and compared it to the simulation time step duration.

For the last time step in each simulation: double resolution:

delta_t * omega_pe(max) = 0.0916
n_max = 1.46e+23 cm^-3

single resolution:

delta_t * omega_pe(max) = 0.1821
n_max = 1.45e+23 cm^-3

Thus the single resolution simulation is not capable of resolving the plasma dynamic thus creating nonphysical results.

Furthermore, the simulation are not equal with respect to the laser time correlation. In laser.param, bot simulations initialize the laser with an offset of 31 cells, but the double resolution simulation would need to initialize the laser at 62 cells to have the same spatial/temporal offset as the single resolution simulation. However this should be just a slight difference since your wavelength is much larger than the simulation box.

PrometheusPi commented 7 years ago

In the worst case the slight difference in position of your laser init plane L_y = 31 * delta_y / lambda_0 = 0.08 leads to a different energy density E ~ v^2 ~ E_x^2 of up to 25% of the maximum value (of course the relative difference can be infinitely large).

Perhaps adjusting the init plane to the same (spatial) value will reduce the difference even though the plasma dynamic is not well resolved.

PrometheusPi commented 7 years ago

Judging from your energy density however, resolving the plasma frequency is the main issue: hc_ion_issue01 (solid lines are maximum values, dashed lines are mean values)

The mean energy evolves smoothly thus it is not an issue of a more or less wide spaced sampling of the phase of the laser and the thus resulting kinetic energy.

NastasiaM commented 7 years ago

Thanks!

So now my question is if the position of the laser init plane leads to a different energy density up to 25% is there any "proper" position of the laser init plane which gives results closer to reality? And why actually this effect occurs?

And just to confirm my understanding... the results with higher temporal resolution are more reliable and it shouldn't happen that if I make the temporal resolution even better results will change significantly? Can you check it (I also can do this but I think you know the faster way)?

PrometheusPi commented 7 years ago

You are welcome.

The slight shift in laser phase (due to the spatial shift of the init plane) could alter the energy dramatically (as seen in a single snap shot in time) if the kinetic energy is mainly caused by the laser field. As an example imagine a single electron in a laser field. It gains kinetic energy and loses it again with respect to the velocity squared. If I now have a second simulation, with a slight shift in phase, the difference between the energies of the electron can be up to 25% of the maximum value for the specific phase shift caused by the init plane difference. For a shift of pi/2 it would be 100%. However in your case: the laser does not penetrate the droplets/spheres thus it is not the main factor of the electron energy. Thus the error in the init plane is not the reason for the problem, (as seen in the above plot as well). The position of the laser init plane is irrelevant as long as it is the same for both simulations. But in both simulations the laser init plane is 31 cells from the y=0 boundary away. For the double resolution 62 cells from y=0 would be equivalent to the single resolution simulation. This offset in space causes a (31 delta_y(single) - 31 delta_y(double)) = -31 delta_y(double) spatial difference, with is equivalent to a -31 delta_y(double) /lambda_0 = -0.08 = -4.6° phase difference.

Yes, higher temporal (and spatial) resolution will improve your numerical accuracy. With increased resolution you will see a trend to a stable solution (hopefully already close to your double resolution setup). Sorry, since I do not know all of your setup parameters I will not be able to rerun your simulation. Please increase the resolution further in your setup, recompile and resubmit. After the better resolved simulation has finished, we can have a look at it together.

n01r commented 7 years ago

Hey @NastasiaM , any news on the matter? :)

NastasiaM commented 7 years ago

Hi everyone! Yes, I tried with twice better resolution and it seems that the simulation result doesn't change too much (less than 20%). So I think now everything is fine. Thanks for help! Nastasia.

PrometheusPi commented 7 years ago

You are welcome. Great that your simulation now show consistent results with increased resolution simulations.