/*The code to turn off all bulbs in a matrix of bulbs */ #include #include #include #define N 499 typedef double ftype; typedef ftype ** MATRIX; typedef ftype * VECTOR; typedef const ftype ** CONST_MATRIX; typedef const ftype * CONST_VECTOR; #define INIT_a 0.0 #define INIT_b 100.0 #define INIT_c 75.0 #define INIT_d 50.0 #define ERROR 0.00001 MATRIX matrix_alloc(int n){ MATRIX x = (MATRIX)malloc(sizeof(ftype *)*n+sizeof(ftype)*n*n); int i; x[0]=(ftype *)(x+n); for(i=1;i0)sum-=in[i-1]; if(i0)sum-=in[i-1][j]; if(i0)sum+=in[i-1][j]; if(i0)sum+=in[i-1]; if(i=1)h[i][i-1]=1.0; h[i][i]=-1.0; if(imax_value){ max_index=j; max_value=v; } } if(max_value>ERROR){ freedom_index[fc]=i; j=max_index; if(j!=fc){ //Swap line j and fc ftype temp; for(k=i;kERROR){ failed=1; fprintf(stderr,"Error value %f\n",ME[i]); } } fprintf(stderr,"Freedom factor = %d\n", n-fc); if(failed){ printf("NO SOLUTION\n"); free(freedom_index); return; } k=n-1; for(i=fc-1;i>=0;i--){ int index=freedom_index[i]; while(k>index){ ME[k--]=0.0; } ME[k]=ME[i]; k--; } while(k>=0){ ME[k--]=0.0; } //Now ME hold's one solution for X[1], //Output one solution { MATRIX x; x=matrix_alloc(n); vector_copy(x[0],ME,n); matrix_mul_vector(x[1],H,x[0],n); vector_rdiff(x[1],init[0],n); for(k=2;k