diff --git a/include/matrix.h b/include/matrix.h index e13a4d29749cc3d9fc55465ef67575c0747c4a63..62d2977b8affaf55dbbd0c346734fb03d283c640 100644 --- a/include/matrix.h +++ b/include/matrix.h @@ -6,20 +6,20 @@ typedef zp* zp_mat; -void matrix_zp_from_int(zp_mat x, int *int_mat, int row, int col); +zp_mat matrix_zp_from_int(int *int_mat, int row, int col); -void matrix_zp_rand(zp_mat x, int row, int col); +zp_mat matrix_zp_rand(int row, int col); -void matrix_identity(zp_mat x, int size); +zp_mat matrix_identity(int size); int matrix_is_identity(zp_mat x, int size); -void matrix_transpose(zp_mat xt, zp_mat x, int row, int col); +zp_mat matrix_transpose(zp_mat x, int row, int col); -void matrix_merge(zp_mat xy, zp_mat x, zp_mat y, int row, int col_x, int col_y); +zp_mat matrix_merge(zp_mat x, zp_mat y, int row, int col_x, int col_y); -void matrix_multiply(zp_mat xy, zp_mat x, zp_mat y, int row_x, int row_y, int col_y); +zp_mat matrix_multiply(zp_mat x, zp_mat y, int row_x, int row_y, int col_y); -void matrix_inverse(zp_mat xi, zp_mat x, int size); +zp_mat matrix_inverse(zp_mat x, int size); #endif //PPANN_MATRIX_H diff --git a/src/matrix.c b/src/matrix.c index ecc3e54ce0515323b67170e7989dc87351d1071c..6396f21619f6dbad938d6fb8431f5ac90435fe9e 100644 --- a/src/matrix.c +++ b/src/matrix.c @@ -1,28 +1,37 @@ #include "matrix.h" -void matrix_zp_from_int(zp_mat x, int *int_mat, int row, int col) { +zp_mat matrix_zp_from_int(int *int_mat, int row, int col) { + zp_mat x; + x = (zp_mat) malloc(sizeof(zp) * row * col); for (int i = 0; i < row; i++) { for (int j = 0; j < col; j++) { zp_from_int(x[i * col + j], int_mat[i * col + j]); } } + return x; } -void matrix_zp_rand(zp_mat x, int row, int col) { +zp_mat matrix_zp_rand(int row, int col) { + zp_mat x; + x = (zp_mat) malloc(sizeof(zp) * row * col); for (int i = 0; i < row; i++) { for (int j = 0; j < col; j++) { rand_zp(x[i * col + j]); } } + return x; } -void matrix_identity(zp_mat x, int size) { +zp_mat matrix_identity(int size) { + zp_mat x; + x = (zp_mat) malloc(sizeof(zp) * size * size); for (int i = 0; i < size; i++) { for (int j = 0; j < size; j++) { if (i == j) zp_from_int(x[i * size + j], 1); else zp_zero(x[i * size + j]); } } + return x; } int matrix_is_identity(zp_mat x, int size) { @@ -35,15 +44,20 @@ int matrix_is_identity(zp_mat x, int size) { return 1; } -void matrix_transpose(zp_mat xt, zp_mat x, int row, int col) { +zp_mat matrix_transpose(zp_mat x, int row, int col) { + zp_mat xt; + xt = (zp_mat) malloc(sizeof(zp) * row * col); for (int i = 0; i < row; i++) { for (int j = 0; j < col; j++) { zp_copy(xt[j * row + i], x[i * col + j]); } } + return xt; } -void matrix_merge(zp_mat xy, zp_mat x, zp_mat y, int row, int col_x, int col_y) { +zp_mat matrix_merge(zp_mat x, zp_mat y, int row, int col_x, int col_y) { + zp_mat xy; + xy = (zp_mat) malloc(sizeof(zp) * row * (col_x + col_y)); for (int i = 0; i < row; i++) { for (int j = 0; j < col_x; j++) { zp_copy(xy[i * (col_x + col_y) + j], x[i * col_x + j]); @@ -52,9 +66,12 @@ void matrix_merge(zp_mat xy, zp_mat x, zp_mat y, int row, int col_x, int col_y) zp_copy(xy[i * (col_x + col_y) + col_x + j], y[i * col_y + j]); } } + return xy; } -void matrix_multiply(zp_mat xy, zp_mat x, zp_mat y, int row_x, int row_y, int col_y) { +zp_mat matrix_multiply(zp_mat x, zp_mat y, int row_x, int row_y, int col_y) { + zp_mat xy; + xy = (zp_mat) malloc(sizeof(zp) * row_x * col_y); zp temp; for (int i = 0; i < row_x; i++) { for (int j = 0; j < col_y; j++) { @@ -65,15 +82,14 @@ void matrix_multiply(zp_mat xy, zp_mat x, zp_mat y, int row_x, int row_y, int co } } } + return xy; } -void matrix_inverse(zp_mat xi, zp_mat x, int size) { +zp_mat matrix_inverse(zp_mat x, int size) { // Declare the row echelon matrix and generate it. - zp *identity, *row_echelon; - identity = (zp *) malloc(size * size * sizeof(zp)); - row_echelon = (zp *) malloc(2 * size * size * sizeof(zp)); - matrix_identity(identity, size); - matrix_merge(row_echelon, x, identity, size, size, size); + zp_mat identity, row_echelon; + identity = matrix_identity(size); + row_echelon = matrix_merge(x, identity, size, size, size); // Declare temp value. zp temp_multiplier, temp_neg; @@ -81,20 +97,20 @@ void matrix_inverse(zp_mat xi, zp_mat x, int size) { // Bottom left half to all zeros. for (int i = 0; i < size; i++) { for (int j = i; j < size; j++) { - if (i == j && fp_cmp_dig(row_echelon[i * 2 * size + j], 1) != RLC_EQ) { + if (i == j && !zp_is_int(row_echelon[i * 2 * size + j], 1)) { zp_inverse(temp_multiplier, row_echelon[i * 2 * size + i]); for (int k = i; k < size * 2; k++) { zp_multiply(row_echelon[j * 2 * size + k], row_echelon[j * 2 * size + k], temp_multiplier); } } - if (i == j && fp_cmp_dig(row_echelon[i * 2 * size + j], 0) == RLC_EQ) break; + if (i == j && zp_is_int(row_echelon[i * 2 * size + j], 0)) break; if (i != j) { zp_copy(temp_multiplier, row_echelon[j * 2 * size + i]); for (int k = i; k < size * 2; k++) { zp_multiply(temp_neg, temp_multiplier, row_echelon[i * 2 * size + k]); - fp_neg(temp_neg, temp_neg); + zp_neg(temp_neg, temp_neg); zp_add(row_echelon[j * 2 * size + k], row_echelon[j * 2 * size + k], temp_neg); } } @@ -107,16 +123,22 @@ void matrix_inverse(zp_mat xi, zp_mat x, int size) { zp_copy(temp_multiplier, row_echelon[j * 2 * size + i]); for (int k = i; k < size * 2; k++) { zp_multiply(temp_neg, temp_multiplier, row_echelon[i * 2 * size + k]); - fp_neg(temp_neg, temp_neg); + zp_neg(temp_neg, temp_neg); zp_add(row_echelon[j * 2 * size + k], row_echelon[j * 2 * size + k], temp_neg); } } } + // Copy over the output. + + zp_mat xi; + xi = (zp_mat) malloc(sizeof(zp) * size * size); + for (int i = 0; i < size; i++) { for (int j = 0; j < size; j++) { zp_copy(xi[i * size + j], row_echelon[i * 2 * size + size + j]); } } + return xi; } diff --git a/tests/test_matrix.c b/tests/test_matrix.c index 2f170742140feccd8db41c22cc6f7ce320dbdff1..24857d3205dc3f5312b243c4077c0bc2aada10c4 100644 --- a/tests/test_matrix.c +++ b/tests/test_matrix.c @@ -1,52 +1,48 @@ #include "matrix.h" int test_zp_from_int() { - int row = 3, col = 3; - zp x[row * col]; - int int_vec[] = {1, 2, 3, 4, 5, 6, 7, 8, 9}; - matrix_zp_from_int(x, int_vec, 3, 3); - return fp_cmp_dig(x[2 * col + 2], 9); + zp_mat x; + int int_mat[] = {1, 2, 3, 4, 5, 6, 7, 8, 9}; + x = matrix_zp_from_int(int_mat, 3, 3); + return zp_is_int(x[8], 9); } int test_transpose() { int row = 3, col = 3; - zp x[row * col], xt[col * row]; - matrix_zp_rand(x, row, col); - matrix_transpose(xt, x, row, col); + zp_mat x, xt; + x = matrix_zp_rand(row, col); + xt = matrix_transpose(x, row, col); return fp_cmp(xt[col - 1], x[2 * row]); } int test_identity() { - int size = 1000; - zp *x; - x = (zp *) malloc(size * size * sizeof(zp)); - matrix_identity(x, size); + int size = 100; + zp_mat x; + x = matrix_identity(size); return matrix_is_identity(x, size); } int test_merge() { int size = 10; - zp xy[(size + size) * size], x[size * size], y[size * size]; - matrix_zp_rand(x, size, size); - matrix_identity(y, size); - matrix_merge(xy, x, y, size, size, size); + zp_mat xy, x, y; + x = matrix_zp_rand(size, size); + y = matrix_identity(size); + xy = matrix_merge(x, y, size, size, size); return fp_cmp(x[2 * size + 1], xy[4 * size + 1]); } int test_multiply_vector() { - int int_x[5] = {1, 2, 3, 4, 5}; - int int_y[15] = {10, 20, 30, + int mat_x[5] = {1, 2, 3, 4, 5}; + int mat_y[15] = {10, 20, 30, 10, 20, 30, 10, 20, 30, 10, 20, 30, 10, 20, 30}; - zp x[5], y[15]; - matrix_zp_from_int(x, int_x, 1, 5); - matrix_zp_from_int(y, int_y, 5, 3); - - zp xy[3]; - matrix_multiply(xy, x, y, 1, 5, 3); + zp_mat x, y, xy; + x = matrix_zp_from_int(mat_x, 1, 5); + y = matrix_zp_from_int(mat_y, 5, 3); + xy = matrix_multiply(x, y, 1, 5, 3); return fp_cmp_dig(xy[1], 300); } @@ -54,24 +50,20 @@ int test_multiply_vector() { int test_inverse() { int size = 100; // Allocate space. - zp *x, *xi, *r; - x = (zp *) malloc(size * size * sizeof(zp)); - xi = (zp *) malloc(size * size * sizeof(zp)); - r = (zp *) malloc(size * size * sizeof(zp)); - matrix_zp_rand(x, size, size); - matrix_inverse(xi, x, size); - matrix_multiply(r, xi, x, size, size, size); + zp_mat x, xi, r; + x = matrix_zp_rand(size, size); + xi = matrix_inverse(x, size); + r = matrix_multiply(xi, x, size, size, size); return matrix_is_identity(r, size); } - int main() { // Init core and setup. core_init(); pc_param_set_any(); // Perform tests. - if (test_zp_from_int() != RLC_EQ) return 1; + if (test_zp_from_int() != 1) return 1; if (test_transpose() != RLC_EQ) return 1; if (test_identity() != 1) return 1; if (test_merge() != RLC_EQ) return 1;