C: Nested functions using structs break the program if order of call is changed-CodePudding

I have ecountered a problem douring a programming exercise and i can not find the source of the problem since the program is compiling with gcc -std=c11 -Wall -Wextra -Wpedantic test.c -o test.exe without errors nor flags. The specific problem arrises when the call of square_matrix_multiplication is put at the end of the triagonal_invert function (FULL CODE AT THE END):

Square_Matrix tridigonal_invert(Square_Matrix matrix){
    int dim = matrix.n;
    double factor;
    Square_Matrix inverse = {dim, calloc(dim*dim,sizeof(double))};
    for(int i = 0;i<dim;i  ){
        inverse.data[i dim*i]=i;
    }
    Square_Matrix lower =   {dim, calloc(dim*dim,sizeof(double))};
    Square_Matrix upper =   {dim, calloc(dim*dim,sizeof(double))};
   // here works:
    print_square_matrix(square_matrix_multiplication(inverse,lower));
   //###############################################################
    lower.data[0]=1;
    upper.data[0]=matrix.data[0];
    upper.data[1]=matrix.data[1];
    for(int i = 1 ; i< dim; i  ){
        lower.data[i i*dim]=1;
        lower.data[i-1 i*dim]=matrix.data[i-1 i*dim]/upper.data[(i-1) (i-1)*dim];
        upper.data[i i*dim]=matrix.data[i i*dim]-lower.data[i-1 i*dim]*matrix.data[i (i-1)*dim];
        upper.data[i 1 i*dim]=matrix.data[i 1 i*dim];
        for(int j =0; j<i;j  ){
            inverse.data[j i*dim]=1;
            for(int k=j;k<i;k  ){
                inverse.data[j i*dim]*=-1*lower.data[k (k 1)*dim];
                
            }
        }
       /* factor=-1*lower.data[i-1 i*dim];
        for(int j=0;j<i;j  ){
            lower.data[j i*dim]=factor*lower.data[j (i-1)*dim];
        }*/
     
    }
    for(int i = dim-1;i>=0;i--){
        upper.data[i dim*i]=1/upper.data[i i*dim];
        factor=-1*upper.data[i dim*i]*upper.data[i 1 i*dim];
        for(int j=i 1;j<dim;j  ){
            upper.data[j i*dim]=factor*upper.data[j (i 1)*dim];            
            }

    }
   // here does not work:
    print_square_matrix(square_matrix_multiplication(inverse,lower));
   //###############################################################
    print_square_matrix(upper);
    print_square_matrix(inverse);
    
}

EDIT To clarify the problem:

If the call print_square_matrix(square_matrix_multiplication(inverse,lower)); is anywere after the second for statement the programs stops at function call exiting without error.

Here follows the full code:

#include<stdio.h>
#include<stdlib.h>
#include<string.h>
#include<math.h>


typedef struct 
{
    int n;
    double *data;
} Square_Matrix;

typedef struct {
    int lenght;
    double *data;
} n_vector;

void print_square_matrix(Square_Matrix matrix){
    int dim= matrix.n;
    printf("%dx%d Matrix\n",dim,dim);
    for(int i = 0;i<dim*dim;i  ){
        printf("%lf ",matrix.data[i]);
        if((i 1)%dim==0){
            printf("\n");
        }
        
    }
}

void print_n_vector(n_vector vector){

    int dim = vector.lenght;

    for( int i=0; i<dim; i  ){
        printf("%lf \n",vector.data[i]);
    }
}

void print_matrix_vector(Square_Matrix matrix,n_vector vector){

     int dim= matrix.n;
     int count  = 0;
    printf("MATRIX|VECTOR\n");
    for(int i = 0;i<dim*dim;i  ){
        printf("%lf ",matrix.data[i]);
        if((i 1)%dim==0){
            printf(" %lf \n",vector.data[count]);
            count =1;
        }
        
    }

}

Square_Matrix square_matrix_multiplication(Square_Matrix matrix1, Square_Matrix matrix2){
    int dim= matrix1.n;
    Square_Matrix result={dim,calloc(dim*dim,sizeof(double))};
    for(int i=0; i<dim;i  ){
        for(int j = 0; j<dim; j  ){
            for(int k=0;k<dim;k  ){
                result.data[j i*dim] =matrix1.data[k i*dim]*matrix2.data[j k*dim];
               
            }
            
        }
    }
    return result;
}


n_vector tridiagonal_solve(Square_Matrix matrix, n_vector vector){
    int dim = matrix.n;
    Square_Matrix lower = {dim, calloc(dim*dim,sizeof(double))};
    Square_Matrix upper = {dim, calloc(dim*dim,sizeof(double))};
    n_vector result = vector;



    lower.data[0]=1;
    upper.data[0]=matrix.data[0];
    upper.data[1]=matrix.data[1];
    


    for(int i = 1 ; i< dim; i  ){
        lower.data[i i*dim]=1;
        lower.data[i-1 i*   dim]=matrix.data[i-1 i*dim]/upper.data[(i-1) (i-1)*dim];
        upper.data[i i*dim]=matrix.data[i i*dim]-lower.data[i-1 i*dim]*matrix.data[i (i-1)*dim];
        upper.data[i 1 i*dim]=matrix.data[i 1 i*dim];
        result.data[i]= result.data[i]-lower.data[i-1 i*dim]*result.data[i-1];
    }

    result.data[dim-1]/= upper.data[dim*dim-1];  
    for (int i = dim-2; i>=0; i--)
    {
        result.data[i]=(result.data[i]-upper.data[i 1 i*dim]*result.data[i 1])/upper.data[i i*dim]; 
    }
    
return result;
}

Square_Matrix tridigonal_invert(Square_Matrix matrix){
    int dim = matrix.n;
    double factor;
    Square_Matrix inverse = {dim, calloc(dim*dim,sizeof(double))};
    for(int i = 0;i<dim;i  ){
        inverse.data[i dim*i]=i;
    }
    Square_Matrix lower =   {dim, calloc(dim*dim,sizeof(double))};
    Square_Matrix upper =   {dim, calloc(dim*dim,sizeof(double))};
    lower.data[0]=1;
    upper.data[0]=matrix.data[0];
    upper.data[1]=matrix.data[1];
    for(int i = 1 ; i< dim; i  ){
        lower.data[i i*dim]=1;
        lower.data[i-1 i*dim]=matrix.data[i-1 i*dim]/upper.data[(i-1) (i-1)*dim];
        upper.data[i i*dim]=matrix.data[i i*dim]-lower.data[i-1 i*dim]*matrix.data[i (i-1)*dim];
        upper.data[i 1 i*dim]=matrix.data[i 1 i*dim];
        for(int j =0; j<i;j  ){
            inverse.data[j i*dim]=1;
            for(int k=j;k<i;k  ){
                inverse.data[j i*dim]*=-1*lower.data[k (k 1)*dim];
                
            }
        }
       /* factor=-1*lower.data[i-1 i*dim];
        for(int j=0;j<i;j  ){
            lower.data[j i*dim]=factor*lower.data[j (i-1)*dim];
        }*/
     
    }
    for(int i = dim-1;i>=0;i--){
        upper.data[i dim*i]=1/upper.data[i i*dim];
        factor=-1*upper.data[i dim*i]*upper.data[i 1 i*dim];
        for(int j=i 1;j<dim;j  ){
            upper.data[j i*dim]=factor*upper.data[j (i 1)*dim];            
            }

    }
     print_square_matrix(square_matrix_multiplication(inverse,lower));
    print_square_matrix(upper);
    print_square_matrix(inverse);
    
}

int main(){
    Square_Matrix newMatrix={4,calloc(16,sizeof(double))};
    Square_Matrix index = {4,calloc(16,sizeof(double))};
    for(int i=0;i<4;i  ){
        index.data[i 4*i]=1;
    }
    n_vector vector={4,calloc(4,sizeof(double))};
    memcpy(newMatrix.data, (double[]){  2,-1, 0, 0, 
                                       -1, 2,-1, 0,
                                        0,-1, 2,-1,
                                        0, 0, -1, 2},16 * sizeof(double));
    memcpy(vector.data, (double[]){-5, 1, 4, 1},4*sizeof(double));

  tridigonal_invert(newMatrix);
  //index=square_matrix_multiplication(newMatrix,index);
  //print_square_matrix(index);
  free(newMatrix.data);
  free(vector.data);
  free(index.data);

  //print_square_matrix(index);
  
  //print_square_matrix(square_matrix_multiplication(newMatrix,index));
  /*vector =  tridiagonal_solve(newMatrix,vector);
  print_n_vector(vector);*/

}

CodePudding user response：

This line in tridigonal_invert

upper.data[i 1 i*dim]=matrix.data[i 1 i*dim];

is out of bounds for both memory regions when i == 3, invoking undefined behavior. dim == 4, both upper and matrix have space for 16 doubles reserved (valid indices of [0, 15]), yet when i == 3 in this loop, i 1 i*dim == 16, out of bounds by 1.

There could be other out-of-bounds occurrences as well, recommend running in valgrind to try to identify them.