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

#define _STR(X) #X
#define STR(X) _STR(X)

#define MAX_FILE_PATH_LENGTH FILENAME_MAX

#define MAX_MATRIX_SIZE 500
#define MAX_MATRIX_IN_FILE_SIZE (MAX_MATRIX_SIZE + 2)

#define MAX_ITERATION_STEPS 100

typedef struct {
	int n;
	double** data;
} Matrix;
Matrix* createMatrix(void);
void freeMatrix(Matrix* matrix);
void createMatrixRows(Matrix* matrix, int rows);
void printMatrix(Matrix* matrix);

typedef struct {
	int n;
	double* data;
} Vector;
Vector* createVector(void);
void initVector(Vector* vector, int size);
void freeVector(Vector* vector);
void printVector(Vector* vector);

typedef enum {
	JACOBI = 0, GAUSS_SEIDEL = 1
} Method;

void flushStdin(void);

bool load(const char* filename, Matrix* A, Vector* b, Vector* x);
Vector* solve(Method method, Matrix* A, Vector* b, Vector* x, double e);

int readMatrixLine(FILE* file, double* matrixLine, int maxCols);

int main(int argc, char* argv[]) {
	char filePath[MAX_FILE_PATH_LENGTH];

	if(argc >= 2) {
		strncpy(filePath, argv[1], MAX_FILE_PATH_LENGTH);
	} else {
		puts("Please enter the path of the file you'd like to open");
		int result = scanf("%s" STR(MAX_FILE_PATH_LENGTH) "[^\n]", filePath);
		if(result == EOF || result == 0 || filePath[0] == 0) {
			return 0;
		}
		if(result != 1) {
			fputs("Couldn't read file path - stopping\n", stderr);
			return 1;
		}
		flushStdin();
	}

	Matrix* matrix = createMatrix();
	Vector* b = createVector();
	Vector* x = createVector();

	int returnCode = 0;

	while(true) {
		if(load(filePath, matrix, b, x)) {

			// Debug outputs
			//puts("Data successfully loaded\nMatrix A:");
			//printMatrix(matrix);
			//puts("Vector b:");
			//printVector(b);
			//puts("Vector x:");
			//printVector(x);

			puts("Please enter the algorithm to use:\n\t0: Jacobi (default)\n\t1: Gauss-Seidel");
			int algorithm;

			while(true) {
				algorithm = getchar();
				if(algorithm == EOF) {
					goto end;
				}
				if(algorithm == '0' || algorithm == '1') {
					algorithm -= '0';
					break;
				}
				if(algorithm == '\n' || algorithm == 0) {
					puts("Defaulting to Jacobi");
					algorithm = 0;
					break;
				}
				fputs("Please enter 0, 1 or leave empty to exit!\n", stderr);
			}

			puts("Please enter the precision to use:");
			double e;
			int scan;
			while(true) {
				scan = scanf("%lf", &e);
				if(scan == EOF) {
					goto end;
				}
				if(scan == 1) {
					break;
				}
				flushStdin();
				fputs("Invalid input - please enter a valid floating point number!\n", stderr);
			}

			Vector* result = solve(algorithm, matrix, b, x, e);
			free(result);
			break;
		} else {
			fputs("Failed to load data from file.\nEnter new file path or leave empty to exit.\n", stderr);

			int result = scanf("%" STR(MAX_FILE_PATH_LENGTH) "[^\n]", filePath);
			if(result == EOF || result == 0 || filePath[0] == 0) {
				goto end;
			}
			if(result != 1) {
				fputs("Couldn't read file path - stopping\n", stderr);
				returnCode = 1;
				goto end;
			}
			flushStdin();
		}
	}

end:
	freeMatrix(matrix);
	freeVector(b);
	freeVector(x);

	printf("\nReturning with code %i\n", returnCode);
	return returnCode;
}

bool load(const char* filename, Matrix* matrix, Vector* b, Vector* x) {
	FILE* file = fopen(filename, "r");
	if(file == NULL) {
		fprintf(stderr, "Failed to open file \"%s\"\n", filename);
		return false;
	} else {
		double* firstLineBuffer = malloc(sizeof(double) * (MAX_MATRIX_IN_FILE_SIZE));

		int colsInFile = readMatrixLine(file, firstLineBuffer, MAX_MATRIX_IN_FILE_SIZE);

		if(colsInFile == -1) {
			fputs("Unexpected input on line 0", stderr);

			free(firstLineBuffer);
			goto failure;
		} else if(colsInFile == -2) {
			fprintf(stderr, "Exceeded maximum matrix size of %i\n", MAX_MATRIX_SIZE);

			free(firstLineBuffer);
			goto failure;
		} else {
			// success
			int cols = colsInFile - 1;
			initVector(b, cols);
			initVector(x, cols);

			b->data[0] = firstLineBuffer[cols];
			createMatrixRows(matrix, cols);
			matrix->data[0] = firstLineBuffer;

			int colsInLine;
			for(int i = 1; i < cols; i++) {
				matrix->data[i] = malloc(sizeof(double) * (colsInFile));
				colsInLine = readMatrixLine(file, matrix->data[i], colsInFile);
				if(colsInLine < 0) {
					if(i == cols - 1) {
						puts("Optional parameters are being used");

						// Matrix is one smaller than assumed
						cols--;
						b->n--;
						x->n--;
						matrix->n--;
						free(matrix->data[i]);

						x->data[0] = 0.4;

						// Copy b to x
						memcpy(x->data, b->data, b->n * sizeof(double));

						// extract b
						for(int j = 0; j < cols; j++) {
							b->data[j] = matrix->data[j][cols];
						}
						goto success;
					} else {
						fprintf(stderr, "Line %i contains illegal formatting - please fix!\n", i + 1);
						goto failure;
					}
				} else if(colsInLine != colsInFile) {
					fprintf(stderr, "Illegal line length found in line %i\n", i + 1);
					goto failure;
				} else {
					b->data[i] = matrix->data[i][cols];
				}
			}
			// successful but with no optional parameter

			// initialize vector x with zeros
			memset(x->data, 0, cols * sizeof(double));

		}
	}
success:
	fclose(file);
	return true;

failure:
	fclose(file);
	return false;
}


int readMatrixLine(FILE* file, double* matrixLine, int maxCols) {
	char nextChar = 0;
	int col = 0;
	double buffer;

	while(true) {
		if(fscanf(file, "%lf", &buffer) != 1) {
			return -1;
		}

		matrixLine[col] = buffer;

		col++;
		if(col > maxCols) {
			return -2;
		}
		nextChar = fgetc(file);
		if(nextChar == EOF || nextChar == '\n') break;
	}
	return col;
}

void flushStdin(void) {
	int c;
	while((c = getchar()) != '\n' && c != EOF && c != 0);
}

Vector* solve(Method method, Matrix* A, Vector* b, Vector* x, double e) {
	Vector* vectors = malloc(sizeof(Vector) * (MAX_ITERATION_STEPS + 1));



	initVector(&vectors[0], b->n);
	memcpy(vectors[0].data, x->data, b->n * sizeof(double));

	int vectorCount = 1;
    double temp;
	double delta = 0.0;

	if (method == 0){

	do {
        delta= 0.0;
        initVector(&vectors[vectorCount],b->n);
        for(int i=0; i<b->n; i++){
        temp=0.0;
            for(int j=0; j<b->n; j++){
                if(j!=i){
                    temp= temp+A->data[i][j]*vectors[vectorCount-1].data[j];
                }
            }
        vectors[vectorCount].data[i]=(b->data[i]-temp)/A->data[i][i];
        delta=fmax(fabs(vectors[vectorCount].data[i]-vectors[vectorCount-1].data[i]),delta);
        }
        printVector(&vectors[vectorCount]);
		vectorCount++;
	} while (delta > e && vectorCount < MAX_ITERATION_STEPS);

	}

	if (method ==1){


        do{
            delta = 0.0;
            if (vectorCount==1){
                initVector(&vectors[vectorCount], b->n);
                memcpy(vectors[vectorCount].data, x->data, b->n * sizeof(double));}
            else {
                    initVector(&vectors[vectorCount], b->n);
                    memcpy(vectors[vectorCount].data, vectors[vectorCount-1].data, b->n * sizeof(double));}

            for (int i = 0; i<b->n; i++){
            temp=0.0;
                for (int j=0; j<b->n; j++){
                    if (j!=i){
                        temp = temp + A->data[i][j]*vectors[vectorCount].data[j];

                    }
                }
            vectors[vectorCount].data[i] = (b->data[i]-temp)/A->data[i][i];
            delta=fmax(fabs(vectors[vectorCount].data[i]-vectors[vectorCount-1].data[i]),delta);

            }
            printVector(&vectors[vectorCount]);
            vectorCount++;
        }while(delta > e && vectorCount < MAX_ITERATION_STEPS);
	}

	// MAX_ITERATION_STEPS enthält die maximal zulässige Anzahl an Iterationsschritten (100)
	// Die einzelnen Vektoren müssen noch mit initVector initialisiert werden


	// HIER kommt der Code hin ;)

	// on success
	// Sei x die Anzahl der durchgeführten Iterationschritte. Dann setzt vectors[x+1].n = 0. Damit weiß das folgende Programm wie viele Schritte getätigt wurden.
	vectors[vectorCount].n=0;

	return vectors;
}

inline Matrix* createMatrix(void) {
	Matrix* matrix = malloc(sizeof(Matrix));
	matrix->n = 0;
	return matrix;
}

inline void createMatrixRows(Matrix* matrix, int rows) {
	matrix->n = rows;
	matrix->data = malloc(sizeof(double*) * rows);
}

inline void freeMatrix(Matrix* matrix) {
	for(int i = 0; i < matrix->n; i++) {
		free(matrix->data[i]);
	}
	free(matrix->data);
	free(matrix);
}

void printMatrix(Matrix* matrix) {
	for(int i = 0; i < matrix->n; i++) {
		for(int j = 0; j < matrix->n; j++) {
			printf("%lf, ", matrix->data[i][j]);
		}
		puts("");
	}
}

inline Vector* createVector(void) {
	Vector* vector = malloc(sizeof(Vector));
	vector->n = 0;
	return vector;
}

inline void initVector(Vector* vector, int size) {
	vector->n = size;
	vector->data = malloc(sizeof(double) * size);
}

inline void freeVector(Vector* vector) {
	free(vector->data);
	free(vector);
}

void printVector(Vector* vector) {
	for(int i = 0; i < vector->n; i++) {
		printf("%lf, ",  vector->data[i]);
	}
	puts("");
}