Segmentation Fault 139 (core dumped) C++ after changing operating system

Question

I had a working program to analyse data in C++ that has produced something like 35 successful data files to date. I was working on Scientific Linux in Code:Blocks when it was working and barring some small errors involving very large grid sizes (1000x1000+) it worked perfectly and produced exactly what I was looking for.

I recently switched to Ubuntu and expected it to work fine and it doesn't. It accepts initial input (the first particle switch) but then crashes immediately with a segmentation fault 139. I have tried to run it in Windows instead with my dual boot but it doesn't seem to recognise the local filing system so I'm forced to ask for help.

It's a long program so I'll reproduce the whole thing. I apologise in advance.

// This program converts the column output of a 1D PIC code into a workable solution

#include <iostream>
#include <fstream>
#include <math.h>
using namespace std;

double calculateXMaximum(double arraymax[], int size)
{
    double maximum = 0;
    for (int k = 1; k < size/2; k++)
    {
        if(arraymax[2*k] > maximum)
        {
            maximum = arraymax[2*k];
        }
    }
    return maximum;
}

double calculateXMinimum(double arraymin[], int size)
{
    double minimum = 0;
    for (int k = 1; k < size/2; k++)
    {
        if(arraymin[2*k] < minimum)
        {
            minimum = arraymin[2*k];
        }
    }
    return minimum;
}

double calculatePXMaximum(double arraymax[], int size)
{
    double maximum = 0;
    for (int k = 1; k < size/2; k++)
    {
        if(arraymax[2*k+1] > maximum)
        {
            maximum = arraymax[2*k+1];
        }
    }
    return maximum;
}

double calculatePXMinimum(double arraymin[], int size)
{
    double minimum = 0;
    for (int k = 1; k < size/2; k++)
    {
        if(arraymin[2*k+1] < minimum)
        {
            minimum = arraymin[2*k+1];
        }
    }
    return minimum;
}

int main()
{
// Variables settable before running program - will set up initialisation later.
double xmin = 0;
double xmax = 0;
double pmin = 0;
double pmax = 0;

int xni = 0;
double xntemp = 0;
double deltax = 0;
int xi; // X interpolates, defined from console for resolution of diagram

int pnj = 0;
double pntemp = 0;
double deltap = 0;
int pi;
int type;
double modifier;

// Determines momentum modifier!

cout << "For particle type, enter 1 (e-) or 2 (p+)" << endl;
cout << "Particle type:  ";
cin >> type;

if (type == 2)
{
    modifier = 1836;
}
else
{
    modifier = 1;
}

ifstream inputFile;
ofstream outputFile;

inputFile.open ("/home/Nick/fi020000.dat");
outputFile.open ("/home/Nick/fi20k.dat");

int dataformat[2];
for(int rd = 0; rd < 2; rd++)
{
    dataformat[rd] = 0;
    inputFile >> dataformat[rd];
}

int records = dataformat[1] + 2;
double data[records];
cout << "Number of particles: " << dataformat[1]/2 << endl;


// Introduction of data from input data file loop.  Produces records.
for (int count = 0; count < records; count++)
{
    inputFile >> data[count];
}

// Calling functions for xmin and xmax.  May streamline later

xmax = calculateXMaximum(data, records) * 1.1;
cout << "Maximum x value: " << xmax << endl;
xmin = calculateXMinimum(data, records) * 1.1;
cout << "Minimum x value: " << xmin << endl;
pmax = calculatePXMaximum(data, records) * 1.1 / modifier;
cout << "Maximum p value: " << pmax << endl;
pmin = calculatePXMinimum(data, records) * 1.1 / modifier;
cout << "Minimum p value: " << pmin << endl;

// Definition of bin size

cout << "Entire desired number of x bins: ";
cin >> xi;
const int xip = xi;
cout << "Enter desired number of p bins: ";
cin >> pi;
const int pip = pi;
cout << "Grid is " << xip << " x " << pip << endl;

// Calculate DELTA X and DELTA P
deltax = (xmax - xmin)/(xip);
deltap = (pmax - pmin)/(pip);

cout << "Resolution of x:  " << deltax << endl;
cout << "Resolution of p:  " << deltap << endl;

int phaseSpace [xip][pip];
for(int i=0; i<xip; i++)
{
    for(int j=0; j<pip; j++)
    {
        phaseSpace[i][j] = 0;
    }
}

for (int phasecount=1; phasecount < (records/2)-1; phasecount++)
{
    xntemp = (data[2*phasecount] - xmin)/deltax;
    xni = floor(xntemp);
    pntemp = ((data[(2*phasecount)+1] / modifier) - pmin)/deltap;
    pnj = floor(pntemp);
    phaseSpace[xni][pnj] = phaseSpace[xni][pnj] + 1;
}

for (int xoutcount = 0; xoutcount < xip; xoutcount++)
{
    for (int poutcount = 0; poutcount < pip; poutcount++)
    {
        outputFile << xmin+((xoutcount+0.5)*deltax) << " " << pmin+((poutcount+0.5)*deltap) << " "<< phaseSpace[xoutcount][poutcount] << endl;
    }
    outputFile << endl;
}


cout << "Program complete" << endl;

return 0;

}

My intent was to finish a 30 page report by this weekend and now the program that I was using to do this has completely fallen apart. I am not a computer scientist - I'm a physicist and I learned C++ less than a month ago. As such, I have no clue what's going on. I know this topic has been seen a lot but I can't really understand the advice.

EDIT: Stack trace is:

#0 0x4010af     ??     ()     (??:??)
#1 0x7ffff7215ea5       __libc_start_main() (/lib/x86_64-linux-gnu/libc.so.6:??)
#2 0x4017f1 ??    () (??:??)

EDIT2: Valgrind results

==4089== Using Valgrind-3.8.1 and LibVEX; rerun with -h for copyright info
==4089== Command: ./Analysis
==4089== 
For particle type, enter 1 (e-) or 2 (p+)
Particle type:  2
==4089== Warning: client switching stacks?  SP change: 0x7fefff9c0 --> 0x7fe6d7118
==4089==          to suppress, use: --max-stackframe=9603240 or greater
==4089== Invalid write of size 8
==4089==    at 0x4010AF: ??? (in /home/paladin/Contour/Analysis/bin/Release/Analysis)
==4089==    by 0x5673EA4: (below main) (libc-start.c:260)
==4089==  Address 0x7fe6d7118 is on thread 1's stack
==4089== 
==4089== 
==4089== Process terminating with default action of signal 11 (SIGSEGV)
==4089==  Access not within mapped region at address 0x7FE6D7118
==4089==    at 0x4010AF: ??? (in /home/paladin/Contour/Analysis/bin/Release/Analysis)
==4089==  If you believe this happened as a result of a stack
==4089==  overflow in your program's main thread (unlikely but
==4089==  possible), you can try to increase the size of the
==4089==  main thread stack using the --main-stacksize= flag.
==4089==  The main thread stack size used in this run was 8388608.
==4089== 
==4089== Process terminating with default action of signal 11 (SIGSEGV)
==4089==  Access not within mapped region at address 0x7FE6D7111
==4089==    at 0x4A256A0: _vgnU_freeres (in /usr/lib/valgrind/vgpreload_core-amd64-linux.so)
==4089==  If you believe this happened as a result of a stack
==4089==  overflow in your program's main thread (unlikely but
==4089==  possible), you can try to increase the size of the
==4089==  main thread stack using the --main-stacksize= flag.
==4089==  The main thread stack size used in this run was 8388608.
==4089== 
==4089== HEAP SUMMARY:
==4089==     in use at exit: 17,520 bytes in 4 blocks
==4089==   total heap usage: 4 allocs, 0 frees, 17,520 bytes allocated
==4089== 
==4089== LEAK SUMMARY:
==4089==    definitely lost: 0 bytes in 0 blocks
==4089==    indirectly lost: 0 bytes in 0 blocks
==4089==      possibly lost: 0 bytes in 0 blocks
==4089==    still reachable: 17,520 bytes in 4 blocks
==4089==         suppressed: 0 bytes in 0 blocks
==4089== Rerun with --leak-check=full to see details of leaked memory
==4089== 
==4089== For counts of detected and suppressed errors, rerun with: -v
==4089== ERROR SUMMARY: 1 errors from 1 contexts (suppressed: 2 from 2)
Segmentation fault (core dumped)

The fault occurs at the ifstream inputFile statement.

EDIT3: As requested, console session:

paladin@paladin:~/Programming/Contour Constructor 2/bin/Debug$ ./Contour\ Constructor\ 2 
For particle type, enter 1 (e-) or 2 (p+)
Particle type:  2
Segmentation fault (core dumped)

There are 1200402 lines in the input file, corresponding to 600200 particles in the PIC code plus 2 descriptive lines.

EDIT4: Complete shot in the dark here but I originally compiled this under GCC 4.4.7 on Scientific Linux. I am now using the latest Ubuntu version (4.8.1). Has something changed in the interim which would invalidate the size of the file I'm using?

Answer 1

The big hint here is that valgrind thinks you're switching threads -- and certainly from the code you've shown, you're not using multiple threads. Since every thread in a multithreaded program has its own stack, valgrind assumes that if the stack pointer changes by more than a certain threshold (see --max-stackframe as mentioned in the valgrind output) that you're switching to a different thread.

In reality, what's happening is you've created a HUGE stack frame, 9603240 bytes to be exact. This exceeds the 8MB or so that you probably are getting by default. You can see the current soft limit by looking in a shell:

$ ulimit -s
8192

In other words, the stack is limited to 8MB. If you exceed that limit, Linux will assume something bad(tm) happened and terminate your process.

One immediate fix you can use is to raise the limit, or just set it to unlimited:

$ ulimit -s 16384 # 16MB stack
$ ulimit -s unlimited # "unlimited" stack

This should stop your process from crashing, and explains why it worked fine on one box and not another (the one where it worked probably had a higher stack size limit set by default).

Now, from a design perspective, it's usually a bad idea to create stack frames this large. For large allocations, you should use heap memory, so something like:

double data[records];

could be replaced by

double *data = new double[records];

// use data ...

delete[] data;

which allocates and frees memory from the heap instead. This will allow you to avoid such problems in the first place. Or, you can always use one of the standard containers like std::vector<>` instead to avoid the issue.

Answer 2

If you see this problem only for large input sizes, chances are you get a stack overflow when creating phaseSpace, since that might be a quite big array. If you use a std::vector instead of the plain array, such a problem should be avoided:

#include <vector>

// ...

// A vector containing vectors containing integers.
// Initialized to contain the appropriate amount of space filled with zeros
std::vector< std::vector<int> > phaseSpace(pip, std::vector<int>(xip, 0));

The rest of your code probably works just the same with the vector.

---

Alternatively, if xni and pni go out of bounds, you will overwrite random memory. You could add output statements to show these values and see if something looks wrong.

Answer 3

I am surprised this program even compiles: you are declaring arrays of size that cannot be determined during compile time, such as below. Could you tell which compiler you are using?

int records = dataformat[1] + 2;
double data[records];

Answer 4

This is what it might look like using std::vector instead of native arrays.

// This program converts the column output of a 1D PIC code into a
// workable solution

#include <cmath>
#include <fstream>
#include <iostream>
#include <limits>
#include <vector>

double calculateXMaximum(const std::vector<double>& arraymax) {
  double maximum = -std::numeric_limits<double>::max();
  for (unsigned k = 1; k < arraymax.size()/2; ++k) {
    maximum = std::max(arraymax[2*k], maximum);
  }
  return maximum;
}

double calculateXMinimum(const std::vector<double>& arraymin) {
  double minimum = std::numeric_limits<double>::max();
  for (unsigned k = 1; k < arraymin.size()/2; ++k) {
    minimum = std::min(arraymin[2*k], minimum);
  }
  return minimum;
}

double calculatePXMaximum(const std::vector<double>& arraymax) {
  double maximum = -std::numeric_limits<double>::max();
  for (unsigned k = 1; k < arraymax.size()/2; ++k) {
    maximum = std::max(arraymax[2*k+1], maximum);
  }
  return maximum;
}

double calculatePXMinimum(const std::vector<double>& arraymin) {
  double minimum = std::numeric_limits<double>::max();
  for (unsigned k = 1; k < arraymin.size()/2; ++k) {
    minimum = std::min(arraymin[2*k+1], minimum);
  }
  return minimum;
}

int main()
{
  // Variables settable before running program
  //  - will set up initialisation later.

  int xni = 0;
  double xntemp = 0;
  int xi; // X interpolates, defined from console for resolution of diagram

  int pnj = 0;
  double pntemp = 0;
  int pi;
  int type;

  // Determines momentum modifier!

  std::cout << "For particle type, enter 1 (e-) or 2 (p+)\n";
  std::cout << "Particle type:  ";
  std::cin >> type;

  const double modifier = (type == 2) ? 1836.0 : 1.0;

  std::ifstream inputFile("fi020000.dat");
  std::ofstream outputFile("fi20k.dat");

  int dataformat[2];
  inputFile >> dataformat[0];
  inputFile >> dataformat[1];

  int records = dataformat[1] + 2;
  std::vector<double> data(records, 0.0);
  std::cout << "Number of particles: " << dataformat[1]/2 << std::endl;


  // Introduction of data from input data file loop.  Produces records.
  for (int count = 0; count < records; ++count) {
    inputFile >> data[count];
  }

  // Calling functions for xmin and xmax.  May streamline later

  const double xmax = calculateXMaximum(data) * 1.1;
  std::cout << "Maximum x value: " << xmax << std::endl;
  const double xmin = calculateXMinimum(data) * 1.1;
  std::cout << "Minimum x value: " << xmin << std::endl;
  const double pmax = calculatePXMaximum(data) * 1.1 / modifier;
  std::cout << "Maximum p value: " << pmax << std::endl;
  const double pmin = calculatePXMinimum(data) * 1.1 / modifier;
  std::cout << "Minimum p value: " << pmin << std::endl;

  // Definition of bin size

  std::cout << "Entire desired number of x bins: ";
  std::cin >> xi;
  const int xip = xi;
  std::cout << "Enter desired number of p bins: ";
  std::cin >> pi;
  const int pip = pi;
  std::cout << "Grid is " << xip << " x " << pip << std::endl;

  // Calculate DELTA X and DELTA P
  const double deltax = (xmax - xmin)/(xip);
  const double deltap = (pmax - pmin)/(pip);

  std::cout << "Resolution of x:  " << deltax << std::endl;
  std::cout << "Resolution of p:  " << deltap << std::endl;

  std::vector< std::vector<int> > phaseSpace(xip, std::vector<int>(pip, 0));

  for (int phasecount=1; phasecount < (records/2)-1; ++phasecount) {
    xntemp = (data[2*phasecount] - xmin)/deltax;
    xni = std::floor(xntemp);
    pntemp = ((data[(2*phasecount)+1] / modifier) - pmin)/deltap;
    pnj = std::floor(pntemp);
    phaseSpace[xni][pnj] = phaseSpace[xni][pnj] + 1;
  }

  for (int xoutcount = 0; xoutcount < xip; ++xoutcount) {
    for (int poutcount = 0; poutcount < pip; ++poutcount) {
      outputFile << xmin+((xoutcount+0.5)*deltax) << " " << pmin+((poutcount+0.5)*deltap) << " " << phaseSpace[xoutcount][poutcount] << "\n";
    }
    outputFile << "\n";
  }

  std::cout << "Program complete" << std::endl;

  return 0;
}

I removed some of the std::endl calls and replaced them with "\n" to avoid unnecessary I/O buffer flushing.

I made a slight change to the way your min's and max's are calculated.

I moved a few variable declarations and made them const where possible.

Important: You should be checking the state of your ifstream as you are reading input. As it stands, you could read past the end of the file and never know it. There also appears to be a misunderstanding of the contents of your data array. It seems as though you think it contains your two leading descriptive lines. However you have already pulled those two numbers from the stream into dataformat, so you do not need to account for them.

A very small set of hand verified test input and output data would be useful.