I have few very large arrays which I want to use to store some calculations. The arrays are of size 10000003,25003*3, 1000000,1000000.I have used malloc to dynamically allocate memory to all the arrays before using any of them.

The code tries to implement an algorithm for a discritized space. I have tried using smaller sizes for my arrays and everything seems to be working fine.

#include<stdio.h>
#include<math.h>
#include<stdlib.h>
double mu_a=0.0;
double mu_b=0.0;
int n_b_d1=50.0;
int n_b_d2=50.0;



double dot(int m,int n,int o,int p,double * array[]){
    m=m-1;
    n=n-1;
    o=o-1;
    p=p-1;
    return ((array[m][0]-array[n][0])*(array[o][0]-array[p][0]))+((array[m][1]-array[n][1])*(array[o][1]-array[p][1]))+((array[m][2]-array[n][2])*(array[o][2]-array[p][2]));
    }



double * linspace(double x1, double x2,int n,double array[]){
    double * x=calloc(n,sizeof(double));
    double step=(x2-x1)/(double)(n-1);
    for (int i=0;i<n;i++){
        array[i]=x1+((double)i*step);
    x[i]=x1+((double)i*step);
    
    }
    return x;
}

/*array length=*(&array+1)-array)*/
void points(double P_1[],double P_2[],double P_3[],double P_4[],double temp[]){
    double dst=400*(((P_3[0]-P_1[0])*((P_2[1]-P_1[1])-(P_4[1]-P_3[1])))-((P_3[1]-P_1[1])*((P_2[0]-P_1[0])-(P_4[0]-P_3[0]))));
    if (dst>0.001){
        temp[3]=999999; return;
        }
    double * array[4]={P_1,P_2,P_3,P_4};  
    mu_a=((dot(1,3,4,3,array)*dot(4,3,2,1,array))-(dot(1,3,2,1,array)*dot(4,3,4,3,array)))/((dot(2,1,2,1,array)*dot(4,3,4,3,array))-(dot(4,3,2,1,array)*dot(4,3,2,1,array)));
    mu_b=(dot(1,3,4,3,array)+(mu_a*dot(4,3,2,1,array)))/dot(4,3,4,3,array);
double P_a[3]={0,0,0};
    double P_b[3]={0,0,0};
    for (int i=0;i<3;i++){
        P_a[i]=P_1[i]+mu_a*(P_2[i]-P_1[i]);
        P_b[i]=P_3[i]+mu_b*(P_4[i]-P_3[i]);
    }
    //printf("%lf,%lf,%lf,n",P_a[0],P_a[1],P_a[2]);
    dst=pow((pow((P_a[0]-P_b[0]),(double)(2))+pow((P_a[1]-P_b[1]),(double)(2))+pow(((P_a[2]-P_b[2])/10),(double)(2))),(0.5));
    //printf("%lf",dst);
    temp[0]=(P_a[0]+P_b[0])/2;
    temp[1]=(P_a[1]+P_b[1])/2;
    temp[2]=(P_a[2]+P_b[2])/2;
    temp[3]=dst;
}




double omega(double a,double b,double d){
    double alpha=a/(double)2/d;
    double beta=b/(double)2/d;
    return (a*b/pow(d,(double)2));
}







double omega_offaxis(double A,double B,double a, double b,double d){
    return (omega(2*(A+a),2*(B+b),d)-omega(2*A,2*(B+b),d)-omega(2*(A+a),2*B,d)+omega(2*A,2*B,d))/(double)4;
}






double calc_prob(double b[],double d1[][3],double d2[][3]){
    double d[2][3];
    double dist;
        if (b[2]<=200){
        /*#tot=omega_4quad(b[0],b[1],50,50,b[2])*/
            dist=400.0-b[2];

        for(int i=0;i<2;i++){
            for (int j=0;j<3;j++){
                d[i][j]=d2[i][j];
            }
        }
    }
    else{
        /*#tot=omega_4quad(b[0],b[1],50,50,400-b[2])*/
        
        //double d[2][3];
    for(int i=0;i<2;i++){
        for (int j=0;j<3;j++){
            d[i][j]=d1[i][j];
        }
    }   
        dist=b[2];
    }    
    double x1=d1[0][0],y1=d1[0][1],z1=d1[0][2];
    double x2=d1[1][0],y2=d1[1][1],z2=d1[1][2];
    double x3=d2[0][0],y3=d2[0][1],z3=d2[0][2];
    double x4=d2[1][0],y4=d2[1][1],z4=d2[1][2];

    double x_c_1=((x3-x1)*b[2]/(double)(400))+x1;
    double x_c_2=((x4-x2)*b[2]/(double)(400))+x2;
    double y_c_1=((y3-y1)*b[2]/(double)(400))+y1;
    double y_c_2=((y4-y2)*b[2]/(double)(400))+y2;
    

    if ((b[0]<x_c_1) || (b[0]>x_c_2) || (b[1]<y_c_1) || (b[1]>y_c_2)){
        return 0;
    }
    double x_a=d[0][0],y_a=d[0][1],z_a=d[0][2];
    double x_b=d[1][0],y_b=d[1][1],z_b=d[1][2];
    double a1=10;
    double b1=10;
    return (omega_offaxis(x_a-b[0],y_a-b[1],a1, b1, dist)/(double)4*M_PI);
}
//double blocks[1000000][3];
void main(){
double ** blocks=(double **)malloc(1000000*sizeof(double *));
for(int i=0;i<1000000;i++)
{
    blocks[i]=(double*)malloc(3*sizeof(double));
}
double ** intersect=(double **)malloc(100000*sizeof(double *));
for(int i=0;i<100000;i++)
{
    intersect[i]=(double*)malloc(3*sizeof(double));
}
//double detectors[2500][3][3]
double *** detectors=(double ***)malloc(2500*sizeof(double **));
for(int i=0;i<2500;i++)
{
    detectors[i]=(double**)malloc(3*sizeof(double*));

for(int j=0;j<3;j++){detectors[i][j]=(double*)malloc(3*sizeof(double));
}}
printf("here1");
printf("%lf",M_PI);
//double blocking[6][6][41];
double x_blocking[51];
double y_blocking[51];
double z_blocking[401];
linspace(0.0,50.0,51,x_blocking);
linspace(0.0,50.0,51,y_blocking);
linspace(0.0,400.0,401,z_blocking);
//printf("%f,%f,%f",blocking[1][4][5],blocking[3][4][5],blocking[1][2][3]);
FILE * filo;
filo=fopen("input.txt","r");
double data[52822][6];
if (filo==NULL){
    printf("%s","noice");
}
printf("%f",filo);
for(int i=0;i<52822;i++){
    //double y=200;
    fscanf(filo,"%lf %lf %lf %lfn",&data[i][0],&data[i][1],&data[i][3],&data[i][4]);
    data[i][2]=0;
    data[i][5]=400;
    
}
printf("here2");
FILE * gilo;

gilo=fopen("output.txt","w");
int kl =0;
for (int i=0;i<52822;i++){
for (int j=i;j<52822;j++){
        double temp[4];
            double P_1[3]={data[i][0],data[i][1],data[i][2]};
            double P_2[3]={data[i][3],data[i][4],data[i][5]};
            double P_3[3]={data[j][0],data[j][1],data[j][2]};
            double P_4[3]={data[j][3],data[j][4],data[j][5]};
         //   """   if project_xy(P_1,P_2,P_3,P_4) && project_xz(P_1,P_2,P_3,P_4) && ////////project_yz(P_1,P_2,P_3,P_4):"""
            points(P_1,P_2,P_3,P_4,temp);
            //printf("%lfn",temp[3]);
          
           if (temp[3]<=0.0001){
               
                //printf("nice||||");
                    //#xs.append((temp[0][0]+temp[1][0])/2)
                    //#ys.append((temp[0][1]+temp[1][1])/2)
                    //#zs.append((temp[0][2]+temp[1][2])/2)
                    intersect[kl][0]=temp[0];
                    intersect[kl][1]=temp[1];
                    intersect[kl][2]=temp[2];
                    kl++;
                    //fprintf(gilo,"%lf,%lf,%lf",intersect[kl][0],intersect[kl][1],intersect[kl][2],intersect[kl][3]);
    }}}
fprintf(gilo,"//////////////////////////////START\\\\\\\\\\\\\\\\\\\");
//printf("%f,%f,%f",data[4][0],data[4][1],data[4][5]);

printf("here3");;

int h=0;
int h1=0;

for(int i=1;i<51;i++){
    for (int j=1;j<51;j++){
        for (int k=1;k<401;k++){
            
            blocks[h][0]=(x_blocking[i-1]+x_blocking[i])/2;
            blocks[h][1]=(y_blocking[j-1]+y_blocking[j])/2;
            blocks[h][2]=(z_blocking[k-1]+z_blocking[k])/2;
            h++;
            if (k==1 || k==400){
                if (k==400) {k+=1;}
                for(int index=-1;index<=0;index++){
                detectors[h1][index+1][0]=x_blocking[i+index];
                detectors[h1][index+1][1]=y_blocking[j+index];
                detectors[h1][index+1][2]=z_blocking[k-1];
                }
                h1++;
                }
           }
}}
/*
for (int i=0;i<h;i++){

printf("i,%lfn",calc_prob(blocks[i],detectors[1],detectors[49]));
}
*/
/*fprintf(gilo,"||||||detectors||||||");
fprintf(gilo,"n");
for (int i=0;i<h1;i++){
fprintf(gilo,"[%lf,%lf,%lf],[%lf,%lf,%lf]",detectors[i][0][0],detectors[i][0][1],detectors[i][0][2],detectors[i][1][0],detectors[i][1][1],detectors[i][1][2]);
fprintf(gilo,"n");
}
fprintf(gilo,"||||||blocks||||||");
fprintf(gilo,"n");
for (int i=0;i<h;i++){
fprintf(gilo,"[%lf,%lf,%lf]",blocks[i][0],blocks[i][1],blocks[i][2]);
fprintf(gilo,"n");
}
*/
///////////////////////////////////////////////////
printf("here4");
int nstar(int i, int j){
int tempu=0;
    for(int k=0;k<h;k++){
                        if( (data[k][0]>detectors[i][0][0]) &&( data[k][0]<detectors[i][1][0]) && (data[k][3]>detectors[i][0][1]) && (data[k][3]<detectors[i][1][1])){
                            if( (data[k][1]>detectors[j][0][0]) &&(data[k][1]<detectors[j][1][0]) && (data[k][4]>detectors[j][0][1] )&&( data[k][4]<detectors[j][1][1])){
                                tempu++;
}}
}

return tempu;
}
//double lambda_0ld[h];
double * lambda_0ld=(double *)malloc(h*sizeof(double ));


//double Lambda[h];
double * Lambda=(double *)malloc(h*sizeof(double ));

void classify(){
    for (int j=0;j<h;j++){
    lambda_0ld[j]=0;
    for (int i=0;i<kl;i++){
        
    if((intersect[i][0]<(blocks[j][0]+0.5))&&(intersect[i][0]>(blocks[j][0]-0.5))&&(intersect[i][1]<(blocks[j][1]+0.5))&&(intersect[i][1]>(blocks[j][1]-0.5))&&(intersect[i][2]<(blocks[j][2]+0.5))&&(intersect[i][2]>(blocks[j][2]-0.5))){
        lambda_0ld[j]++;};
    }}
}
classify();
printf("here5"); 
double lambda(int b){double sum=0;
    classify();
    double var;
    for(int i=0;i<h1;i++){
        for(int j=i;j<h1;j++){
            if (detectors[i][2]==detectors[j][2]) continue;
            double temp=0.0;
    
            for (int p=0;p<h;p++){var=0.0;
                //if (detectors[i][2]!=detectors[j][2]){
                    var=(lambda_0ld[p]*calc_prob(blocks[p],detectors[i],detectors[j]));
                    //printf("n%lf,%lf,%lfn",calc_prob(blocks[p],detectors[i],detectors[j]),lambda_0ld[p],var);
    temp=temp+var;  
    }// if (temp!=temp) printf("popdsfoopds");
    sum=sum+((double)(nstar(i,j)*calc_prob(blocks[b],detectors[i],detectors[j]))/(double)temp);}}
    return (double)lambda_0ld[b]*(sum);
}

//printf("%lf",lambda_0ld[j];)
fprintf(gilo,"//////////////////////////////MLEM\\\\\\\\\\\\\\\\\\\n");
//printf("%lf,%lfn",lambda(20),lambda_0ld[20]);
printf("here6");
for (int u=0;u<h;u++){

    Lambda[u]=lambda(u);
    fprintf(gilo,"%d-%lf,%lf,%lf-%lfn",u,blocks[u][0],blocks[u][1],blocks[u][2],Lambda[u]);

//printf("%lfn",lambda_0ld[u]);
}       
//printf();

fclose(filo);


/*
for (int i=0;i<6;i++){
for (int j=0;j<6;j++){
for(int k=0;k<41;k++){
fprintf(gilo,"%lf,",blocking[i][j][k]);
}fprintf(gilo,"n");}
fprintf(gilo,"new");}*/
fclose(gilo);}

All the help I can find on the internet and stack overflow says dynamic allocation should fix the segmentation fault. But it does not in my case. I am returning back to C after a lot of years. Please point out if I am making any logical error while handling pointers.

PS Please ignore all of the comments. Most of them were either some unimplemented stuff or used for debugging.