External electromagnetic field

[naokikatsura]

Hello, everyone,

Recently, I learned about the existence of the cylindrical coordinate system in EPOCH, compiled the latest version, and started using it. Is it possible to apply an external magnetic field in epoch_axial? I have tried to apply a uniform magnetic field in the x-axis direction, similar to epoch2d, but it seems that the external magnetic field is not being generated correctly.

Here is my input.deck.

[naokikatsura]

input.txt

[Status-Mirror]

Hi Naokikatsura,

Welcome to the cylindrical EPOCH project! It's not currently possible to load an external field in the code. The cylindrical field-solver uses complex field modes instead of real fields, which our input deck does not yet support.

Even if you could initialise a magnetic field, EPOCH's field solver will update it along with the other fields for Maxwell's equations. This is the same issue you would face with epoch2d too.

If you don't mind the fields being updated, then there is a way you can hack the initial conditions (at least for Bx). If you go to housekeeping/setup.F90, you'll see a subroutine called set_initial_values. Here, you can directly set the initialised values of the complex field modes (see DOCUMENTATION.pdf for more detail on these). For your uniform $B_x$, you only need to change the $m=0$ mode, so you can add lines like:

bxm(:,:,0) = 20.0_num
bxm_old(:,:,0) = 20.0_num

Hope this helps, Stuart

[naokikatsura]

Dear Stuart,

Thank you for giving your insight into this issue. I tried to initialize magnetic field directly by edditing 'set_initial_values' in housekeeping/setup.F90. As you told me, the magnetic field is updated, which is undesireble for me.

So far, I have applied a uniform magnetic field, but my ultimate goal is to observe the interaction between the magnetic field and the plasma due to the coil current. I would like to set the external magnetic field driven by the coil current as shown in the figure below.

The current should be steady, and keep flowing during the simulation. I would be grateful for any advice.

Best regards, Naoki

[Status-Mirror]

Hi Naoki,

What you want to do is possible, but it requires quite a bit of tinkering with the source-code. I'll walk you through a quick example, where I'll assume you're using the default macro-particle shape.

The scripts responsible for calculating the magnetic field on a particle are present in epoch_axial/src/include/triangle/b_part.inc. Here we calculate the cylindrical components of the magnetic field $(B_x, Br, B\theta)$, which we convert to Cartesian components $(B_x, B_y, B_z)$ for the default particle pusher. If you wanted a uniform, unchanging magnetic field on the particles, you could just add bx_part = bx_part + 20.0_num to the end of this file.

If you want a magnetic field which varied with position, then I think you're able to use the variables part_x, part_y, part_z and part_r in this file. These refer to the position of the current particle, so you could write:

bx_part = bx_part + part_x

To have a linearly rising external magnetic field in the $x$ direction. Note that part_x is measured in metres, and bx_part in Tesla.

I think this would work, but you could hit problems in the field ionisation or qed packages. Let me know if you run into problems with this.

Hope this helps, Stuart

[naokikatsura]

Dear Stuart,

Thank you for your reply. This way I can set the value of the magnetic field acting on the particle, but the electromagnetic fields are not updated. So why not add bxm_external(ix, ir, im) in the subroutine updata_b_field in src/fields.F90 as following?

bxm(ix, ir, im) = bxm(ix, ir, im) - (im_fac_x * erm(ix, ir, im) + 0.5_num * (etm(ix, ir, im) + etm(ix, ir-1, im)) / r_d + (etm(ix, ir, im) - etm(ix, ir-1, im)) / dy) * 0.5_num * dt + bxm_external(ix, ir, im)

Best regards, Naoki

[Status-Mirror]

Hi Naoki,

I thought you didn't want the external magnetic field to update? Your solution will add a magnetic field to the field-solver, and any gradients in this field will drive an electric field. This will become unphysical, as your driven electric field will not be able to affect the external magnetic field which produces it, and the fields could blow up. I'm not sure if this is the behaviour you want.

Adding a line like:

bx_part = bx_part + part_x

will keep the original magnetic field updated by the code bx_part unchanged, and add an additional term to represent your constant external field (here this is the part_x term).

If you did want an external field which affects the field-solver in the code, this would be a difficult thing to implement. Your array bxm_external would need to be added to the load-balancer, as the MPI domains in EPOCH can shift to follow the particles around. Modifying your external field from the input deck would also require some code extensions, so it would be best to initialise it in the source-code itself. If you do decide to go down this road, try and find all the references to bxm_old in the code, as your external field will need equivalent lines too.

Hope this helps, Stuart

[naokikatsura]

Dear Stuart,

Thank you for your comments. I apologize for my previous misunderstandings regarding the simulation codes and the calculation of electromagnetic fields. I want the electromagnetic field to be updated just like setting fields block in cartesian epoch2d. In other words, I supose it would be more suitable for my work to initialise the field in housekeeping/setup.F90, as you told me in your first comment.

However, I am encountering an issue where an electromagnetic field is unexpectedly excited from the left-hand (xmin) boundary after initializing the magnetic field.

Here are the initialization settings I used in housekeeping/setup.F90:

bxm(:,:,0) = 4.0_num brm(:,:,0) = 4.0_num bxm_old(:,:,0) = 4.0_num brm_old(:,:,0) = 4.0_num

I also reverted my changes in src/fields.F90 and src/include/triangle/b_part.inc. After running the code, I calculated the electromagnetic field strength at θ=0 by summing up the real parts of modes 0 and 1, and created 2D contour figures. Black lines show the stream plot of the magentic field lines in the figure of magnetic flux density strength (the bottom one).

It seems this fields excitation occurs when I initialized brm. I confirmed there is no fields excitation when I only initialized bxm. I am unsure how this electromagnetic field was generated. Could you advise on any possible causes?

Best regards, Naoki

[Status-Mirror]

Hey Naoki,

Do you have the input deck which generated this simulation?

Cheers, Stuart

[naokikatsura]

Dear Stuart,

Attached is the input.deck file for this simulation: input_deck_20240806.txt Since I injected no particles in the domain, most of the part for species block are commented out.

I have examined the time variations of the electromagnetic field components Ex, Er, Et, Bx, Br, and Bt. It appears that the Et, Bx, and Br components are showing unreasonable behavior, while the others remain steady.

This may be because of the boundary condition of xmin and xmax. In this simulation, these boudary conditions are set to be simple_laser. However, when I set periodic boundary condition, the electromagnetic field showed steady.

Best regards, Naoki

[Status-Mirror]

Hey Naoki,

I think this is likely the reason you're seeing this error. I've been looking into the fields myself, and I think the boundary conditions may be overwriting your set fields. In the ghost-cells surrounding the simulation window (where fields are estimated around the edges so that particles passing the boundary see a continuous field before they exit) the fields can be clamped such that they pass through 0 on the simulation edge. This sets up a sharp gradient in your otherwise uniform field, and field gradients drive further fields in Maxwell's equations.

I'm surprised these errors are less pronounced when you have periodic fields in $x$ though. There is special boundary treatment in the $r$ direction too, and you cannot switch to periodic boundaries in this direction.

I had a tricky time setting the boundary conditions which are currently in the code, I'm not sure I can provide much help for your custom setups. For a quick fix, it might be worth working in geometries where your custom fields drop to zero on the simulation edges, so there's nothing to interfere with the boundary field solver. You'd also avoid these problems if you just added an extra b_part to your particles as discussed before, although I think you did want some interaction with the field solver.

Cheers, Stuart

[naokikatsura]

Dear Stuart,

I sincerely apologize for the delayed response due to personal issues. Thank you also for your valuable comments.

I believe that the rewriting of field components due to boundary conditions has led to sharp gradients and unphysical electromagnetic fields as you mentioned. Have you considered implementing other open boudary condition such as PML (Perfectly Matched Layer) or CPML (Convolutional Perfectly Matched Layer) boundaries? If you have no plans to implement them, I would like to personally introduce CPML boundary conditions to this code. I would appreciate any advice you could provide on this.

Best regards, Naoki

[Status-Mirror]

Hi Naoki,

I'm aware the cartesian forms of EPOCH support CPML boundaries, but I could never find a good reference for these which explains how they work, so I've always avoided them. You can check the cartesian epoch2d source-code for an example CPML implementation on a 2D grid, and the DOCUMENTATION.pdf file on this repository explains the the existing field-solver.

I don't have plans to implement CPML myself, but could you point me to a good reference which explains it if you have one? I probably should learn how it works at some point.

Cheers, Stuart