Segmentation fault when adding a PDE to a crypt cell population

[Kgeorge1729]

I've recently been experimenting with adding an Elliptic PDE to a simulation of an intestinal crypt. I've modified the code detailed in, https://chaste.cs.ox.ac.uk/trac/wiki/PaperTutorials/CryptProliferationDistribution. More specifically, I've modified the file TestCryptProliferationDistributionLiteratePaper.hpp. I've hard coded the simulation to only use one set of parameters:

unsigned cell_proliferation_model =2;

bool contact_inhibition = true;

bool wnt_dependent_ccd = false;
	// Sim and sweep Params
	// TS_ASSERT(CommandLineArguments::Instance()->OptionExists("-end_time"));
	// double end_time = atof(CommandLineArguments::Instance()->GetStringCorrespondingToOption("-end_time").c_str());
double end_time = 100;
	// TS_ASSERT(CommandLineArguments::Instance()->OptionExists("-min_CI"));
//   double min_CI_param = atof(CommandLineArguments::Instance()->GetStringCorrespondingToOption("-min_CI").c_str());
double min_CI_param = 0;
  // TS_ASSERT(CommandLineArguments::Instance()->OptionExists("-max_CI"));
  // double max_CI_param = atof(CommandLineArguments::Instance()->GetStringCorrespondingToOption("-max_CI").c_str());
double max_CI_param = 1;
	// TS_ASSERT(CommandLineArguments::Instance()->OptionExists("-num_CI_sweeps"));
	// double num_CI_sweeps = atof(CommandLineArguments::Instance()->GetStringCorrespondingToOption("-num_CI_sweeps").c
double num_CI_sweeps = 10;
  // Crypt Setup
  double cell_radius = 3.5;
  double crypt_length = 140; //70
  double crypt_radius = 16.0/M_PI*6.0;                 // Choosing same dimensions as for halted migration paper
  // For this size domain there are about 75 cells in the bottom hemisphere so this makes about 20% Paneth cells
	unsigned num_paneth_cells = 15;//15;
	unsigned num_stem_cells = 60;
	unsigned num_cells = num_paneth_cells + num_stem_cells;
  double stem_retainer_force_magnitude = 7.5*10;
  double paneth_retainer_force_magnitude = 7.5*10;
  double param =1.0;
	// Extra code to calculate the generation if Model 1
	unsigned upper_max_transit_generations = (unsigned)ceil(param);
	unsigned lower_max_transit_generations = (unsigned)floor(param);
	double prob_of_upper_max_transit_generations = param - floor(param);
	//Sweep over CI param
  unsigned CIindex = 10;
  double CIparam = min_CI_param + (double)CIindex / double(num_CI_sweeps) * (max_CI_param-min_CI_param);`

I believe this is the only relevant code for defining the parameters of the model. If I run the code with these parameters, the simulation works and does completes without any segmentation fault occurring. 

When I add the following code to the simulation:

`// Create PDE and boundary condition objects

MAKE_PTR_ARGS(CellwiseSourceEllipticPde<3>, p_pde, (crypt, 0.03));

MAKE_PTR_ARGS(ConstBoundaryCondition<3>, p_bc, (0.0));

MAKE_PTR_ARGS(EllipticGrowingDomainPdeModifier<3>, p_pde_modifier, (p_pde, p_bc, false));

p_pde_modifier->SetDependentVariableName("hhog");

simulator.AddSimulationModifier(p_pde_modifier);

I get the following error:

 Entering Test3DCrypt
547: [0]PETSC ERROR: ------------------------------------------------------------------------
547: [0]PETSC ERROR: Caught signal number 11 SEGV: Segmentation Violation, probably memory access out of range
547: [0]PETSC ERROR: Try option -start_in_debugger or -on_error_attach_debugger
547: [0]PETSC ERROR: or see https://petsc.org/release/faq/#valgrind
547: [0]PETSC ERROR: or try http://valgrind.org on GNU/linux and Apple Mac OS X to find memory corruption errors
547: [0]PETSC ERROR: configure using --with-debugging=yes, recompile, link, and run 
547: [0]PETSC ERROR: to get more information on the crash.
547: [0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------
547: [0]PETSC ERROR: Signal received
547: [0]PETSC ERROR: See https://petsc.org/release/faq/ for trouble shooting.
547: [0]PETSC ERROR: Petsc Release Version 3.15.5, Sep 29, 2021 
547: [0]PETSC ERROR: /home/chaste/lib/projects/CryptProliferationTutorial/test/TestCryptBuild on a  named a231b79934d2 by chaste Fri Sep 15 15:24:09 2023
547: [0]PETSC ERROR: Configure options --build=x86_64-linux-gnu --prefix=/usr --includedir=${prefix}/include --mandir=${prefix}/share/man --infodir=${prefix}/share/info --sysconfdir=/etc --localstatedir=/var --with-option-checking=0 --with-silent-rules=0 --libdir=${prefix}/lib/x86_64-linux-gnu --runstatedir=/run --with-maintainer-mode=0 --with-dependency-tracking=0 --with-debugging=0 --shared-library-extension=_real --with-shared-libraries --with-pic=1 --with-cc=mpicc --with-cxx=mpicxx --with-fc=mpif90 --with-cxx-dialect=C++11 --with-opencl=1 --with-blas-lib=-lblas --with-lapack-lib=-llapack --with-scalapack=1 --with-scalapack-lib=-lscalapack-openmpi --with-ptscotch=1 --with-ptscotch-include=/usr/include/scotch --with-ptscotch-lib="-lptesmumps -lptscotch -lptscotcherr" --with-fftw=1 --with-fftw-include="[]" --with-fftw-lib="-lfftw3 -lfftw3_mpi" --with-yaml=1 --with-hdf5-include=/usr/include/hdf5/openmpi --with-hdf5-lib="-L/usr/lib/x86_64-linux-gnu/hdf5/openmpi -L/usr/lib/x86_64-linux-gnu/openmpi/lib -lhdf5 -lmpi" --CXX_LINKER_FLAGS=-Wl,--no-as-needed --with-hypre=1 --with-hypre-include=/usr/include/hypre --with-hypre-lib=-lHYPRE --with-mumps=1 --with-mumps-include="[]" --with-mumps-lib="-ldmumps -lzmumps -lsmumps -lcmumps -lmumps_common -lpord" --with-suitesparse=1 --with-suitesparse-include=/usr/include/suitesparse --with-suitesparse-lib="-lumfpack -lamd -lcholmod -lklu" --with-superlu=1 --with-superlu-include=/usr/include/superlu --with-superlu-lib=-lsuperlu --with-superlu_dist=1 --with-superlu_dist-include=/usr/include/superlu-dist --with-superlu_dist-lib=-lsuperlu_dist --with-ml=1 --with-ml-include=/usr/include/trilinos --with-ml-lib=-ltrilinos_ml --with-zoltan=1 --with-zoltan-include=/usr/include/trilinos --with-zoltan-lib=-ltrilinos_zoltan --prefix=/usr/lib/petscdir/petsc3.15/x86_64-linux-gnu-real --PETSC_ARCH=x86_64-linux-gnu-real CFLAGS="-g -O2 -ffile-prefix-map=/build/petsc-tFq6Xk/petsc-3.15.5+dfsg1=. -flto=auto -ffat-lto-objects -flto=auto -ffat-lto-objects -fstack-protector-strong -Wformat -Werror=format-security -fPIC" CXXFLAGS="-g -O2 -ffile-prefix-map=/build/petsc-tFq6Xk/petsc-3.15.5+dfsg1=. -flto=auto -ffat-lto-objects -flto=auto -ffat-lto-objects -fstack-protector-strong -Wformat -Werror=format-security -fPIC" FCFLAGS="-g -O2 -ffile-prefix-map=/build/petsc-tFq6Xk/petsc-3.15.5+dfsg1=. -flto=auto -ffat-lto-objects -flto=auto -ffat-lto-objects -fstack-protector-strong -fPIC -ffree-line-length-0" FFLAGS="-g -O2 -ffile-prefix-map=/build/petsc-tFq6Xk/petsc-3.15.5+dfsg1=. -flto=auto -ffat-lto-objects -flto=auto -ffat-lto-objects -fstack-protector-strong -fPIC -ffree-line-length-0" CPPFLAGS="-Wdate-time -D_FORTIFY_SOURCE=2" LDFLAGS="-Wl,-Bsymbolic-functions -flto=auto -ffat-lto-objects -flto=auto -Wl,-z,relro -fPIC" MAKEFLAGS=w
547: [0]PETSC ERROR: #1 User provided function() at unknown file:0
547: [0]PETSC ERROR: Run with -malloc_debug to check if memory corruption is causing the crash.
547: --------------------------------------------------------------------------
547: MPI_ABORT was invoked on rank 0 in communicator MPI_COMM_WORLD
547: with errorcode 59.
547: 
547: NOTE: invoking MPI_ABORT causes Open MPI to kill all MPI processes.
547: You may or may not see output from other processes, depending on
547: exactly when Open MPI kills them.
547: --------------------------------------------------------------------------
1/1 Test #547: TestCryptBuild ...................***Failed    1.76 sec

0% tests passed, 1 tests failed out of 1

Is there something about how the crypt cell population is defined which prevents the addition of a PDE to the simulation? Any help or suggestions with this problem are welcome.

[AlexFletcher]

Hi George,

@jmosborne may be better placed to comment on this as an author on that paper, but just a quick question to help understand what you're trying to do - is your intention to solve an elliptic PDE numerically on the mesh generated by NodeBasedCellPopulation<3>::GetTetrahedralMeshForPdeModifier(), i.e. the mesh with 3D elements comprising the crypt 'volume', or is your intention to solve an elliptic PDE numerically on a mesh with 2D elements comprising by the crypt 'surface'?

Best wishes,
Alex

[Kgeorge1729]

Hi Alex,

Thanks for your response. The intention is to solve an elliptic PDE numerically on the mesh generated by a Node based cell population. I'm not entirely sure if I'm using the function GetTetrahedralMeshForPdeModifier.

Perhaps it might be good for me to provide a little more detail. I first modified the following UserTutorial, https://chaste.cs.ox.ac.uk/trac/wiki/UserTutorials/RunningTumourSpheroidSimulations, to use a node based cell population as opposed to a mesh based cell population. To do so, all I needed to do was replace the honeycomb mesh with a mesh based on the nodes using the NodesOnlyMesh class. This ran fine and didn't cause a segmentation fault. It was even possible to convert this simulation into a 3d simulation.

I don't know if this means that Chaste then uses GetTetrahedralMeshForPdeMoifier() under the hood but I haven't explicitly used it anywhere. The intention is to solve a PDE on a mesh which occupies the same volume as the crypt.

I hope that clarifies what I've done a little better. I'll have a look again and see if Chaste is using GetTetrahedralMeshForPdeMoifier().

Thanks for you help

Kind regards,

George

[Kgeorge1729]

I have an update on this issue. I've managed to stop the segmentation from occurring but I'm not sure why my fix works. The cells of the crypt are initially in a ring, if I run the simulation without the PDE modifier they form a crypt shape. I've found that applying the PDE modifier after the cells are in a crypt shape prevents the segmentation fault. I'm happy to close this discussion, but does anyone have a suggestion on why this solves the issue?

[jmosborne]

So the issue is that initially the cells are in a plane so there is no 3d object to solve the PDE on wherea once cells lie on the surface you can make a 3d mesh. Whether that 3d mesh is the one you want is another matter.

From the code it looks like you are solving the PDE on the domain enclosed by the surface defined by the nodes which lie in the shape of a crypt. The mesh from this example will not be a good mesh to solve a pde on as thee will be very poor aspect ratios. The GrowingDomain pde examples are meant for domains with cells evenly distributed among them so a monlayer in 2d or a spheroid in 3d. It wont work (it may compile but the results will be junk) for a sparse setup like the crytp example you're using.

If you want to solve a 3d pde in the domain contain the cells you should consider using the FixedDomain Pdes