I can't allocate memory for bigger matrices then ~200x200 Lapack Visual Studio in C Here is my code [duplicate]

Question

This question already has answers here:

What and where are the stack and heap? (32 answers)

Closed 8 years ago.

I am bench-marking the execution time of solving matrices , and i cant get more then ~200x200, i should go probably 1500x1500 or close to that. I am running this on VS .

#include <stdlib.h>
#include <stdio.h>
#include "lapacke.h"
#include "lapacke_config.h"
#include <time.h>

/* Auxiliary routines prototypes */
extern void print_matrix(lapack_complex_double *a, int m, int n);
extern void generate_matrix(lapack_complex_double * matrix, int w, int h);

/* Parameters */
#define N 179

/* Main program */
int main() {
    clock_t t;
    /* Locals */
    lapack_int n = N,info;
    lapack_int ipiv[N];
    lapack_complex_double a[N*N];
    lapack_complex_double b[N*N];

    FILE* fp1 = fopen("complexNumsLog.txt", "w");

    int w = 1, h = 1, h2 = 1, i = 0, j = 0; 


    for(i = 1; i <= N; i++){
        for(j = 1; j <= N; j++){
            w = i;
            h = i;
            h2 = j;

            generate_matrix(a, w, h);
            generate_matrix(b, w, h2);

            // print_matrix(a, w, h);
            // print_matrix(b, w, h2);

            // getchar(); 


            t = clock();
            info = LAPACKE_zgesv(LAPACK_ROW_MAJOR, w, h2, a, h, ipiv, b, h2);
            t = clock() - t;

            fprintf(fp1, "Matrix A: %3dx%3d ", w, h);
            fprintf(fp1, "Matrix B: %3dx%3d ", w, h2);
            fprintf(fp1, "%3d milliseconds %2.3f seconds\n",t,((float)t)/CLOCKS_PER_SEC);

            /* Check for the exact singularity */
            if( info > 0 ) {
                printf( "The diagonal element of the triangular factor of A,\n" );
                printf( "U(%i,%i) is zero, so that A is singular;\n", info, info );
                printf( "the solution could not be computed.\n" );
                getchar();
                exit( 1 );
            }
        }
        printf("%d\n", i);
    }
    getchar();
    exit( 0 );
} 

void print_matrix(lapack_complex_double* a, int m, int n) {
    int i, j;
    for( i = 0; i < m; i++ ) {
        for( j = 0; j < n; j++ )
            printf( " (%6.2f,%6.2f)", a[i*m+j].real, a[i*m+j].imag );
        printf( "\n" );
    }
    printf( "********************************************\n" );
}



void generate_matrix(lapack_complex_double * matrix, int w, int h){
    int i,j;
    double r;
    for(i = 0; i < w; i++){
        for(j = 0; j < h; j++){
            r = (rand()%1000 - 500)/100.0;
            matrix[i*w+j].real = r;
            r = (rand()%1000 - 500)/100.0;
            matrix[i*w+j].imag = r;;
        }
    }
}

Answer 1

You're blowing your stack. Your program is typically only allocated a few MB of stack space, and you're trying to allocate all of your data on the stack. When you get above about 200x200, you very quickly run out of stack space.

To fix this, you need to either allocate the memory at global scope or on the heap, where the size is only limited by your virtual address space and/or total available physical memory. Since the size is known at compile time, it's easiest to allocate it at global scope:

/* Parameters */
#define N 179
lapack_int ipiv[N];
lapack_complex_double a[N*N];
lapack_complex_double b[N*N];    

/* Main program */
int main() {
    ...
}

Answer 2

That's quite normal. You have this block of declarations in your code:

lapack_int ipiv[N];
lapack_complex_double a[N*N];
lapack_complex_double b[N*N];

which, I guess, is the dominant memory allocation. I suppose that a Lapack complex double is represented by a couple of double and a Lapack int as int, but even if this is not the case, the final result won't change that much.

Let's take some standard values (in Bytes) for sizeof (double), 8, and sizeof (int), 4.

Ok, we are ready to calculate how much memory you need

Memory occupation = 2*(8*2)*N² + 4*N = 32N² + 4N Bytes

What does this mean?

|    N   | Memory occupation |
|--------|-------------------|
|    50  |       80kb        |
|   100  |      320kb        |
|   150  |      721kb        |
|   500  |        8MB        |
|  1000  |       32MB        |
|  1500  |       72MB        |

Since those are local variables, they are allocated on the stack and... wow! 72MB of stack! That's quite a lot, consider that standard stack size is 1MB on Windows with MSVC.

You could allocate them dynamically with malloc:

lapack_int *ipiv = malloc(N * sizeof (lapack_int));
lapack_complex_double *a = malloc(N * N * sizeof (lapack_complex_double);
lapack_complex_double *b = malloc(N * N * sizeof (lapack_complex_double);

but remember to release resources before leaving the function, or you will leak (a lot of) memory. Do:

free(ipiv);
free(a);
free(b);

just before return.