Easy Program, Wont Work With Large N

[breenmachine]

This is just a small chunk of code from something larger but it will run with n<834, anything greater and I get a segmentation fault and it doesnt even print out my "Test" message I have no idea why. I'm guessing maybe I'm not allowed to use as much memory as I need on the system but I don't know, does it work for anyone else or does anyone know the problem? Thanks

#include <stdio.h>
#include <math.h>
#include <time.h>

main(){

        int n=4000,iter,maxit,rcheck,s;
        int i,j;
        double xnew[n],xold[n];
        double b[n],a[n*n],aii,sum;
        int r[n+1],c[(n*n)],l,m,k;
        double toler,diff,maxdif;

        clock_t t1;

        t1 = clock();

        r[0] = 2;

        printf("test");

        for (i=1;i<=n-1;i++){
                r[n+1] = 2;
                for (i=1;i<=n-1;i++){
                        c[i]=3;
                        a[i]=4;
                }
        }

        printf("Time: %f",clock()-t1);

}

[robatino]

  r[n+1] = 2;

That index is always out of bounds. BTW, in C90 variable-length arrays aren't allowed. Try using a macro for the size instead.

Edit: It also doesn't make sense to let the size of a[] and c[] be n*n when you're only using the first n elements.

Edit: And you have two nested loops, using the same variable i. That must be wrong.

[breenmachine]

Oops sorry, made a typo, it is supposed to be

r[i+1] = 2;

strange it did run with n+1 when n was under 834. How would you use a Macro for the array size? I'm pretty new to C

[robatino]

Some compilers such as gcc will let you use variable-length arrays in C90 as a nonstandard extension. If you compile with the flags "-W -Wall -ansi -pedantic", it will warn about this sort of thing. You could use

#define N 4000

at the top of your code, just under the #include's (it's conventional to use all caps for macro names). Then use N instead of n in your code. The preprocessor will substitute 4000 for N everywhere in your code before compiling.

[breenmachine]

Yea, I am using gcc, thats probably why it was working. The macro is kind of cool but it didn't change the segmentation fault i get with numbers above 834, do you see anything else wrong with it? Because I've been looking a few hours and it shouldnt be that hard, thanks.

[robatino]

Did you fix the problem with the two nested loops using the same variable? Also, the index in a[i] and c[i] probably needs changing, since currently you're using only the first N values even though you allocated N*N. There's also the possibility that the segmentation fault is due to code you haven't posted. Could you post the code you have so far? (First the part corresponding to what you posted before, and if that looks okay, post more later.)

[breenmachine]

Oh my god... sorry another typo, ill just post the whole code that isnt typo ridden, those 2 loops dont really use the same variable.

/****************************************************************
* Sequential program to solve Ax=b using Jacobi's method *
****************************************************************
* A             -       Compact Storage Vector
* C             -       Column vector for A
* R 		-       Row position vector for A
* x(n)          -       solution vector *
* b(n)          -       RHS vector *
****************************************************************/

#include <stdio.h>
#include <math.h>
#include <time.h>
#define n 4000

main()
{
	int iter,maxit,rcheck,s;
	int i,j;
	double xnew[n],xold[n];
	double b[n],a[n*n],aii,sum;
	int r[n+1],c[(n*n)],l,m,k;
	double toler,diff,maxdif;

	
	clock_t t1;
	toler = .000001;
	iter=0;
	maxit=50;
	maxdif=1.0;
	s=1;
	rcheck=1;

/* generate sample data for compact storage vector A and right- hand-side-vector (RHS) b */

	for (i = 0;i<=n-1;i++){
		rcheck = 1;
		for(j=0;j<=n-1;j++){
				if (i==j) { 
					a[s] = 4.0;
		  			c[s] = j;
					s=s+1;
                 			if (rcheck==1){ 
						r[i] = s-1;
						rcheck = 0;
					}

				}
                  
				if (i==j-1||i==j+1){
					a[s] = -1.0;
                                	c[s] = j; 
					s=s+1;                                
                               		if (rcheck==1){
                                        	r[i] = s-1;
                                        	rcheck = 0;
					}
                                

		    	  	}
		}
            b[i] = 2.0;

	}

/*FOR TESTING THE COMPACT STORAGE SCHEME*/
/*	for(i=1;i<=s;i++){
		printf("\n&#37;4.2f    %d    %d\n",a[i],c[i],r[i]);
        }
*/


	b[0] = 3.0;
	b[n-1] = 3.0;

/* The N+1 element in R must be set to the size of A+1 to give the sum loop
work to do on the last iteration*/
      	
	r[n] = s+1;

/* time section of code to solve Ax = b*/
	t1 = clock();


//set initial solution vector to RHS Vector b

      for (i = 0;i<=n-1;i++){
         xold[i] = b[i];
      }

// iterate until |xold - xnew| < toler
// or max number of iterations exceeded
// R(I) represents the first non zero in a row
// so it works to sum 1 row at a time, C(J) is the 
// corresponding column that element was in in the matrix

      while ((iter<maxit)&&(maxdif>=toler)){
		for (i=0;i<=n-1;i++){
			sum = b[i];		
			for (j = r[i];j<=r[i+1]-1;j++){
				if (c[j]!=i) sum = sum + (-a[j]*xold[c[j]]);
              			
			}



//Finds the A(I,I) element in the Compact Storage Vector
                	m=r[i+1]-1;
                	l=r[i];
                	k=l;
                	if (c[k]==i) aii=a[k];
                                
                	while ((c[k]!=i)&&(k<=m)){
                       		k=k+1;
                       		if (c[k]==i) aii=a[k];
                	}
                                


	      

		xnew[i] = sum / aii; 
	}

// determine maximum vector norm

              maxdif = fabs(xold[0]-xnew[0]);
              for (i=1;i<=n-1;i++){
              		diff = fabs(xold[i]-xnew[i]);
              		if (diff>maxdif) maxdif = diff;   
	      }

// reset recent approximation to the solution vector x to the newly generated x values

           for(i=0;i<=n-1;i++){
                xold[i]=xnew[i];
           }

// keep track of the number of iterations required to reach the given tolerance

           iter=iter+1;
	}

      printf("\nElapsed time = %f",clock()-t1);
      printf("\nSystem size = %d x %d",n,n);
      printf("\nNumber of iterations = %d",iter);
      printf("\nMaximum Error = %f\n",maxdif);

// used for testing solution
//           printf("\n");
//           for(i=0;i<=n-1;i++){
//                printf("%f\n",xold[i]);
//           }
 
}

[robatino]

That's a lot of code. Try putting a printf() statement in the code, and then move it around to narrow down where it's crashing. You might have to use fflush() right afterwards - I'm not sure whether the output stream is flushed after a segfault.

[breenmachine]

Even with a printf directly after the variable declarations nothing gets printed, it gives the Seg fault immediatly, no matter where I put a print I get nothing.

[robatino]

I tried compiling and running the code you posted originally (letting gcc use nonstandard extensions) and I get the same sort of behavior, although not with the exact same maximum for n. I think the problem is that you're trying to allocate too much space using static arrays, which explains why it crashes during the variable definitions. You commonly can't allocate more than a few MB, sometimes not more than one, this way. Try allocating the arrays dynamically, using malloc() instead. This also has the advantage of letting the array sizes be determined at runtime, so you don't need the macro anymore.

[breenmachine]

ok cool, ill google malloc() and see if i can figure it out thanks