#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <errno.h>
#include "matrix.h"

/* Create a new matrix structure of the
 * given number of rows and columns.  Data
 * is allocated as an array of pointers
 * rather than one contiguous block of
 * memory.  This is because the for large
 * matrix sizes the multiplication to
 * compute address offset can result in
 * overflow of a long data type.  Also
 * instead of one multiplication and one
 * addition to reference an entry there
 * are two pointers to follow.  Return
 * the matrix_t object or NULL if an
 * error occurred.
 */

matrix_t * create_matrix(size_t rows, size_t cols)
{
	matrix_t *matrix;
	size_t i;

	if (rows <= 0 || cols <= 0)
	{
		fprintf(stderr, "Error in mcreate: matrix dimensions " \
				"must not be zero or negative.\n");
	}
	
	matrix = malloc(sizeof(matrix_t));
	if (errno == ENOMEM)
	{
		fprintf(stderr, "Error in mcreate: insufficient memory " \
				"to allocate %d bytes for a new matrix " \
				"structure.\n", sizeof(matrix_t));
		return NULL;
	}

	matrix->rows = rows;
	matrix->cols = cols;
	
	matrix->data = malloc(rows * sizeof(MATRIXTYPE *));
	if (errno == ENOMEM)
	{
		fprintf(stderr, "Error in mcreate: insufficient memory " \
				"to allocate %d bytes for matrix row pointers.\n",
				rows * sizeof(MATRIXTYPE *));
		return NULL;
	}

	for (i = 0; i < rows; i++)
	{
		matrix->data[i] = malloc(cols * sizeof(MATRIXTYPE));
		if (errno == ENOMEM)
		{
			fprintf(stderr, "Error in mcreate: insufficient memory " \
					"to allocate %d bytes for matrix row %d.\n",
					cols * sizeof(MATRIXTYPE), i);
			return NULL;
		}
	}

	return matrix;
}

/* Deallocate all memory used by the
 * matrix object.
 */

void destroy_matrix(matrix_t *matrix)
{
	size_t i;

	if (matrix == NULL)
	{
		fprintf(stderr, "Error in mdestroy: null matrix " \
				"parameter passed to the function.  Possibly " \
				"you forgot to first call mcreate.\n");
		return;
	}

	for (i = 0; i < matrix->rows; i++)
		free(matrix->data[i]);	
	free(matrix->data);
	free(matrix);
	matrix = NULL;
}

/* Set all data elements of a matrix
 * to zero.  Return 0 if an error
 * occurred or 1 otherwise.
 */

int zero_matrix(matrix_t *matrix)
{
	size_t i;

	if (matrix == NULL)
	{
		fprintf(stderr, "Error in mzero: null matrix " \
				"parameter passed to the function.  Possibly " \
				"you forgot to first call mcreate.\n");
		return 0;
	}

	for (i = 0; i < matrix->rows; i++)
		memset(matrix->data[i], 0, matrix->cols * sizeof(MATRIXTYPE));

	return 1;
}

/* For the given matrix overwrite its
 * data to make it an identity matrix.
 * This only works on square matrices.
 * Return 0 if an error occurred or 1
 * otherwise.
 */

int identity_matrix(matrix_t *matrix) 
{
	size_t i;

	if (matrix == NULL)
	{
		fprintf(stderr, "Error in midentity: null matrix " \
				"parameter passed to the function.  Possibly " \
				"you forgot to first call mcreate.\n");
		return 0;
	}

	if (matrix->rows != matrix->cols)
	{
		fprintf(stderr, "midentity: source must be a square matrix.\n");
		return 0;
	}

	// Fill the entire matrix with '0'	
	for (i = 0; i < matrix->rows; i++)
	{
		memset(matrix->data[i], 0, matrix->cols * sizeof(MATRIXTYPE));
	}

	// Make '1' on the diagonal
	for (i = 0; i < matrix->rows; i++)
	{
		matrix->data[i][i] = 1;
	}

	return 1;
}

/* For the given matrix overwrite its
 * data to make it an identity matrix
 * where the diagonal values equal the
 * scale value instead of the normal 1.
 * This only works on square matrices.
 * Return 0 if an error occurred or 1
 * otherwise.
 */

int scaled_identity_matrix(matrix_t *matrix, double scale) 
{
	size_t i;

	if (matrix == NULL)
	{
		fprintf(stderr, "Error in midentity: null matrix " \
				"parameter passed to the function.  Possibly " \
				"you forgot to first call mcreate.\n");
		return 0;
	}

	if (matrix->rows != matrix->cols)
	{
		fprintf(stderr, "midentity: source must be a square matrix.\n");
		return 0;
	}

	// Fill the entire matrix with '0'	
	for (i = 0; i < matrix->rows; i++)
	{
		memset(matrix->data[i], 0, matrix->cols * sizeof(MATRIXTYPE));
	}

	// Make scale value on the diagonal
	for (i = 0; i < matrix->rows; i++)
	{
		matrix->data[i][i] = scale;
	}

	return 1;
}

/* Negate the matrix data.  This is
 * equivalent to muls_matrix by -1 but
 * the operation is performed inline.
 * Return 0 if an error occurred or 1
 * otherwise.
 */

int neg_matrix(matrix_t *matrix)
{
	size_t i, j;

	if (matrix == NULL)
	{
		fprintf(stderr, "Error in neg_matrix: null matrix " \
				"parameter passed to the function.  Possibly " \
				"you forgot to first call mcreate.\n");
		return 0;
	}

	for (i = 0; i < matrix->rows; i++)
	{
		for (j = 0; j < matrix->cols; j++)
		{
			matrix->data[i][j] *= -1;
		}
	}

	return 1;
}

/* Add two matrices and return the
 * result in the destination matrix.
 * The destination object must exist.
 * Inline operation are supported
 * where the destination is the same
 * as either the first or second
 * source matrix.  Return 0 if an
 * error occurred or 1 otherwise.
 */

int add_matrix(matrix_t *src1, matrix_t *src2, matrix_t *dest)
{
	size_t i, j;

	if (src1 == NULL || src2 == NULL || dest == NULL)
	{
		fprintf(stderr, "Error in add_matrix: null matrix " \
				"parameter passed to the function.  Possibly " \
				"you forgot to first call mcreate.\n");
		return 0;
	}

	if (src1->rows != src2->rows || src1->rows != dest->rows)
	{
		fprintf(stderr, "add_matrix: source and destination " \
			"matrices must have the same number of rows.\n");
		return 0;
	}

	if (src1->cols != src2->cols || src1->cols != dest->cols)
	{
		fprintf(stderr, "add_matrix: source and destination " \
			"matrices must have the same number of columns.\n");
		return 0;
	}

	for (i = 0; i < src1->rows; i++)
	{
		for (j = 0; j < src2->cols; j++)
		{
			dest->data[i][j] = src1->data[i][j] + src2->data[i][j];
		}
	}

	return 1;
}

/* Subtract two matrices and return the
 * result in the destination matrix.
 * The destination object must exist.
 * Inline operation are supported
 * where the destination is the same
 * as either the first or second
 * source matrix.  Return 0 if an
 * error occurred or 1 otherwise.
 */

int sub_matrix(matrix_t *src1, matrix_t *src2, matrix_t *dest)
{
	size_t i, j;

	if (src1 == NULL || src2 == NULL || dest == NULL)
	{
		fprintf(stderr, "Error in sub_matrix: null matrix " \
				"parameter passed to the function.  Possibly " \
				"you forgot to first call mcreate.\n");
		return 0;
	}

	if (src1->rows != src2->rows || src1->rows != dest->rows)
	{
		fprintf(stderr, "sub_matrix: source and destination " \
			"matrices must have the same number of rows.\n");
		return 0;
	}

	if (src1->cols != src2->cols || src1->cols != dest->cols)
	{
		fprintf(stderr, "sub_matrix: source and destination " \
			"matrices must have the same number of columns.\n");
		return 0;
	}

	for (i = 0; i < src1->rows; i++)
	{
		for (j = 0; j < src2->cols; j++)
		{
			dest->data[i][j] = src1->data[i][j] - src2->data[i][j];
		}
	}

	return 1;
}

/* Multiply two matrices and return the
 * result in the destination matrix.
 * The destination object must exist.
 * Inline operation are supported
 * only in the case where matrix size
 * criteria is met.  Return 0 if an
 * error occurred or 1 otherwise.
 */

int mul_matrix(matrix_t *src1, matrix_t *src2, matrix_t *dest)
{
	size_t i, j, k;
	MATRIXTYPE value;

	if (src1 == NULL || src2 == NULL || dest == NULL)
	{
		fprintf(stderr, "Error in mul_matrix: null matrix " \
				"parameter passed to the function.  Possibly " \
				"you forgot to first call mcreate.\n");
		return 0;
	}

	if (src1->cols != src2->rows)
	{
		fprintf(stderr, "mul_matrix: the number of columns in the first " \
			"source matrix must equal the number of rows in the second " \
			"source matrix.\n");
		fprintf(stderr, "columns in the first = %d; rows in the " \
			"second = %d\n", src1->cols, src2->rows);
		return 0;
	}
	else if (src1->rows != dest->rows)
	{
		fprintf(stderr, "mul_matrix: the number of rows in the first " \
			"source matrix must equal the number of rows in the " \
			"destination matrix.\n");
		fprintf(stderr, "rows in the first = %d; rows in the " \
			"destination = %d\n", src1->rows, dest->rows);
		return 0;
	}

	for (i = 0; i < src1->rows; i++)
	{
		for (j = 0; j < src2->cols; j++)
		{
			value = 0;
			for (k = 0; k < src1->cols; k++)
			{
				value += src1->data[i][k] *
						 src2->data[k][j];
			}
			dest->data[i][j] = value;
		}
	}

	return 1;
}

/* Add a matrix and scalar value and
 * return the result in the destination
 * matrix. The destination object must
 * exist. Inline operation are supported
 * where the destination is the same
 * as the source matrix.  Return 0 if
 * an error occurred or 1 otherwise.
 */

int adds_matrix(matrix_t *src, MATRIXTYPE scalar, matrix_t *dest)
{
	size_t i, j;

	if (src == NULL || dest == NULL)
	{
		fprintf(stderr, "Error in madds: null matrix " \
				"parameter passed to the function.  Possibly " \
				"you forgot to first call mcreate.\n");
		return 0;
	}

	if (src->rows != dest->rows)
	{
		fprintf(stderr, "madds: source and destination " \
			"matrices must have the same number of rows.\n");
		return 0;
	}

	if (src->cols != dest->cols)
	{
		fprintf(stderr, "madds: source and destination " \
			"matrices must have the same number of columns.\n");
		return 0;
	}

	for (i = 0; i < src->rows; i++)
	{
		for (j = 0; j < src->cols; j++)
		{
			dest->data[i][j] = src->data[i][j] + scalar;
		}
	}

	return 1;
}

/* Subtract a matrix and scalar value and
 * return the result in the destination
 * matrix. The destination object must
 * exist. Inline operation are supported
 * where the destination is the same
 * as the source matrix.  Return 0 if
 * an error occurred or 1 otherwise.
 */

int subs_matrix(matrix_t *src, MATRIXTYPE scalar, matrix_t *dest)
{
	size_t i, j;

	if (src == NULL || dest == NULL)
	{
		fprintf(stderr, "Error in subs_matrix: null matrix " \
				"parameter passed to the function.  Possibly " \
				"you forgot to first call mcreate.\n");
		return 0;
	}

	if (src->rows != dest->rows)
	{
		fprintf(stderr, "subs_matrix: source and destination " \
			"matrices must have the same number of rows.\n");
		return 0;
	}

	if (src->cols != dest->cols)
	{
		fprintf(stderr, "subs_matrix: source and destination " \
			"matrices must have the same number of columns.\n");
		return 0;
	}

	for (i = 0; i < src->rows; i++)
	{
		for (j = 0; j < src->cols; j++)
		{
			dest->data[i][j] = src->data[i][j] - scalar;
		}
	}

	return 1;
}

/* Multiply a matrix and scalar value and
 * return the result in the destination
 * matrix. The destination object must
 * exist. Inline operation are supported
 * where the destination is the same
 * as the source matrix.  Return 0 if
 * an error occurred or 1 otherwise.
 */

int muls_matrix(matrix_t *src, MATRIXTYPE scalar, matrix_t *dest)
{
	size_t i, j;

	if (src == NULL || dest == NULL)
	{
		fprintf(stderr, "Error in muls_matrix: null matrix " \
				"parameter passed to the function.  Possibly " \
				"you forgot to first call mcreate.\n");
		return 0;
	}

	if (src->rows != dest->rows)
	{
		fprintf(stderr, "muls_matrix: source and destination " \
			"matrices must have the same number of rows.\n");
		return 0;
	}

	if (src->cols != dest->cols)
	{
		fprintf(stderr, "muls_matrix: source and destination " \
			"matrices must have the same number of columns.\n");
		return 0;
	}

	for (i = 0; i < src->rows; i++)
	{
		for (j = 0; j < src->cols; j++)
		{
			dest->data[i][j] = src->data[i][j] / scalar;
		}
	}

	return 1;
}

/* Divide a matrix and scalar value and
 * return the result in the destination
 * matrix. The destination object must
 * exist. Inline operation are supported
 * where the destination is the same
 * as the source matrix.  Return 0 if
 * an error occurred or 1 otherwise.
 */

int divs_matrix(matrix_t *src, MATRIXTYPE scalar, matrix_t *dest)
{
	size_t i, j;

	if (src == NULL || dest == NULL)
	{
		fprintf(stderr, "Error in divs_matrix: null matrix " \
				"parameter passed to the function.  Possibly " \
				"you forgot to first call mcreate.\n");
		return 0;
	}

	if (src->rows != dest->rows)
	{
		fprintf(stderr, "divs_matrix: source and destination " \
			"matrices must have the same number of rows.\n");
		return 0;
	}

	if (src->cols != dest->cols)
	{
		fprintf(stderr, "divs_matrix: source and destination " \
			"matrices must have the same number of columns.\n");
		return 0;
	}

	for (i = 0; i < src->rows; i++)
	{
		for (j = 0; j < src->cols; j++)
		{
			dest->data[i][j] = src->data[i][j] / scalar;
		}
	}

	return 1;
}

/* Transpose the source matrix and return
 * a new matrix object.  Memory is
 * allocated and it is the responsibility
 * of the programmer to deallocate any
 * transposed matrix object.  Return NULL
 * if an error occurred or a pointer to
 * the transposed matrix otherwise.
 */

matrix_t * trans_matrix(matrix_t *source)
{
	size_t r, c;
	matrix_t *dest;

	if ((dest = create_matrix(source->cols, source->rows)) == NULL)
		return NULL;

	for (r = 0; r < source->rows; r++)
	{
		for (c = 0; c < source->cols; c++)
		{
			dest->data[c][r] = source->data[r][c];
		}
	}
	
	return dest;
}