Include temperature deposition shape factor

[emarchet]

This PR includes the possibility to select the order of the plasma temperature deposition shape factor, when applicable.

hipace.temp_depos_order = 2

• Small enough (< few 100s of lines), otherwise it should probably be split into smaller PRs
• Tested (describe the tests in the PR description)
• Runs on GPU (basic: the code compiles and run well with the new module)
• Contains an automated test (checksum and/or comparison with theory)
• Documented: all elements (classes and their members, functions, namespaces, etc.) are documented
• Constified (All that can be const is const)
• Code is clean (no unwanted comments, )
• Style and code conventions are respected at the bottom of https://github.com/Hi-PACE/hipace
• Proper label and GitHub project, if applicable

[emarchet]

Input file

max_step = 0
amr.n_cell = 64 64 128

geometry.prob_lo = -0.39e-3 -0.39e-3 -1.e-3
geometry.prob_hi = 0.39e-3 0.39e-3 0.1e-3  #1.1 mm long plasma
boundary.field = Dirichlet
boundary.particle = Periodic

amr.max_level = 0
hipace.verbose = 1

my_constants.n0 = 8.e21  #m^-3

beams.names = no_beam

hipace.dt = 1.e-12

hipace.do_device_synchronize = 0

plasmas.names = electrons ions
plasmas.neutralize_background = 0

electrons.element = electron
electrons.density(x,y,z) = n0
electrons.ppc = 4 4
electrons.radius = 0.39e-3
electrons.temperature_in_ev = 1

ions.element = H  
ions.density(x,y,z) = n0
ions.ppc = 4 4
ions.radius = 0.39e-3
ions.initial_ion_level = 1
ions.temperature_in_ev = 1

diagnostic.output_period = 1
diagnostic.diag_type = xyz

# Temperature deposition settings
hipace.deposit_temp_individual = 1
hipace.temp_depos_order = 2

Figure script

import numpy as np
import scipy.constants as const

import matplotlib.pyplot as plt

# Set global font sizes for the plots
plt.rcParams['font.size'] = 10          # global default

from openpmd_viewer import OpenPMDTimeSeries

ts_0 = OpenPMDTimeSeries('/data/dust/user/emarchet/runs/branch_Tdep/shapef_0/input_0_3d/diags/hdf5/', check_all_files=False)
ts_1 = OpenPMDTimeSeries('/data/dust/user/emarchet/runs/branch_Tdep/shapef_1/input_0_3d/diags/hdf5/', check_all_files=False)
ts_2 = OpenPMDTimeSeries('/data/dust/user/emarchet/runs/branch_Tdep/shapef_2/input_0_3d/diags/hdf5/', check_all_files=False)
ts_3 = OpenPMDTimeSeries('/data/dust/user/emarchet/runs/branch_Tdep/shapef_3/input_0_3d/diags/hdf5/', check_all_files=False)

iteration = 0

def get_u (ts):
    
    print(ts.avail_fields, '\n')

    ux_sq_e, _ = ts.get_field(field='ux^2_electrons', iteration=iteration, slice_across=None)
    ux_e, info_ux_e = ts.get_field(field='ux_electrons', iteration=iteration, slice_across=None)

    #print(ux_sq_e.ndim)
    #print(ux_e.ndim)

    uy_sq_e, _ = ts.get_field(field='uy^2_electrons', iteration=iteration, slice_across=None)
    uy_e, _ = ts.get_field(field='uy_electrons', iteration=iteration, slice_across=None)

    #ts.get_field?
    #info_ux_e?

    x = info_ux_e.x
    y = info_ux_e.y
    z = info_ux_e.z

    # Calculate the plasma temperatures
    Tx_e = (ux_sq_e - ux_e**2) * const.m_e * const.c**2 / const.e
    Ty_e = (uy_sq_e - uy_e**2) * const.m_e * const.c**2 / const.e

    print(Tx_e.shape, '\n')
    #print(Ty_e.shape)

    return Tx_e, Ty_e, x, y, z


Tx_e_0, Ty_e_0, x, y, z = get_u(ts_0)
Tx_e_1, Ty_e_1, *_, = get_u(ts_1)
Tx_e_2, Ty_e_2, *_, = get_u(ts_2)
Tx_e_3, Ty_e_3, *_, = get_u(ts_3)

# Get the shape of the array (128, 64, 64)
shape = Tx_e_0.shape
x_mid = shape[2]//2     # Middle of the x-dimension
y_mid = shape[1]//2     # Middle of the y-dimension
z_mid = shape[0]//2     # Middle of the z-dimension

# Plot
Tx_list = [Tx_e_0, Tx_e_1, Tx_e_2, Tx_e_3]
Tx2d_list = [Tx[z_mid, :, :] for Tx in Tx_list]
titles = [f"Shapef {i}" for i in range(4)]

fig, axes = plt.subplots(2, 2, figsize=(10, 8))

for i, ax in enumerate(axes.ravel()):
    m = ax.pcolormesh(x, y, Tx2d_list[i].T, shading="auto", cmap='plasma', vmin=0, vmax=2)
    fig.colorbar(m, ax=ax)
    ax.set_title(titles[i])
    ax.set_xlabel("x (m)")
    ax.set_ylabel("y (m)")
    ax.ticklabel_format(axis="both", style="sci", scilimits=(0, 0))
    ax.grid(True)

plt.suptitle("Electron Temperature Tx_e (eV)")
plt.tight_layout()
plt.show()

[MaxThevenet]

Thanks for this PR! Could you fix conflicts with dev and make sure CI passes?

[MaxThevenet]

Thanks for this PR!