Завершение процесса (код ошибки -1073741571)

Поэтому я работаю над этой симуляцией, для которой я должен написать код.
Первоначально я получал код ошибки 255 (0xff), который, как я думал, мог быть связан с добавлением Eigen в мои блоки кода. Однако, как только я переместил свои собственные Eigen-заголовки в папку include блоков кода, я теперь получаю новую ошибку с сообщением:

Процесс завершен со статусом -1073741571 (0 минут, 2 секунд)

Я посмотрел эту ошибку, и она говорит, что она возникает, когда ошибка вызывает переполнение стека. Я запустил мой отладчик и получил сообщение:

Зарегистрирован новый тип: wxString

Зарегистрирован новый тип: STL String

Зарегистрирован новый тип: STL Vector

Установка точек останова Имя и версия отладчика: GNU GDB (GDB) 7.6.1

PID дочернего процесса: 13284

Программа получила сигнал SIGSEGV, Ошибка сегментации.

Таким образом, я знал, что есть ошибка с одной из моих функций.
Вот мой код (он довольно большой)

#include <iostream>
#include <fstream>
#include <conio.h>
#include <random>
#include<math.h>
#include <Eigen/Dense>
#include<time.h>
#define pi atan(1) * 4
using namespace std;
using namespace Eigen;
int Poisson(int N)
{
int n=1;
float x,y,k=0;
do
{
x=(float)rand()/RAND_MAX;
y = -log(x);
k=k+y;
n=n+1;     //n increases for each random number generated
}
while(k<N);
return n;
}
int sumrow(MatrixXd A,int row, int col_len)
{
int sum=0;
for(int i=0;i<col_len;i++)
{
sum=sum+A(row,i);
}
return sum;
}
void generate_D(MatrixXd A, int dim, MatrixXd D)
{
for(int x=0;x<dim;x++)
{
for(int y=0;y<dim;y++)
{
if(x==y)
{
D(x,y)=1/sqrt(sumrow(A,x,y));
}
else
{ D(x,y) =0;  }
}
}
}
void sqmatrix_multiply(MatrixXd A,MatrixXd B,int dim, MatrixXd C)
{
for(int i=0;i<dim;i++)
{
for(int j=0;j<dim;j++)
{
for(int k=0;k<dim;k++)
{
C(i,j) = C(i,j) + (A(i,k)*B(k,j));
}
}
}
}
void generate_B(MatrixXd A, int dim, MatrixXd B)
{
MatrixXd D;
MatrixXd C, inter1, inter2;
generate_D(A,dim,D);
sqmatrix_multiply(A,D,dim,inter1);
sqmatrix_multiply(D,C,dim,inter2);
for(int i=0;i<dim;i++)
{
for(int j=0;j<dim;j++)
{
if(i==j)
{
B(i,j) = 1 - inter2(i,j);
}
else
{
B(i,j) = 0 - inter2(i,j);
}

}
}
}
int main()
{
srand(time(NULL));
double X,Y,x,y,a,b,L[5000][5000],M[5000][5000];
MatrixXd A, B;
for(int i=0; i<5000;i++)
{   L[i][i]=0;  M[i][i]=0;   }
int N, I;
cout<<"Enter your number and iterations";
cin>>N>>I;
int P[5000];
for(int i=0; i<I; i++)
{
P[i]=Poisson(N);
}
for(int k=0; k<I; k++)
{
for (int j=0;j<P[k];j++)
{
A:
X=(float)rand()/RAND_MAX;
Y=(float)rand()/RAND_MAX;
x=sqrt(X);
y = 2*(pi)*Y;
a=x*cos(y);
b=x*sin(y);
if((a*a)+(b*b)>1)
{
goto A;
}
L[k][j] = a; M[k][j] = b;
}
}
for(int a=0;a<I;a++)
{
for(int a1=0;a1<P[a];a1++)
{
for(int a2=0;a2<P[a];a2++)
{
if(a1==a2)
{
A(a1,a2)=0;
}
else
{
if(((L[a][a1]-L[a][a2])*(L[a][a1]-L[a][a2]))+((M[a][a1]-M[a][a2])*M[a][a1]-M[a][a2]) < 1)  //distance formula
{
A(a1,a2)=1;
}
else
{
A(a1,a2)=0;
}
}
}
}
generate_B(A,P[a],B);
EigenSolver<MatrixXd> b(B,false);
cout<<a+1<<" Eigen values are: \n "<<b.eigenvalues()<<"\n";
}
return 0;
}

Я опробовал несколько функций по отдельности в разных программах. Проблема не связана с заголовками Eigen, так как sumrow () работал отлично.
Однако функция Пуассона () выдала мне ту же ошибку (-1073741571). Я попытался заменить do-while простым, но это не сработало.

Ранее, когда я использовал Poisson () отдельно, он работал нормально. Но теперь, это не работает независимо друг от друга.

Извините за длинный вопрос, но я был на нем уже больше недели и не могу понять это. Если бы вы могли сказать мне, что, возможно, не так с моим кодом, я был бы очень благодарен.
Заранее большое спасибо!