-
Notifications
You must be signed in to change notification settings - Fork 120
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Segmentation fault (core dumped) on running for D-D fusion simulation #615
Comments
attachment |
The function of uploading attachment seems to get wrong. I will post text directly here. namelist# ----------------------------------------------------------------------------------------
# SIMULATION PARAMETERS FOR THE PIC-CODE SMILEI
# ----------------------------------------------------------------------------------------
import math
# Reference
l0 = 2.0 * math.pi # laser wavelength 800nm/1um (2*pi*lr)
t0 = l0 # optical cicle 2.66851fs/3.333fs
# Simulation parameters for 800nu
# lm = l0 * 5.0 / 4.0 # 1um
# lr = 8e-7 / (2.0 * math.pi) # show the value of unit of length in smilei as SI unit
# tm = 5.0 * math.pi / 2.0 # 3.33fs
# Simulation parameters for 1um
lm = l0 # 1um
lr = 1e-6 / (2.0 * math.pi) # show the value of unit of length in smilei as SI unit
tm = t0 # 3.33fs
Lsim = [4.0 * lm, 2.0 * lm] # length of the simulation 4um × 2um
# Lsim = [5.0 * lm, 15.0 * lm]
# Lcell = [0.01*l0, 0.01*l0]
resx = 128.0 # nb of cells in on a lm
resy = 128.0
Tsim = 200.0 * t0 # duration of the simulation for 10 laser wavelength (about 48 tm)
# rest = 60. # time of timestep in one optical cycle
# Laser
omega = 1.0 # = 2PI/t0
a0 = 5.0 # 8.4 、10 、11.85
# particles in per cell
ppc = 16
def preprocess():
# density = 22.2203 # density of Carbon atom on 800nm
density = 34.6643 # density of Carbon atom on 1um
if smilei_mpi_rank == 0:
import numpy as np
from numpy import linspace, meshgrid, ceil
import h5py
# 输出模拟尺寸 及 网格大小(用作delta)
with open("Box.txt", "w") as f:
f.write(f"{Lsim[0]} {Lsim[1]} {2.0*lm/resx}")
vertices = [
(0.2 * Lsim[0], 0.2 * Lsim[1]), # (0.8 * Lsim[0], 0.1 * Lsim[1]),
(0.2 * Lsim[0], 0.8 * Lsim[1]), # (0.8 * Lsim[0], 0.9 * Lsim[1]),
(0.8 * Lsim[0], 0.8 * Lsim[1]), # (0.9 * Lsim[0], 0.9 * Lsim[1]),
(0.8 * Lsim[0], 0.2 * Lsim[1]), # (0.9 * Lsim[0], 0.1 * Lsim[1]),
] # 多边形顶点数据
with open("vertices.txt", "w") as f:
for vertex in vertices:
f.write(f"{vertex[0]} {vertex[1]}\n")
vertices = np.loadtxt("vertices.txt")
p = linspace(0.0, Lsim[0], int(resx * ceil(Lsim[0] / lm)) + 1)
q = linspace(0.0, Lsim[1], int(resy * ceil(Lsim[1] / lm)) + 1)
xx, yy = meshgrid(p, q)
points = np.vstack((xx.ravel(), yy.ravel())).T
def is_point_in_polygon(point, polygon):
"""
判断一个点是否在多边形内部
Args:
point: tuple, 待判断的点坐标
polygon: ndarray, 多边形的顶点坐标, 形状为(n, 2)
Returns:
bool, True表示在多边形内部, False表示在外部
"""
# 将多边形的顶点按顺时针排序
x, y = polygon[:, 0], polygon[:, 1]
cx, cy = x.mean(), y.mean()
angles = np.arctan2(y - cy, x - cx)
sorted_idx = np.argsort(angles)
polygon = polygon[sorted_idx]
# 计算射线与多边形的交点个数
x, y = polygon[:, 0], polygon[:, 1]
n = len(polygon)
count = 0
for i in range(n):
x1, y1 = x[i], y[i]
x2, y2 = x[(i + 1) % n], y[(i + 1) % n]
if (point[0] == x1 and point[1] == y1) or (
point[0] == x2 and point[1] == y2
):
return True # 点在多边形上
if (y1 > point[1]) != (y2 > point[1]) and point[0] < x1 + (x2 - x1) * (
point[1] - y1
) / (y2 - y1):
count += 1
return count % 2 == 1
is_inside = np.array([is_point_in_polygon(p, vertices) for p in points])
lenx = p.size
leny = q.size
data_C = np.zeros((lenx, leny)) # for Carbon atom
data_D = np.zeros((lenx, leny)) # for Deuterium atom
# data_E = np.zeros((lenx, leny)) # for Electron
for j in range(leny):
for i in range(lenx):
data_C[i, j] = density if (is_inside[j * lenx + i]) else 0.0
data_D[i, j] = density * 2.0 if (is_inside[j * lenx + i]) else 0.0
# data_E[i, j] = density*8. if (is_inside[j*lenx+i]) else 0.0
with h5py.File("CD2-profile.h5", "w") as f:
g = f.create_group("density_group")
g.create_dataset("C_profile", data=data_C)
g.create_dataset("D_profile", data=data_D)
# g.create_dataset("E_profile", data=data_E)
Main(
geometry="2Dcartesian",
interpolation_order=2,
# reference_angular_frequency_SI=2*math.pi*3e8/1.e-6/0.8, # Lr=0.8/2pi
reference_angular_frequency_SI=3e8 / lr, # ωr = c/lr
grid_length=Lsim,
cell_length=[lm / resx, lm / resx], # [l0/resx, l0/resx],
number_of_patches=[8, 8],
# timestep=t0/rest,
timestep_over_CFL=0.9,
simulation_time=Tsim,
EM_boundary_conditions=[
["silver-muller"],
["periodic"],
],
# print_every = 1,
cluster_width=1,
)
LaserPlanar1D(
box_side="xmin",
a0=a0,
omega=omega,
polarization_phi=0.0,
ellipticity=0.0,
# time_envelope=tconstant(),
time_envelope=tgaussian(
start=0.0, duration=33 * t0, fwhm=15.0 * t0, center=15.0 * t0
),
)
# Deuterium
Species(
name="Deuterium",
position_initialization="regular",
momentum_initialization="mj",
ionization_model="tunnel",
ionization_electrons="Electron",
temperature=[5.0e-8],
particles_per_cell=ppc,
# c_part_max=1.0,
mass=3670.03, # 3675.18,
charge=0.0,
atomic_number=1,
number_density="CD2-profile.h5/density_group/D_profile",
# mean_velocity=[0.5, 0., 0.],
# time_frozen=0.0,
boundary_conditions=[
["stop", "stop"],
["periodic", "periodic"],
],
)
# Carbon
Species(
name="Carbon",
position_initialization="regular",
momentum_initialization="mj",
ionization_model="tunnel",
ionization_electrons="Electron",
temperature=[5.0e-8],
particles_per_cell=ppc,
# c_part_max=1.0,
mass=21866.04,
charge=0.0,
atomic_number=6,
number_density="CD2-profile.h5/density_group/C_profile",
# mean_velocity=[-0.5, 0., 0.],
# time_frozen=0.0,
boundary_conditions=[
["stop", "stop"],
["periodic", "periodic"],
],
)
# Helium-3
Species(
name="He3",
position_initialization="regular",
momentum_initialization="cold",
# ionization_model='tunnel',
# ionization_electrons='Electron',
# temperature=[5.0e-8],
particles_per_cell=ppc,
# c_part_max=1.0,
mass=5495.66,
charge=0.0,
atomic_number=2,
number_density=0.0, # 'CD2-profile.h5/density_group/C_profile',
# mean_velocity=[-0.5, 0., 0.],
# time_frozen=0.0,
boundary_conditions=[
["stop", "stop"],
["periodic", "periodic"],
],
)
# Neutron
Species(
name="Neutron",
position_initialization="regular",
momentum_initialization="mj",
temperature=[5.0e-8],
particles_per_cell=ppc,
# c_part_max=1.0,
mass=1822.17,
charge=0.0,
atomic_number=0,
number_density=0.0, # 'CD2-profile.h5/density_group/C_profile',
# mean_velocity=[-0.5, 0., 0.],
# time_frozen=0.0,
boundary_conditions=[
["stop", "stop"],
["periodic", "periodic"],
],
)
# Electron
Species(
name="Electron",
position_initialization="regular",
momentum_initialization="mj",
temperature=[5.0e-8],
particles_per_cell=ppc,
# c_part_max=1.0,
mass=1.0,
charge=-1.0,
number_density=0.0, # 'CD2-profile.h5/density_group/E_profile',
# mean_velocity=[0.5, 0., 0.],
# time_frozen=0.0,
boundary_conditions=[
["stop", "stop"],
["periodic", "periodic"],
],
)
Collisions(
species1=["Deuterium"],
species2=["Deuterium"],
# coulomb_log=0.001,
nuclear_reaction=["He3", "Neutron"],
debug_every=10,
)
every = 400
DiagParticleBinning(
deposited_quantity="weight",
every=every,
flush_every=5 * every,
species=["Neutron"],
axes=[["ekin", 0.0, 2.0, 10000, "edge_inclusive"]], # 0keV-1MeV
# axes=[["ekin", 0.002, 20, 1000, "logscale", "edge_inclusive"]], # 1keV-10MeV
)
DiagParticleBinning(
deposited_quantity="weight",
every=every,
flush_every=5 * every,
species=["Electron"],
axes=[["ekin", 0.0, 2.0, 10000, "edge_inclusive"]], # 0keV-1MeV
# axes=[["ekin", 0.002, 20, 1000, "logscale", "edge_inclusive"]], # 1keV-10MeV
)
DiagParticleBinning(
deposited_quantity="weight",
every=every,
flush_every=5 * every,
species=["Deuterium"],
axes=[["ekin", 0.0, 2.0, 10000, "edge_inclusive"]], # 0keV-1MeV
# axes=[["ekin", 0.002, 20, 1000, "logscale", "edge_inclusive"]], # 1keV-10MeV
)
DiagScalar(
every=every,
vars=["Dens_Neutron", "Ntot_Neutron"],
precision=5,
) outlog
|
Are you able to reproduce the bug with a smaller simulation? |
Please also provide an input that does not rely on external resources ( |
Here is a new inputfile which has smaller scale and don't rely on any external resources. But the same error still occurs! |
For the record I paste the input here
|
At which timestep do you see the error occur? Is this reproducible? |
I was able to run your reduced simulation without error. Please specify your version of smilei (and where you got it) and of MPI ( |
After simplifying the code overhead, I ran it again and encountered the same error as before when it reached 30%. Strangely, now when I run it again, it seems to be working properly again. I have been using Smilei-4.7.tar.gz downloaded from GitHub, and the version of MPI is mpirun (Open MPI) 4.1.1. I will make further attempts. |
I have carried out simulations of two sizes and repeated each three times. The same error occurred each time, but the timing of the error was not consistent. Attached are my namelist and outlog files. |
Can you try using less threads? Like 16? You are using way too many compared to the number of patches. Your pushtime is catastrophic because of this, and I wonder whether it could cause memory issues. |
Additionally, I don't think you have provided the exact same input file as that you have used because the basic data I obtain at the start of the log is different from yours |
I have tried repeatedly to run the program using I'm sorry, but I have also been puzzled by this bug. I am no longer sure if it is caused by the fusion part, as this error can occur with changes in simulation scale, grid division, and even density distribution. The most lethal thing is that even without any changes, you can repeatedly run it and get different results. |
The only thing that remains unchanged is the error message.
|
Could you try with pure MPI (no openmp), using OMP_NUM_THREADS set to 1? If it works, try to recompile with the no_mpi_tm option:
Then run again with openmp activated |
Thank you a lot! According to your previous advice, I carried out the operation (let OMP_NUM_THREADS=1) and the previous error magically disappeared. I repeatedly ran it in a loop and modified the initial conditions, and now I am confident that it can run normally. So what caused this? |
Have you tried what I indicated above? no_mpi_tm If it works again then the problem is that mpi was not compiled with the THREAD_MULTIPLE option |
I just tried with pure MPI, then it worked well. So I will continue to push forward with my work in this way for now. After I'm done with this busy period, I will come back and test the recompile option. |
Let me close the issue. Please reopen if the solution did not work in the end |
Recently, I have try Worst of all, the error mentioned above will still occur. It seems to work well only with pure MPI. |
Is this still happening with the same previous input files? And did you |
Yes, I did |
Description
Hi, friends! Some error occur when I use SMILEI for D-D fusion. Code will be Interrupted with "Segmentation fault (core dumped)" for some cell number(eg. resx=resy=128), but run well for others(eg. resx=resy=64). And it also run well while I removed the fusion section. Any suggestion for me? Thanks a lot!
Parameters
Distributor ID: Ubuntu
Description: Ubuntu 20.04.5 LTS
Release: 20.04
Codename: focal
g++ --version
)g++ (Ubuntu 9.4.0-1ubuntu1~20.04.1) 9.4.0
mpic++ --version
,mpic++ -show
)g++ (Ubuntu 9.4.0-1ubuntu1~20.04.1) 9.4.0
g++ -I/usr/local/openmpi/include -pthread -Wl,-rpath -Wl,/usr/local/openmpi/lib -Wl,--enable-new-dtags -L/usr/local/openmpi/lib -lmpi
h5cc --version
)g++ (Ubuntu 9.4.0-1ubuntu1~20.04.1) 9.4.0
python --version
)Python 3.8.10
Namelist and outlog are attached.
The text was updated successfully, but these errors were encountered: