ECP-WarpX / WarpX

WarpX is an advanced electromagnetic & electrostatic Particle-In-Cell code.
https://ecp-warpx.github.io
Other
310 stars 196 forks source link

Issue with Limiting External Electric Field to Specific Boundary in 3D Simulation #5046

Open jinze-liu opened 4 months ago

jinze-liu commented 4 months ago

For the three-dimensional case, when adding an external electric field, I want to add it on the line in the x-direction where yz=0 (i.e., the boundary). According to the documentation, I used upper and lower to limit the region of application. However, I found that it didn't work (I tried both AnalyticInitialField and AnalyticAppliedField ), and the value of the external electric field was still added to the entire simulation area. My code is as follows:

from pywarpx import callbacks, fields, libwarpx, picmi
import numpy as np

# Define simulation domain
nx, ny, nz = 256, 256, 256
Lx, Ly, Lz = 0.04188, 0.04188, 0.04188  # domain size in mete

# Define the simulation grid
grid = picmi.Cartesian3DGrid(
    number_of_cells=[nx, ny, nz],
    lower_bound=[0, 0, 0],
    upper_bound=[Lx, Ly, Lz],
    lower_boundary_conditions=['periodic', 'periodic', 'periodic'],
    upper_boundary_conditions=['periodic', 'periodic', 'periodic']
)

# Electromagnetic solver
solver = picmi.ElectromagneticSolver(
    grid=grid,
    cfl=1.0
)

# Plasma parameters
density = 7.9e18  # in m^-3
temperature = 0  # in Kelvin

# Define electron species
electrons = picmi.Species(
    particle_type='electron',
    name='electrons',
    initial_distribution=picmi.UniformDistribution(
        density=density,
        upper_bound=[None, None, None],
        directed_velocity=[0., 0., 0.]
    )
)
# Initialize the simulation
sim = picmi.Simulation(
    solver=solver,
    max_steps=6000,
    verbose=1,
    particle_shape=3,
    warpx_current_deposition_algo='esirkepov',
)

# Add species to the simulation
sim.add_species(electrons, layout=picmi.GriddedLayout(grid=grid, n_macroparticle_per_cell=[1, 1, 1]))

# Define and add the external electric field
external_Efield = picmi.AnalyticInitialField(
    Ex_expression="10.0",   ##Ex_expression= 10.0*cos(4000*t)
    Ey_expression="0.0",
    Ez_expression="0.0",
    lower_bound=[0.0, 0.0, 0.0],
    upper_bound=[Lx, 0.0, 0.0]
)

sim.add_applied_field(external_Efield)  # Correct method to add applied fields
sim.step()

Where did I go wrong? Are there any other methods to add an external electromagnetic field at a fixed plane or fixed point?

dpgrote commented 4 months ago

WarpX does not currently support the lower_bound and upper_bound options. However, you can get the same affect by adding the limits to the expression. Ex_expression="10.0*(y==0)*(z==0)"

jinze-liu commented 4 months ago

Thank you very much! So, if I want to add a time-varying external electromagnetic field in a specific region (AnalyticAppliedField), following your example, can I use Ex_expression="10.0*cos(ωt)*(y==0)*(z==0)?

dpgrote commented 4 months ago

Yes, that would work, though it would be Ex_expression="10.0*cos(w*t)*(y==0)*(z==0) and you wold have to define w.

jinze-liu commented 4 months ago

Thank you again! Your answer saved me a lot of time! Can we consider modifying the lower_bound and upper_bound in the Applied fields section of the Parameters: Python (PICMI) documentation? I initially thought this was feasible while reading it, but currently, it appears that these two are not supported.

ax3l commented 4 months ago

Great, glad this could be solved.

For action items we need to address for WarpX' PICMI implementation:

  1. We should error out if defined PICMI parameters are silently ignored.
  2. In the special case here, we could in the backend of WarpX wrap the passed analytic field with if(..., ..., if( .., ..., ...)) (six if: per direction and up/lo) for the lower and upper bound to support the PICMI feature in AnalyticInitialField and AnalyticAppliedField without changing any backend - or we could change the C++ backend implementation, of course
jinze-liu commented 4 months ago

@ax3l Hi! Using the same simulation area and parameters, settingEx_expression="10.0*cos(4000*t)*(y==0)*(z==0)" will result in all diagnosed data being zero. When using an if-statement you suggested, such as Ex_expression="10*cos(4000*t)* if (y == 0 and z == 0) else 0" , a syntax error is shown. This might be due to my misunderstanding of how to use conditional statements. Could you explain the usage of conditional statements in more detail?

My simulation space length in the xyz directions is 0.000134, so when I try using Ex_expression="10.0*cos(4000*t)*(y<1e-4)*(z<1e-4)", the diagnosed data is not zero. Is this behavior due to numerical errors, meaning that strict conditions like y=0 and z=0 are not achievable?

Thanks!

dpgrote commented 4 months ago

The conditional expression would be if (y==0 and z==0, 10*cos(4000*t), 0.). Though this should give the same result as what I suggested. In your grid setup, are there grid cells that lie at exactly y and z equal to zero?

jinze-liu commented 4 months ago

My grid points do start from (0, 0, 0), my grid is set to

grid = picmi.Cartesian3DGrid(
    number_of_cells=[nx, ny, nz],
    lower_bound=[0, 0, 0],
    upper_bound=[Lx, Ly, Lz],
    lower_boundary_conditions=['periodic', 'periodic', 'periodic'],
    upper_boundary_conditions=['periodic', 'periodic', 'periodic'])

Therefore, I'm puzzled why setting y==0 and z==0 results in all data being 0.

dpgrote commented 3 months ago

I had a chance to look at this again. I gave the wrong advice above. The expressions in AnalyticInitialField don't allow for time dependence, only spatial dependence. With the electromagnetic solver, the fields are only applied to the grid during initialization, which is at t=0.

jinze-liu commented 3 months ago

@dpgrote Hi!In fact, if I want to define a time-dependent condition, I use the AnalyticAppliedField, which should work. The main issue is still with the boundary conditions that were not correctly recognized and were silently ignored.