NO

Author Topic: Matrix limitation ...  (Read 7143 times)

delper

  • Guest
Matrix limitation ...
« on: July 22, 2008, 06:07:46 AM »
Hi C Wizards ...

     Although I posted a similar question before, this is slightly different.  It looks like the C compiler is self-limiting itself to small matrices...  I need to perform matrix operations that should be 10000x10000 (I perform quite a bit of matrices operation, but multiplication is a good example).  However, depending in how I run the matrix multiplication code I get different maximum matrix sizes ....

     For example if I run this code:

Code: [Select]
#include <stdio.h>
#include <string.h>
#include <math.h>
#include <stdlib.h>
#include <time.h>
#include <stdlib.h>
#define MaxMatrix 20000

// Protos ...
float *MatrixMultiply (int r1, int c1, float Array1[][c1], int c2, float Array2[][c2], float Result[r1][c2]);   
float *RandomMatrix (int Size, float Matrix[][Size]);

int main (void)
{
int Num;
for (Num = 1; Num < MaxMatrix; Num ++)
{
printf ("\n Matrix %ix%i \n", Num, Num);
float Random_M[Num][Num];
float Random_M2[Num][Num];
RandomMatrix (Num, Random_M);
MatrixMultiply (Num, Num, Random_M, Num, Random_M, Random_M2);
}
return 0;
}

float *MatrixMultiply (int r1, int c1, float Array1[][c1], int c2, float Array2[][c2], float Result[r1][c2])   
{
int i, j, k;
for (i = 0; i < r1; i++)
{
for (j = 0; j < c2; j++)
{
Result[i][j] = 0.0;
for (k = 0; k < c1; k++)  
{
Result[i][j] += Array1[i][k] * Array2[k][j];
}
}
}
return Result[0];
}


float *RandomMatrix (int Size, float Matrix[][Size])
{
int i, j;
srand((unsigned)time(NULL) );
if (Size > 1)
{
for (i = 0; i < Size; i++)
{
for (j = 0; j < Size; j++)
{

Matrix [i][j] = ((float) (rand() / ((float)RAND_MAX + 1.0)))*3;
}
}
}
return Matrix[0];
}

With this I can multiply matrices up to 73x73 and then afterwards I get an error message saying that:

====================
Program Error:
Program.exe has generated errors and will be closed by Windows.  You will need to restart the program.

An error log is being generated.
====================

     However, if I set up the matrix size to a specific size then I can multiply matrices that are 359x359.   For example:

Code: [Select]
#include <stdio.h>
#include <string.h>
#include <math.h>
#include <stdlib.h>
#include <time.h>
#include <stdlib.h>
#define MaxMatrix 20000

// Protos ...
float *MatrixMultiply (int r1, int c1, float Array1[][c1], int c2, float Array2[][c2], float Result[r1][c2]);   
float *RandomMatrix (int Size, float Matrix[][Size]);

int main (void)
{
int Num = 359;
  float Random_M[Num][Num];
  float Random_M2[Num][Num];
RandomMatrix (Num, Random_M);
MatrixMultiply (Num, Num, Random_M, Num, Random_M, Random_M2);
return 0;
}

float *MatrixMultiply (int r1, int c1, float Array1[][c1], int c2, float Array2[][c2], float Result[r1][c2])   
{
int i, j, k;
for (i = 0; i < r1; i++)
{
for (j = 0; j < c2; j++)
{
Result[i][j] = 0.0;
for (k = 0; k < c1; k++)
{
Result[i][j] += Array1[i][k] * Array2[k][j];
}
}
}
return Result[0];
}


float *RandomMatrix (int Size, float Matrix[][Size])
{
int i, j;
srand((unsigned)time(NULL) );
if (Size > 1)
{
for (i = 0; i < Size; i++)
{
for (j = 0; j < Size; j++)
{
Matrix [i][j] = ((float) (rand() / ((float)RAND_MAX + 1.0)))*3;
}
}
}
return Matrix[0];
}

If I attempt to multiply matrices that are bigger than 359 then I get the same error than before, but as you can see this time the matrix was bigger (however not enough for my case).

     Does anyone of you know why I get this limitation in matrix handling?   

     Is there any way to increase the matrix handling capacity in my code?   

     Is this limitation due to hardware (RAM limitation)?

Any help in this matter is greatly appreciated.

    NOTE:  I compiled the application as a Win32 Console application (if this have to be compiled in a different way, please let me know).

Regards,

Delper

JohnF

  • Guest
Re: Matrix limitation ...
« Reply #1 on: July 22, 2008, 08:09:47 AM »
You are running out of stack space.

Have to calculated how much ram you will require for a 10000x10000 matrix?

John

delper

  • Guest
Re: Matrix limitation ...
« Reply #2 on: July 23, 2008, 05:27:07 AM »
Hi JohnF,

     No, I haven't calculated the ram requirement.   How do I do that?  (probably 10000x10000x4 =  400000000 bytes = 400 Megabytes).  In my deep ignorance, shouldn't virtual memory kick in?  My matrices will have variable sizes ranging from 100x100 to almost 10000x10000 (probably 99% of my cases will be in the lower 1000x1000).   

     However, if your theory is right, why I get two different results when running "almost" the same code?   It looks like the one running in a loop does not return the memory to the system.  Is there anything I should do in order to "release" the memory used?    Anyway, it does not matter if I am running the program alone or with other applications open, the compiled PellesC program end at the same point (in both cases).   Unless, you have gotten a different result, I'd like to know what kind of (super)computer you have!  ;-)

     On the other hand, I tested the same code with a different compiler (in the same computer), and somehow I get bigger matrices (not yet what I need, but matrices bigger than 359x359).   Why that happen?   I am wondering if there is any flag or compiler option that I should be using.

     Another stupid question: is there anyway to bypass the stack limitation? 

     Your help is greatly appreciated.

Regards,

JohnF

  • Guest
Re: Matrix limitation ...
« Reply #3 on: July 23, 2008, 07:42:40 AM »
Hi JohnF,

     No, I haven't calculated the ram requirement.   How do I do that?  (probably 10000x10000x4 =  400000000 bytes = 400 Megabytes).  In my deep ignorance, shouldn't virtual memory kick in?  My matrices will have variable sizes ranging from 100x100 to almost 10000x10000 (probably 99% of my cases will be in the lower 1000x1000).   

     However, if your theory is right, why I get two different results when running "almost" the same code?   It looks like the one running in a loop does not return the memory to the system.  Is there anything I should do in order to "release" the memory used?    Anyway, it does not matter if I am running the program alone or with other applications open, the compiled PellesC program end at the same point (in both cases).   Unless, you have gotten a different result, I'd like to know what kind of (super)computer you have!  ;-)

     On the other hand, I tested the same code with a different compiler (in the same computer), and somehow I get bigger matrices (not yet what I need, but matrices bigger than 359x359).   Why that happen?   I am wondering if there is any flag or compiler option that I should be using.

     Another stupid question: is there anyway to bypass the stack limitation? 

     Your help is greatly appreciated.

Regards,

Virtual memory will not 'kick' in when using the stack frame for the buffer.
     
Those two codes were not alike in that with the for loop, you were allocating new stack memory each loop.
     
Getting slightly bigger matrices with other compiles will not solve the problem as you will still over-run the limit. You can increase the stack size but I doubt it will allow 400 megs.
     
You could more easily use global arrays as these are not allocated on the stack.
     
My guess is that you will have to use allocated memory, using malloc for example. I would also suggest you ask about this on a forum that deals with these problems, say a 3D modeling forum or something like that.
   
Unless someone here has experience with these large arrays?

EDIT:

To make better use of the stack: VLA's are only de-allocated when the leaving the function.

Code: [Select]
void DoMat(int Num)
{
float Random_M[Num][Num];
  float Random_M2[Num][Num];
RandomMatrix (Num, Random_M);
MatrixMultiply (Num, Num, Random_M, Num, Random_M, Random_M2);
}

int main (void)
{
int Num;
for (Num = 1; Num < MaxMatrix; Num ++)
{
printf ("\n Matrix %ix%i \n", Num, Num);
DoMat(Num);
}
return 0;
}
   

John
« Last Edit: July 23, 2008, 08:27:47 AM by JohnF »

delper

  • Guest
Re: Matrix limitation ...
« Reply #4 on: July 24, 2008, 01:25:35 AM »
Hi JohnF

     Thank you a lot for your suggestion and all the useful information.  I've been naive when writing the matrix code expecting it to tackle any size, when the stack limitation is a critical issue.   

    Your suggestion on using malloc sounds great, but do you mind helping me construct a matrix using mallocs?  ... And a basic operation as addition, so I can understand how to transform my entire math library to "malloc" base matrices....

    Boy, I just hope someone else can enlighten me on how to do this with the less pain possible....  This is already looking uphill ... ;-(

Regards,

JohnF

  • Guest
Re: Matrix limitation ...
« Reply #5 on: July 24, 2008, 08:36:21 AM »
Pain is not always a bad thing. :)

I think you should attempt it yourself - that's the way to learn.

timovjl already supplied you with the function MatrixMultiply3x3p which should spare you some of the pain.

John
« Last Edit: July 24, 2008, 09:09:01 AM by JohnF »

JohnF

  • Guest
Re: Matrix limitation ...
« Reply #6 on: July 26, 2008, 08:17:14 PM »
I'm not brushing you off, come back if you get stuck.

John

delper

  • Guest
Re: Matrix limitation ...
« Reply #7 on: August 04, 2008, 05:30:16 AM »
Hi JohnF (and any other C-Experts),

     Taking your word, I hope you can help me out on this quite puzzling (and very simple) example:

Code: [Select]
#include <stdio.h>
#include <stdlib.h>

#define ISize 10

int main(void)
{
float *M2_ptr;
int i;
float p = 6.6f;

M2_ptr = (float *)malloc((ISize)*sizeof(float));

// Loading some data into the reserved space.....
for (i = 0; i< ISize; i++, M2_ptr++)
{
M2_ptr = &p;
p += 0.01f;
}
for (i = 0; i< ISize; i++, M2_ptr--);

// Checking content in memory....
for (i = 0; i < ISize; i++, M2_ptr++)
{
printf ("%f, ", &M2_ptr); 
}
free (M2_ptr);
return 0;
}     

As you can see, I am trying to use "malloc" to reserve the memory space.  I was planning to use malloc as a way to overcome the limitation in the heap you described before (please let me know if there is any limitation here).

However, I am doing something wrong, since no data seems to be written in the memory allocated.   What should I do to make this work?

JohnF

  • Guest
Re: Matrix limitation ...
« Reply #8 on: August 04, 2008, 08:09:05 AM »
Here is the corrected code. I have marked the lines that were wrong with 'HERE'.

Code: [Select]
#include <stdio.h>
#include <stdlib.h>

#define ISize 10

int main(void)
{
float * M2_ptr;
int i;
float p = 6.6f;

M2_ptr = (float *)malloc((ISize)*sizeof(float));

// Loading some data into the reserved space.....
for (i = 0; i< ISize; i++, M2_ptr++)
{
*M2_ptr = p; // HERE
p += 0.11f;
}
for (i = 0; i< ISize; i++, M2_ptr--);

// Checking content in memory....
for (i = 0; i < ISize; i++, M2_ptr++)
{
printf ("%f, ", *M2_ptr);   // HERE
}
free (M2_ptr);
return 0;
}     

As you seem to be having trouble using pointers I suggest you use the code I've attached. It uses the following method to allocate memory.

Code: [Select]
float (*Random_M)[MATSIZE] = malloc(...);

which allows using the allocated memory like you would an array.

John

delper

  • Guest
Re: Matrix limitation ...
« Reply #9 on: August 05, 2008, 06:20:37 AM »
Hi JohnF,

      Thanks a lot for the correction.   Wow that was a pretty simple star!  ;-)

      Your suggestion also makes a lot of sense, since my data comes in rows and then I'd be able to handle matrix data in a more "physical sense".   

      The code now works beautifully up to a point.   I guess I am hitting another "resource limitation" of some short (maybe two).   Please take a look at the code:

Code: [Select]
#include <stdio.h>
#include <stdlib.h>

#define ISize 635
#define Start 0.0001f
#define Increment 0.0001f

int main(void)
{
float *M2_ptr;
int i;
float p = Start;

M2_ptr = (float *)malloc((ISize)*sizeof(float));

// Loading some data .....
for (i = 0; i< ISize; i++, M2_ptr++)
{
*M2_ptr = p;
p += Increment;
}
for (i = 0; i< ISize; i++, M2_ptr--);
for (i = 0; i < ISize; i++, M2_ptr++)
{
printf ("%f, ", *M2_ptr); 
}
free (M2_ptr);
return 0;
}

At least in my computer, I get the numbers starting from 0.000100 and correctly up to 0.063400 (as you can see I am incrementing the float number by 0.0001), however, the next numbers are somehow wrong: 0.063401.   How come this last number increased by 0.000001    (the following numbers keep this weird increment!).   

      If Start and Increment are 0.001 then the results are wrong in a different location (315).   Is this some short of floating point limitation?   

      Another weird and probably more devastating issue is that the array cannot be as large as I want it to be.   It seems that (at least in my computer) I can only let ISize to be as big as 6349, if I try to go higher (6350) then I get an error message "Program.exe has stopped working .... bla bla bla".    Now, what I have done wrong????    Oooh boy, this is not my week! 

      I hope you can chime in and enlighten this "missing in combat" boy to the freedom of computer resources limitation.... ;-)

JohnF

  • Guest
Re: Matrix limitation ...
« Reply #10 on: August 05, 2008, 07:52:52 AM »
Yes I think you are hitting a floating point limitation.

EDIT:
Using multiplication would produce less floating point errors.
Code: [Select]
// Loading some data .....
for (i = 0; i < ISIZE; i++, M2_ptr++)
{
sum = INCREMENT * i;
*M2_ptr = sum;
}
ENDEDIT:

What error message are you getting?

Take a look here, always check to see if malloc was successful. Make a copy of M2_ptr so that you can set it back to its original value.


Code: [Select]
#include <stdio.h>
#include <stdlib.h>

#define ISize 635
#define Start 0.0001f
#define Increment 0.0001f

int main(void)
{
float *M2_ptr, *origptr;
int i;
float p = Start;

M2_ptr = malloc((ISize) * sizeof(float));
if(!M2_ptr)
{
printf("Error Allocating memory!\n");
return 0;
}
origptr = M2_ptr; // Save pointer

// Loading some data .....
for (i = 0; i < ISize; i++, M2_ptr++)
{
*M2_ptr = p;
p += Increment;
}

M2_ptr = origptr; // Set pointer back to original
for (i = 0; i < ISize; i++, M2_ptr++)
{
printf("%f, ", *M2_ptr);
}

M2_ptr = origptr; // Set pointer back to original
free(M2_ptr);
return 0;
}

John
« Last Edit: August 05, 2008, 04:14:59 PM by JohnF »