CRPropa / CRPropa3

CRPropa is a public astrophysical simulation framework for propagating extraterrestrial ultra-high energy particles. https://crpropa.github.io/CRPropa3/
https://crpropa.desy.de
GNU General Public License v3.0
65 stars 66 forks source link

Segmentation fault by using transformModelVector function #397

Closed JassimarS closed 5 months ago

JassimarS commented 2 years ago

Hi all, I am new to CRPropa but I was just playing around and trying to replicate some of the example plots on the documentation website and I had an issue with one on this page: https://crpropa.github.io/CRPropa3/pages/example_notebooks/galactic_lensing/lensing_maps.v4.html

The healpy dipole plot works correctly but then applying a magnetic lens using the transformModelVector command (as done in the documentation) produces an output:

Segmentation fault (core dumped)

I am not even using a custom lens, I am using an example lens downloaded from: https://web.physik.rwth-aachen.de/Auger_MagneticFields/PARSEC/downloads.php. This lens does work with the plots on the previous page so I do not think that is the issue.

I am quite unsure what's going on. Assuming there is nothing wrong on my end, I looked at the source code for the function to look for oddities and I found a function 'prod_up' inside the definition of 'transformModelVector' which I then traced back to the ModelMatrix.cpp file. Since this is a seg fault, I think there may be an issue with the copy-pasting of model data into the different variables. I have not dealt with Eigen before so maybe even there.

My code just in case:

import crpropa
import healpy
import matplotlib.pyplot as plt
import numpy as np
import copy
# The lens can be downloaded here: https://crpropa.github.io/CRPropa3/ --> Additional Resources
lens = crpropa.MagneticLens('JF12regular_Gamale/lens.cfg')

#generate random input map
#including a monopole and a dipole
l = [1,1,0,0]
inputMap = abs(healpy.synfast(l, nside=64))

#plot map
healpy.mollview(map=inputMap, title='Unlensed')
plt.savefig('random_unlensed_map.png')

#copy input data as we need it later
outputMap = copy.copy(inputMap)
print(type(inputMap))     # Debugging attempt
print(type(outputMap))    # Debugging attempt

#transform map using a rigidity of 5 EeV
lens.transformModelVector(outputMap, 5 * crpropa.EeV)

#plot transformed map
healpy.mollview(map=outputMap, title='Lensed')
plt.savefig('random_lensed_map.png')
rafaelab commented 2 years ago

Hi @JassimarS I confirm the bug with the current development version of CRPropa on OSX, using both python 3.9 and 3.10. However, using CRPropa 3.1.7 on a Linux system with python 3.7 the example ran just fine.

QuinFin commented 1 year ago

Hey, having the same exact problem when running python 3.7.13 and CRPropa 3.1.7 on a Linux system. After applying the changes suggested here https://github.com/tdwiser/CRPropa3/commit/a9830a9e, https://github.com/tdwiser/CRPropa3/commit/a303407 I pass all tests but get Segmentation fault (core dumped) when trying to call "lens.transformModelVector()". Sadly running out of ideas what to do now.

JassimarS commented 1 year ago

Apologies for returning to this so late but this issue is still not resolved for me. I started working on other things and forgot about this, but having returned to this, the issue still persists. I tested this code on the interactive online demo a while ago on the VISPA cluster (link) which uses a python version of 3.8.10 and it displayed the same segmentation fault.

Note: For my original post uses a machine which runs python version 3.6.8 and CRPropa 3.1.7 on a Linux system.

QuinFin commented 1 year ago

It is running for me fine when going back to python version 2.7 (including downgrading all required packages to a compatible version), but I get the segfault error with all python3 versions that I have tried. This is a little suboptimal as python version 2.7 days have more or less past by now.

lukasmerten commented 5 months ago

@JassimarS and @QuinFin, should be fixed with #466. Please reopen if the problem persists.