diff --git a/ext/coo/coo_def.c b/ext/coo/coo_def.c index 563d1b9..b8ea941 100644 --- a/ext/coo/coo_def.c +++ b/ext/coo/coo_def.c @@ -1,6 +1,6 @@ VALUE coo_init(int argc, VALUE* argv, VALUE self) { coo_matrix* mat; - Data_Get_Struct(self, coo_matrix, mat); + TypedData_Get_Struct(self, coo_matrix, &coo_data_type, mat); if(argc > 0){ mat->dtype = sp_float64; @@ -31,19 +31,34 @@ VALUE coo_init(int argc, VALUE* argv, VALUE self) { VALUE coo_alloc(VALUE klass) { coo_matrix* mat = ALLOC(coo_matrix); - return Data_Wrap_Struct(klass, NULL, coo_free, mat); + return TypedData_Wrap_Struct(klass, &coo_data_type, mat); } -void coo_free(coo_matrix* mat) { +void coo_free(void* ptr) { + coo_matrix *mat = (coo_matrix*)ptr; + if (mat->shape) xfree(mat->shape); + if (mat->elements) xfree(mat->elements); + if (mat->ia) xfree(mat->ia); + if (mat->ja) xfree(mat->ja); xfree(mat); } +size_t coo_memsize(const void* ptr) { + coo_matrix *mat = (coo_matrix*)ptr; + size_t size = sizeof(mat); + if (mat->shape) size += mat->ndims; + if (mat->elements) size += mat->count; + if (mat->ia) size += mat->count; + if (mat->ja) size += mat->count; + return size; +} + VALUE coo_get_elements(VALUE self) { coo_matrix* input; - Data_Get_Struct(self, coo_matrix, input); + TypedData_Get_Struct(self, coo_matrix, &coo_data_type, input); VALUE* array = ALLOC_N(VALUE, input->count); for (size_t index = 0; index < input->count; index++) { @@ -57,7 +72,7 @@ VALUE coo_get_elements(VALUE self) { VALUE coo_get_coords(VALUE self) { coo_matrix* input; - Data_Get_Struct(self, coo_matrix, input); + TypedData_Get_Struct(self, coo_matrix, &coo_data_type, input); VALUE* array_ia = ALLOC_N(VALUE, input->count); for (size_t index = 0; index < input->count; index++) { @@ -78,7 +93,7 @@ VALUE coo_get_coords(VALUE self) { VALUE coo_get_count(VALUE self) { coo_matrix* input; - Data_Get_Struct(self, coo_matrix, input); + TypedData_Get_Struct(self, coo_matrix, &coo_data_type, input); return SIZET2NUM(input->count); } @@ -87,7 +102,7 @@ VALUE coo_get_count(VALUE self) { VALUE coo_get_ndims(VALUE self) { coo_matrix* input; - Data_Get_Struct(self, coo_matrix, input); + TypedData_Get_Struct(self, coo_matrix, &coo_data_type, input); return SIZET2NUM(input->ndims); } @@ -96,7 +111,7 @@ VALUE coo_get_ndims(VALUE self) { VALUE coo_get_dtype(VALUE self) { coo_matrix* input; - Data_Get_Struct(self, coo_matrix, input); + TypedData_Get_Struct(self, coo_matrix, &coo_data_type, input); return ID2SYM(rb_intern(DTYPE_NAMES[input->dtype])); } @@ -105,7 +120,7 @@ VALUE coo_get_dtype(VALUE self) { VALUE coo_get_shape(VALUE self) { coo_matrix* input; - Data_Get_Struct(self, coo_matrix, input); + TypedData_Get_Struct(self, coo_matrix, &coo_data_type, input); VALUE* array = ALLOC_N(VALUE, input->ndims); for (size_t index = 0; index < input->ndims; index++) { @@ -115,5 +130,102 @@ VALUE coo_get_shape(VALUE self) { return rb_ary_new4(input->ndims, array); } +VALUE coo_to_csr(VALUE self) { + coo_matrix* left; + TypedData_Get_Struct(self, coo_matrix, &coo_data_type, left); + + csr_matrix* result = ALLOC(csr_matrix); + result->dtype = left->dtype; + result->count = left->count; + result->ndims = left->ndims; + result->shape = ALLOC_N(size_t, result->ndims); + + for (size_t index = 0; index < result->ndims; index++) { + result->shape[index] = left->shape[index]; + } + + result->elements = ALLOC_N(double, left->count); + result->ip = ALLOC_N(size_t, left->shape[0] + 1); + result->ja = ALLOC_N(size_t, left->count); + + for(size_t index = 0; index <= result->shape[0]; index++) result->ip[index] = 0; + for(size_t index = 0; index < result->count; index++) { + size_t row_index_val = left->ia[index]; + result->ip[row_index_val] += 1; + } + for(size_t index = 0, cumsum = 0; index < result->shape[0]; index++) { + size_t temp = result->ip[index]; + result->ip[index] = cumsum; + cumsum += temp; + } + result->ip[result->shape[0]] = result->count; + for(size_t index = 0; index < result->count; index++) { + size_t row = left->ia[index]; + size_t dest = result->ip[row]; + + result->ja[dest] = left->ja[index]; + result->elements[dest] = left->elements[index]; + + result->ip[row] += 1; + } + + for(size_t index = 0, last = 0; index <= result->shape[0]; index++) { + size_t temp = result->ip[index]; + result->ip[index] = last; + last = temp; + } + + return TypedData_Wrap_Struct(CSR, &csr_data_type, result); +} + +VALUE coo_to_csc(VALUE self) { + coo_matrix* left; + TypedData_Get_Struct(self, coo_matrix, &coo_data_type, left); + + csc_matrix* result = ALLOC(csc_matrix); + result->dtype = left->dtype; + result->count = left->count; + result->ndims = left->ndims; + result->shape = ALLOC_N(size_t, result->ndims); + + for (size_t index = 0; index < result->ndims; index++) { + result->shape[index] = left->shape[index]; + } + + result->elements = ALLOC_N(double, left->count); + result->ia = ALLOC_N(size_t, left->count); + result->jp = ALLOC_N(size_t, left->shape[1] + 1); + + for(size_t index = 0; index <= result->shape[1]; index++) result->jp[index] = 0; + for(size_t index = 0; index < result->count; index++) { + size_t col_index_val = left->ja[index]; + result->jp[col_index_val] += 1; + } + for(size_t index = 0, cumsum = 0; index < result->shape[1]; index++) { + size_t temp = result->jp[index]; + result->jp[index] = cumsum; + cumsum += temp; + } + result->jp[result->shape[1]] = result->count; + for(size_t index = 0; index < result->count; index++) { + size_t col = left->ja[index]; + size_t dest = result->jp[col]; + + result->ia[dest] = left->ia[index]; + result->elements[dest] = left->elements[index]; + + result->jp[col] += 1; + } + + for(size_t index = 0, last = 0; index <= result->shape[1]; index++) { + size_t temp = result->jp[index]; + result->jp[index] = last; + last = temp; + } + + return TypedData_Wrap_Struct(CSC, &csc_data_type, result); +} + -#include "elementwise.c" \ No newline at end of file +#include "elementwise.c" +#include "iteration.c" \ No newline at end of file diff --git a/ext/coo/elementwise.c b/ext/coo/elementwise.c index 55818f8..7e6559a 100644 --- a/ext/coo/elementwise.c +++ b/ext/coo/elementwise.c @@ -1,3 +1,7 @@ +#include "string.h" +#include "math.h" + + double coo_perform_oper(double val_a, double val_b, char oper) { switch(oper) { case '+': @@ -18,10 +22,11 @@ double coo_perform_oper(double val_a, double val_b, char oper) { VALUE coo_elementwise_binary(VALUE self, VALUE another, char oper) { coo_matrix* left; coo_matrix* right; - Data_Get_Struct(self, coo_matrix, left); - Data_Get_Struct(another, coo_matrix, right); + TypedData_Get_Struct(self, coo_matrix, &coo_data_type, left); + TypedData_Get_Struct(another, coo_matrix, &coo_data_type, right); coo_matrix* result = ALLOC(coo_matrix); + result->dtype = left->dtype; result->count = 0; result->ndims = left->ndims; result->shape = ALLOC_N(size_t, result->ndims); @@ -79,7 +84,12 @@ VALUE coo_elementwise_binary(VALUE self, VALUE another, char oper) { } while (left_index < left->count) { - result->elements[result_index] = left->elements[left_index]; + double result_val = coo_perform_oper(left->elements[left_index], 0.0, oper); + if(fabs(result_val) < 1e-6) { //near to zero + left_index++; + continue; //skip current result value + } + result->elements[result_index] = result_val; result->ia[result_index] = left->ia[left_index]; result->ja[result_index] = left->ja[left_index]; @@ -88,7 +98,12 @@ VALUE coo_elementwise_binary(VALUE self, VALUE another, char oper) { } while (right_index < right->count) { - result->elements[result_index] = right->elements[right_index]; + double result_val = coo_perform_oper(0.0, right->elements[right_index], oper); + if(fabs(result_val) < 1e-6) { //near to zero + right_index++; + continue; //skip current result value + } + result->elements[result_index] = result_val; result->ia[result_index] = right->ia[right_index]; result->ja[result_index] = right->ja[right_index]; @@ -96,7 +111,7 @@ VALUE coo_elementwise_binary(VALUE self, VALUE another, char oper) { result->count++; } - return Data_Wrap_Struct(COO, NULL, coo_free, result); + return TypedData_Wrap_Struct(COO, &coo_data_type, result); } VALUE coo_add(VALUE self, VALUE another) { @@ -109,4 +124,57 @@ VALUE coo_sub(VALUE self, VALUE another) { VALUE coo_mul(VALUE self, VALUE another) { return coo_elementwise_binary(self, another, '*'); -} \ No newline at end of file +} + +double coo_unary_oper(double val, const char* oper) { + if (strcmp(oper, "sin") == 0) + return sin(val); + else if (strcmp(oper, "cos") == 0) + return cos(val); + else if (strcmp(oper, "tan") == 0) + return tan(val); + else + return 0.00; +} + +/* + Takes the matrix and performs unary operator elementwise +*/ +VALUE coo_elementwise_unary(VALUE self, const char* oper) { + coo_matrix* left; + TypedData_Get_Struct(self, coo_matrix, &coo_data_type, left); + + coo_matrix* result = ALLOC(coo_matrix); + result->dtype = left->dtype; + result->count = left->count; + result->ndims = left->ndims; + result->shape = ALLOC_N(size_t, result->ndims); + + for (size_t index = 0; index < result->ndims; index++) { + result->shape[index] = left->shape[index]; + } + + result->elements = ALLOC_N(double, left->count); + result->ia = ALLOC_N(size_t, left->count); + result->ja = ALLOC_N(size_t, left->count); + + for(size_t index = 0; index < result->count; index++) { + result->elements[index] = coo_unary_oper(left->elements[index], oper); + result->ia[index] = left->ia[index]; + result->ja[index] = left->ja[index]; + } + + return TypedData_Wrap_Struct(COO, &coo_data_type, result); +} + +VALUE coo_sin(VALUE self) { + return coo_elementwise_unary(self, "sin"); +} + +VALUE coo_cos(VALUE self) { + return coo_elementwise_unary(self, "cos"); +} + +VALUE coo_tan(VALUE self) { + return coo_elementwise_unary(self, "tan"); +} diff --git a/ext/coo/iteration.c b/ext/coo/iteration.c new file mode 100644 index 0000000..fab96ee --- /dev/null +++ b/ext/coo/iteration.c @@ -0,0 +1,27 @@ +VALUE coo_each(VALUE self) { + coo_matrix* input; + + TypedData_Get_Struct(self, coo_matrix, &coo_data_type, input); + + for (size_t index = 0; index < input->count; index++) { + rb_yield(DBL2NUM(input->elements[index])); + } + + return self; +} + +VALUE coo_each_with_indices(VALUE self) { + coo_matrix* input; + + TypedData_Get_Struct(self, coo_matrix, &coo_data_type, input); + + VALUE* result_array = ALLOC_N(VALUE, input->ndims + 1); + for (size_t index = 0; index < input->count; index++) { + result_array[0] = DBL2NUM(input->elements[index]); + result_array[1] = INT2NUM(input->ia[index]); + result_array[2] = INT2NUM(input->ja[index]); + rb_yield(rb_ary_new4(input->ndims + 1, result_array)); + } + + return self; +} \ No newline at end of file diff --git a/ext/csc/csc_def.c b/ext/csc/csc_def.c new file mode 100644 index 0000000..4d41287 --- /dev/null +++ b/ext/csc/csc_def.c @@ -0,0 +1,189 @@ +VALUE csc_init(int argc, VALUE* argv, VALUE self) { + csc_matrix* mat; + TypedData_Get_Struct(self, csc_matrix, &csc_data_type, mat); + + if(argc > 0){ + mat->dtype = sp_float64; + mat->ndims = 2; + mat->count = (size_t)RARRAY_LEN(argv[1]); + mat->shape = ALLOC_N(size_t, mat->ndims); + for(size_t index = 0; index < mat->ndims; index++) { + mat->shape[index] = (size_t)FIX2LONG(RARRAY_AREF(argv[0], index)); + } + mat->elements = ALLOC_N(double, mat->count); + mat->ia = ALLOC_N(size_t, mat->count); + mat->jp = ALLOC_N(size_t, mat->shape[1] + 1); + for(size_t index = 0; index <= mat->shape[1]; index++) mat->jp[index] = 0; + for(size_t index = 0; index < mat->count; index++) { + size_t col_index_val = (size_t)NUM2SIZET(RARRAY_AREF(argv[3], index)); + mat->jp[col_index_val] += 1; + } + for(size_t index = 0, cumsum = 0; index < mat->shape[1]; index++) { + size_t temp = mat->jp[index]; + mat->jp[index] = cumsum; + cumsum += temp; + } + mat->jp[mat->shape[1]] = mat->count; + for(size_t index = 0; index < mat->count; index++) { + size_t col = (size_t)NUM2SIZET(RARRAY_AREF(argv[3], index)); + size_t dest = mat->jp[col]; + + mat->ia[dest] = (size_t)NUM2SIZET(RARRAY_AREF(argv[2], index)); + mat->elements[dest] = (double)NUM2DBL(RARRAY_AREF(argv[1], index)); + + mat->jp[col] += 1; + } + + for(size_t index = 0, last = 0; index <= mat->shape[1]; index++) { + size_t temp = mat->jp[index]; + mat->jp[index] = last; + last = temp; + } + } + + return self; +} + + +VALUE csc_alloc(VALUE klass) { + csc_matrix* mat = ALLOC(csc_matrix); + + return TypedData_Wrap_Struct(klass, &csc_data_type, mat); +} + + +void csc_free(void* ptr) { + csc_matrix *mat = (csc_matrix*)ptr; + if (mat->shape) xfree(mat->shape); + if (mat->elements) xfree(mat->elements); + if (mat->ia) xfree(mat->ia); + if (mat->jp) xfree(mat->jp); + xfree(mat); +} + +size_t csc_memsize(const void* ptr) { + csc_matrix *mat = (csc_matrix*)ptr; + size_t size = sizeof(mat); + if (mat->shape) size += mat->ndims; + if (mat->elements) size += mat->count; + if (mat->ia) size += mat->count; + if (mat->jp) size += mat->shape[1]+1; + return size; +} + + +VALUE csc_get_elements(VALUE self) { + csc_matrix* input; + + TypedData_Get_Struct(self, csc_matrix, &csc_data_type, input); + + VALUE* array = ALLOC_N(VALUE, input->count); + for (size_t index = 0; index < input->count; index++) { + array[index] = DBL2NUM(input->elements[index]); + } + + return rb_ary_new4(input->count, array); +} + + +VALUE csc_get_indices(VALUE self) { + csc_matrix* input; + + TypedData_Get_Struct(self, csc_matrix, &csc_data_type, input); + + VALUE* array_ia = ALLOC_N(VALUE, input->count); + for (size_t index = 0; index < input->count; index++) { + array_ia[index] = SIZET2NUM(input->ia[index]); + } + + return rb_ary_new4(input->count, array_ia); +} + + +VALUE csc_get_indptr(VALUE self) { + csc_matrix* input; + + TypedData_Get_Struct(self, csc_matrix, &csc_data_type, input); + + VALUE* array_jp = ALLOC_N(VALUE, input->shape[1] + 1); + for (size_t index = 0; index <= input->shape[1]; index++) { + array_jp[index] = SIZET2NUM(input->jp[index]); + } + + return rb_ary_new4(input->shape[1] + 1, array_jp); +} + + +VALUE csc_get_count(VALUE self) { + csc_matrix* input; + + TypedData_Get_Struct(self, csc_matrix, &csc_data_type, input); + + return SIZET2NUM(input->count); +} + + +VALUE csc_get_ndims(VALUE self) { + csc_matrix* input; + + TypedData_Get_Struct(self, csc_matrix, &csc_data_type, input); + + return SIZET2NUM(input->ndims); +} + + +VALUE csc_get_dtype(VALUE self) { + csc_matrix* input; + + TypedData_Get_Struct(self, csc_matrix, &csc_data_type, input); + + return ID2SYM(rb_intern(DTYPE_NAMES[input->dtype])); +} + + +VALUE csc_get_shape(VALUE self) { + csc_matrix* input; + + TypedData_Get_Struct(self, csc_matrix, &csc_data_type, input); + + VALUE* array = ALLOC_N(VALUE, input->ndims); + for (size_t index = 0; index < input->ndims; index++) { + array[index] = SIZET2NUM(input->shape[index]); + } + + return rb_ary_new4(input->ndims, array); +} + +VALUE csc_to_coo(VALUE self) { + csc_matrix* left; + TypedData_Get_Struct(self, csc_matrix, &csc_data_type, left); + + coo_matrix* result = ALLOC(coo_matrix); + result->dtype = left->dtype; + result->count = left->count; + result->ndims = left->ndims; + result->shape = ALLOC_N(size_t, result->ndims); + + for (size_t index = 0; index < result->ndims; index++) { + result->shape[index] = left->shape[index]; + } + + result->elements = ALLOC_N(double, left->count); + result->ia = ALLOC_N(size_t, left->count); + result->ja = ALLOC_N(size_t, left->count); + + for(size_t j = 0; j < left->shape[1]; j++) { + for(size_t i = left->jp[j]; i < left->jp[j + 1]; i++) { + // size_t index = (left->shape[0] * j) + left->ia[i]; + result->elements[i] = left->elements[i]; + result->ia[i] = left->ia[i]; + result->ja[i] = j; + } + } + + return TypedData_Wrap_Struct(COO, &coo_data_type, result); +} + + +#include "elementwise.c" +#include "iteration.c" \ No newline at end of file diff --git a/ext/csc/elementwise.c b/ext/csc/elementwise.c new file mode 100644 index 0000000..a0f82e4 --- /dev/null +++ b/ext/csc/elementwise.c @@ -0,0 +1,192 @@ +#include "string.h" +#include "math.h" + + +double csc_perform_oper(double val_a, double val_b, char oper) { + switch(oper) { + case '+': + return (val_a + val_b); + case '-': + return (val_a - val_b); + case '*': + return (val_a * val_b); + default: //unknown operator + return 0.00; + } +} + + +/* + Takes two matrices and performs operator elementwise +*/ +VALUE csc_elementwise_binary(VALUE self, VALUE another, char oper) { + csc_matrix* left; + csc_matrix* right; + TypedData_Get_Struct(self, csc_matrix, &csc_data_type, left); + TypedData_Get_Struct(another, csc_matrix, &csc_data_type, right); + + csc_matrix* result = ALLOC(csc_matrix); + result->count = 0; + result->ndims = left->ndims; + result->shape = ALLOC_N(size_t, result->ndims); + + for (size_t index = 0; index < result->ndims; index++) { + result->shape[index] = left->shape[index]; + } + + result->elements = ALLOC_N(double, 2 * left->count); + result->ia = ALLOC_N(size_t, 2 * left->count); + result->jp = ALLOC_N(size_t, left->shape[1] + 1); + + for(size_t index = 0; index <= result->shape[1]; index++) result->jp[index] = 0; + + size_t left_index = 0, right_index = 0, result_index = 0; + for(size_t index = 1; index <= result->shape[1]; index++) { + size_t left_col_ele_count = left->jp[index]; + size_t right_col_ele_count = right->jp[index]; + + while (left_index < left_col_ele_count && right_index < right_col_ele_count) { + if (left->ia[left_index] == right->ia[right_index]) { //left and right indices equal + double result_val = csc_perform_oper(left->elements[left_index], right->elements[right_index], oper); + if(fabs(result_val) < 1e-6) { //near to zero + left_index++, right_index++; + continue; //skip current result value + } + result->elements[result_index] = result_val; + // result->ia[result_index] = left->ia[left_index]; + result->ia[result_index] = left->ia[left_index]; + + left_index++, right_index++, result_index++; + } + else if (left->ia[left_index] < right->ia[right_index]) { //left index smaller + double result_val = csc_perform_oper(left->elements[left_index], 0.0, oper); + if(fabs(result_val) < 1e-6) { //near to zero + left_index++; + continue; //skip current result value + } + result->elements[result_index] = result_val; + // result->ia[result_index] = left->ia[left_index]; + result->ia[result_index] = left->ia[left_index]; + + left_index++, result_index++; + } + else { //right index smaller + double result_val = csc_perform_oper(0.0, right->elements[right_index], oper); + if(fabs(result_val) < 1e-6) { //near to zero + right_index++; + continue; //skip current result value + } + result->elements[result_index] = result_val; + // result->ia[result_index] = right->ia[right_index]; + result->ia[result_index] = right->ia[right_index]; + + right_index++, result_index++; + } + + result->count++; + result->jp[index] += 1; + } + + while (left_index < left_col_ele_count) { + double result_val = csc_perform_oper(left->elements[left_index], 0.0, oper); + if(fabs(result_val) < 1e-6) { //near to zero + left_index++; + continue; //skip current result value + } + result->elements[result_index] = result_val; + // result->ia[result_index] = left->ia[left_index]; + result->ia[result_index] = left->ia[left_index]; + + left_index++, result_index++; + result->count++; + result->jp[index] += 1; + } + + while (right_index < right_col_ele_count) { + double result_val = csc_perform_oper(0.0, right->elements[right_index], oper); + if(fabs(result_val) < 1e-6) { //near to zero + right_index++; + continue; //skip current result value + } + result->elements[result_index] = result_val; + // result->ia[result_index] = right->ia[right_index]; + result->ia[result_index] = right->ia[right_index]; + + right_index++, result_index++; + result->count++; + result->jp[index] += 1; + } + + result->jp[index] += result->jp[index - 1]; + } + + return TypedData_Wrap_Struct(CSC, &csc_data_type, result); +} + +VALUE csc_add(VALUE self, VALUE another) { + return csc_elementwise_binary(self, another, '+'); +} + +VALUE csc_sub(VALUE self, VALUE another) { + return csc_elementwise_binary(self, another, '-'); +} + +VALUE csc_mul(VALUE self, VALUE another) { + return csc_elementwise_binary(self, another, '*'); +} + +double csc_unary_oper(double val, const char* oper) { + if (strcmp(oper, "sin") == 0) + return sin(val); + else if (strcmp(oper, "cos") == 0) + return cos(val); + else if (strcmp(oper, "tan") == 0) + return tan(val); + else + return 0.00; +} + +/* + Takes the matrix and performs unary operator elementwise +*/ +VALUE csc_elementwise_unary(VALUE self, const char* oper) { + csc_matrix* left; + TypedData_Get_Struct(self, csc_matrix, &csc_data_type, left); + + csc_matrix* result = ALLOC(csc_matrix); + result->dtype = left->dtype; + result->count = left->count; + result->ndims = left->ndims; + result->shape = ALLOC_N(size_t, result->ndims); + + for (size_t index = 0; index < result->ndims; index++) { + result->shape[index] = left->shape[index]; + } + + result->elements = ALLOC_N(double, left->count); + result->ia = ALLOC_N(size_t, left->count); + result->jp = ALLOC_N(size_t, left->shape[1] + 1); + + for(size_t index = 0; index < result->count; index++) { + result->elements[index] = csc_unary_oper(left->elements[index], oper); + result->ia[index] = left->ia[index]; + } + + for(size_t index = 0; index <= left->shape[1]; index++) { + result->jp[index] = left->jp[index]; + } + + return TypedData_Wrap_Struct(CSC, &csc_data_type, result); +} + +VALUE csc_sin(VALUE self) { + return csc_elementwise_unary(self, "sin"); +} + +VALUE csc_cos(VALUE self) { + return csc_elementwise_unary(self, "cos"); +} + +VALUE csc_tan(VALUE self) { + return csc_elementwise_unary(self, "tan"); +} diff --git a/ext/csc/iteration.c b/ext/csc/iteration.c new file mode 100644 index 0000000..a3af78e --- /dev/null +++ b/ext/csc/iteration.c @@ -0,0 +1,53 @@ +VALUE csc_each(VALUE self) { + csc_matrix* input; + + TypedData_Get_Struct(self, csc_matrix, &csc_data_type, input); + + for (size_t index = 0; index < input->count; index++) { + rb_yield(DBL2NUM(input->elements[index])); + } + + return self; +} + +VALUE csc_each_column(VALUE self) { + csc_matrix* input; + + TypedData_Get_Struct(self, csc_matrix, &csc_data_type, input); + + for (size_t j = 0; j < input->shape[1]; j++) { + size_t col_len = input->jp[j + 1] - input->jp[j]; + VALUE* col_array = ALLOC_N(VALUE, col_len); + + for(size_t i = input->jp[j]; i < input->jp[j + 1]; i++) { + col_array[i - input->jp[j]] = DBL2NUM(input->elements[i]); + } + + rb_yield(rb_ary_new4(col_len, col_array)); + } + + return self; +} + +VALUE csc_each_with_indices(VALUE self) { + csc_matrix* input; + + TypedData_Get_Struct(self, csc_matrix, &csc_data_type, input); + + VALUE* result_array = ALLOC_N(VALUE, input->ndims + 1); + + for (size_t j = 0; j < input->shape[1]; j++) { + size_t col_index = j; + for(size_t i = input->jp[j]; i < input->jp[j + 1]; i++) { + size_t row_index = input->ia[i]; + size_t index = i; + result_array[0] = DBL2NUM(input->elements[index]); + result_array[1] = INT2NUM(row_index); + result_array[2] = INT2NUM(col_index); + rb_yield(rb_ary_new4(input->ndims + 1, result_array)); + + } + } + + return self; +} \ No newline at end of file diff --git a/ext/csr/csr_def.c b/ext/csr/csr_def.c new file mode 100644 index 0000000..5d32ebe --- /dev/null +++ b/ext/csr/csr_def.c @@ -0,0 +1,189 @@ +VALUE csr_init(int argc, VALUE* argv, VALUE self) { + csr_matrix* mat; + TypedData_Get_Struct(self, csr_matrix, &csr_data_type, mat); + + if(argc > 0){ + mat->dtype = sp_float64; + mat->ndims = 2; + mat->count = (size_t)RARRAY_LEN(argv[1]); + mat->shape = ALLOC_N(size_t, mat->ndims); + for(size_t index = 0; index < mat->ndims; index++) { + mat->shape[index] = (size_t)FIX2LONG(RARRAY_AREF(argv[0], index)); + } + mat->elements = ALLOC_N(double, mat->count); + mat->ip = ALLOC_N(size_t, mat->shape[0] + 1); + mat->ja = ALLOC_N(size_t, mat->count); + for(size_t index = 0; index <= mat->shape[0]; index++) mat->ip[index] = 0; + for(size_t index = 0; index < mat->count; index++) { + size_t row_index_val = (size_t)NUM2SIZET(RARRAY_AREF(argv[2], index)); + mat->ip[row_index_val] += 1; + } + for(size_t index = 0, cumsum = 0; index < mat->shape[0]; index++) { + size_t temp = mat->ip[index]; + mat->ip[index] = cumsum; + cumsum += temp; + } + mat->ip[mat->shape[0]] = mat->count; + for(size_t index = 0; index < mat->count; index++) { + size_t row = (size_t)NUM2SIZET(RARRAY_AREF(argv[2], index)); + size_t dest = mat->ip[row]; + + mat->ja[dest] = (size_t)NUM2SIZET(RARRAY_AREF(argv[3], index)); + mat->elements[dest] = (double)NUM2DBL(RARRAY_AREF(argv[1], index)); + + mat->ip[row] += 1; + } + + for(size_t index = 0, last = 0; index <= mat->shape[0]; index++) { + size_t temp = mat->ip[index]; + mat->ip[index] = last; + last = temp; + } + } + + return self; +} + + +VALUE csr_alloc(VALUE klass) { + csr_matrix* mat = ALLOC(csr_matrix); + + return TypedData_Wrap_Struct(klass, &csr_data_type, mat); +} + + +void csr_free(void* ptr) { + csr_matrix *mat = (csr_matrix*)ptr; + if (mat->shape) xfree(mat->shape); + if (mat->elements) xfree(mat->elements); + if (mat->ip) xfree(mat->ip); + if (mat->ja) xfree(mat->ja); + xfree(mat); +} + +size_t csr_memsize(const void* ptr) { + csr_matrix *mat = (csr_matrix*)ptr; + size_t size = sizeof(mat); + if (mat->shape) size += mat->ndims; + if (mat->elements) size += mat->count; + if (mat->ip) size += mat->shape[0]+1; + if (mat->ja) size += mat->count; + return size; +} + + +VALUE csr_get_elements(VALUE self) { + csr_matrix* input; + + TypedData_Get_Struct(self, csr_matrix, &csr_data_type, input); + + VALUE* array = ALLOC_N(VALUE, input->count); + for (size_t index = 0; index < input->count; index++) { + array[index] = DBL2NUM(input->elements[index]); + } + + return rb_ary_new4(input->count, array); +} + + +VALUE csr_get_indices(VALUE self) { + csr_matrix* input; + + TypedData_Get_Struct(self, csr_matrix, &csr_data_type, input); + + VALUE* array_ja = ALLOC_N(VALUE, input->count); + for (size_t index = 0; index < input->count; index++) { + array_ja[index] = SIZET2NUM(input->ja[index]); + } + + return rb_ary_new4(input->count, array_ja); +} + + +VALUE csr_get_indptr(VALUE self) { + csr_matrix* input; + + TypedData_Get_Struct(self, csr_matrix, &csr_data_type, input); + + VALUE* array_ip = ALLOC_N(VALUE, input->shape[0] + 1); + for (size_t index = 0; index <= input->shape[0]; index++) { + array_ip[index] = SIZET2NUM(input->ip[index]); + } + + return rb_ary_new4(input->shape[0] + 1, array_ip); +} + + +VALUE csr_get_count(VALUE self) { + csr_matrix* input; + + TypedData_Get_Struct(self, csr_matrix, &csr_data_type, input); + + return SIZET2NUM(input->count); +} + + +VALUE csr_get_ndims(VALUE self) { + csr_matrix* input; + + TypedData_Get_Struct(self, csr_matrix, &csr_data_type, input); + + return SIZET2NUM(input->ndims); +} + + +VALUE csr_get_dtype(VALUE self) { + csr_matrix* input; + + TypedData_Get_Struct(self, csr_matrix, &csr_data_type, input); + + return ID2SYM(rb_intern(DTYPE_NAMES[input->dtype])); +} + + +VALUE csr_get_shape(VALUE self) { + csr_matrix* input; + + TypedData_Get_Struct(self, csr_matrix, &csr_data_type, input); + + VALUE* array = ALLOC_N(VALUE, input->ndims); + for (size_t index = 0; index < input->ndims; index++) { + array[index] = SIZET2NUM(input->shape[index]); + } + + return rb_ary_new4(input->ndims, array); +} + +VALUE csr_to_coo(VALUE self) { + csr_matrix* left; + TypedData_Get_Struct(self, csr_matrix, &csr_data_type, left); + + coo_matrix* result = ALLOC(coo_matrix); + result->dtype = left->dtype; + result->count = left->count; + result->ndims = left->ndims; + result->shape = ALLOC_N(size_t, result->ndims); + + for (size_t index = 0; index < result->ndims; index++) { + result->shape[index] = left->shape[index]; + } + + result->elements = ALLOC_N(double, left->count); + result->ia = ALLOC_N(size_t, left->count); + result->ja = ALLOC_N(size_t, left->count); + + for(size_t i = 0; i < left->shape[0]; i++) { + for(size_t j = left->ip[i]; j < left->ip[i + 1]; j++) { + // size_t index = (left->shape[1] * i) + left->ja[j]; + result->elements[j] = left->elements[j]; + result->ia[j] = i; + result->ja[j] = left->ja[j]; + } + } + + return TypedData_Wrap_Struct(COO, &coo_data_type, result); +} + + +#include "elementwise.c" +#include "iteration.c" \ No newline at end of file diff --git a/ext/csr/elementwise.c b/ext/csr/elementwise.c new file mode 100644 index 0000000..24c20f4 --- /dev/null +++ b/ext/csr/elementwise.c @@ -0,0 +1,192 @@ +#include "string.h" +#include "math.h" + + +double csr_perform_oper(double val_a, double val_b, char oper) { + switch(oper) { + case '+': + return (val_a + val_b); + case '-': + return (val_a - val_b); + case '*': + return (val_a * val_b); + default: //unknown operator + return 0.00; + } +} + + +/* + Takes two matrices and performs operator elementwise +*/ +VALUE csr_elementwise_binary(VALUE self, VALUE another, char oper) { + csr_matrix* left; + csr_matrix* right; + TypedData_Get_Struct(self, csr_matrix, &csr_data_type, left); + TypedData_Get_Struct(another, csr_matrix, &csr_data_type, right); + + csr_matrix* result = ALLOC(csr_matrix); + result->count = 0; + result->ndims = left->ndims; + result->shape = ALLOC_N(size_t, result->ndims); + + for (size_t index = 0; index < result->ndims; index++) { + result->shape[index] = left->shape[index]; + } + + result->elements = ALLOC_N(double, 2 * left->count); + result->ip = ALLOC_N(size_t, left->shape[0] + 1); + result->ja = ALLOC_N(size_t, 2 * left->count); + + for(size_t index = 0; index <= result->shape[0]; index++) result->ip[index] = 0; + + size_t left_index = 0, right_index = 0, result_index = 0; + for(size_t index = 1; index <= result->shape[0]; index++) { + size_t left_row_ele_count = left->ip[index]; + size_t right_row_ele_count = right->ip[index]; + + while (left_index < left_row_ele_count && right_index < right_row_ele_count) { + if (left->ja[left_index] == right->ja[right_index]) { //left and right indices equal + double result_val = csr_perform_oper(left->elements[left_index], right->elements[right_index], oper); + if(fabs(result_val) < 1e-6) { //near to zero + left_index++, right_index++; + continue; //skip current result value + } + result->elements[result_index] = result_val; + // result->ia[result_index] = left->ia[left_index]; + result->ja[result_index] = left->ja[left_index]; + + left_index++, right_index++, result_index++; + } + else if (left->ja[left_index] < right->ja[right_index]) { //left index smaller + double result_val = csr_perform_oper(left->elements[left_index], 0.0, oper); + if(fabs(result_val) < 1e-6) { //near to zero + left_index++; + continue; //skip current result value + } + result->elements[result_index] = result_val; + // result->ia[result_index] = left->ia[left_index]; + result->ja[result_index] = left->ja[left_index]; + + left_index++, result_index++; + } + else { //right index smaller + double result_val = csr_perform_oper(0.0, right->elements[right_index], oper); + if(fabs(result_val) < 1e-6) { //near to zero + right_index++; + continue; //skip current result value + } + result->elements[result_index] = result_val; + // result->ia[result_index] = right->ia[right_index]; + result->ja[result_index] = right->ja[right_index]; + + right_index++, result_index++; + } + + result->count++; + result->ip[index] += 1; + } + + while (left_index < left_row_ele_count) { + double result_val = csr_perform_oper(left->elements[left_index], 0.0, oper); + if(fabs(result_val) < 1e-6) { //near to zero + left_index++; + continue; //skip current result value + } + result->elements[result_index] = result_val; + // result->ia[result_index] = left->ia[left_index]; + result->ja[result_index] = left->ja[left_index]; + + left_index++, result_index++; + result->count++; + result->ip[index] += 1; + } + + while (right_index < right_row_ele_count) { + double result_val = csr_perform_oper(0.0, right->elements[right_index], oper); + if(fabs(result_val) < 1e-6) { //near to zero + right_index++; + continue; //skip current result value + } + result->elements[result_index] = result_val; + // result->ia[result_index] = right->ia[right_index]; + result->ja[result_index] = right->ja[right_index]; + + right_index++, result_index++; + result->count++; + result->ip[index] += 1; + } + + result->ip[index] += result->ip[index - 1]; + } + + return TypedData_Wrap_Struct(CSR, &csr_data_type, result); +} + +VALUE csr_add(VALUE self, VALUE another) { + return csr_elementwise_binary(self, another, '+'); +} + +VALUE csr_sub(VALUE self, VALUE another) { + return csr_elementwise_binary(self, another, '-'); +} + +VALUE csr_mul(VALUE self, VALUE another) { + return csr_elementwise_binary(self, another, '*'); +} + +double csr_unary_oper(double val, const char* oper) { + if (strcmp(oper, "sin") == 0) + return sin(val); + else if (strcmp(oper, "cos") == 0) + return cos(val); + else if (strcmp(oper, "tan") == 0) + return tan(val); + else + return 0.00; +} + +/* + Takes the matrix and performs unary operator elementwise +*/ +VALUE csr_elementwise_unary(VALUE self, const char* oper) { + csr_matrix* left; + TypedData_Get_Struct(self, csr_matrix, &csr_data_type, left); + + csr_matrix* result = ALLOC(csr_matrix); + result->dtype = left->dtype; + result->count = left->count; + result->ndims = left->ndims; + result->shape = ALLOC_N(size_t, result->ndims); + + for (size_t index = 0; index < result->ndims; index++) { + result->shape[index] = left->shape[index]; + } + + result->elements = ALLOC_N(double, left->count); + result->ip = ALLOC_N(size_t, left->shape[0] + 1); + result->ja = ALLOC_N(size_t, left->count); + + for(size_t index = 0; index < result->count; index++) { + result->elements[index] = csr_unary_oper(left->elements[index], oper); + result->ja[index] = left->ja[index]; + } + + for(size_t index = 0; index <= left->shape[0]; index++) { + result->ip[index] = left->ip[index]; + } + + return TypedData_Wrap_Struct(CSR, &csr_data_type, result); +} + +VALUE csr_sin(VALUE self) { + return csr_elementwise_unary(self, "sin"); +} + +VALUE csr_cos(VALUE self) { + return csr_elementwise_unary(self, "cos"); +} + +VALUE csr_tan(VALUE self) { + return csr_elementwise_unary(self, "tan"); +} diff --git a/ext/csr/iteration.c b/ext/csr/iteration.c new file mode 100644 index 0000000..c1a4781 --- /dev/null +++ b/ext/csr/iteration.c @@ -0,0 +1,53 @@ +VALUE csr_each(VALUE self) { + csr_matrix* input; + + TypedData_Get_Struct(self, csr_matrix, &csr_data_type, input); + + for (size_t index = 0; index < input->count; index++) { + rb_yield(DBL2NUM(input->elements[index])); + } + + return self; +} + +VALUE csr_each_row(VALUE self) { + csr_matrix* input; + + TypedData_Get_Struct(self, csr_matrix, &csr_data_type, input); + + for (size_t i = 0; i < input->shape[0]; i++) { + size_t row_len = input->ip[i + 1] - input->ip[i]; + VALUE* row_array = ALLOC_N(VALUE, row_len); + + for(size_t j = input->ip[i]; j < input->ip[i + 1]; j++) { + row_array[j - input->ip[i]] = DBL2NUM(input->elements[j]); + } + + rb_yield(rb_ary_new4(row_len, row_array)); + } + + return self; +} + +VALUE csr_each_with_indices(VALUE self) { + csr_matrix* input; + + TypedData_Get_Struct(self, csr_matrix, &csr_data_type, input); + + VALUE* result_array = ALLOC_N(VALUE, input->ndims + 1); + + for (size_t i = 0; i < input->shape[0]; i++) { + size_t row_index = i; + for(size_t j = input->ip[i]; j < input->ip[i + 1]; j++) { + size_t col_index = input->ja[j]; + size_t index = j; + result_array[0] = DBL2NUM(input->elements[index]); + result_array[1] = INT2NUM(row_index); + result_array[2] = INT2NUM(col_index); + rb_yield(rb_ary_new4(input->ndims + 1, result_array)); + + } + } + + return self; +} \ No newline at end of file diff --git a/ext/dia/dia_def.c b/ext/dia/dia_def.c index 3bd4072..3b590df 100644 --- a/ext/dia/dia_def.c +++ b/ext/dia/dia_def.c @@ -1,6 +1,6 @@ VALUE dia_init(int argc, VALUE* argv, VALUE self) { dia_matrix* mat; - Data_Get_Struct(self, dia_matrix, mat); + TypedData_Get_Struct(self, dia_matrix, &dia_data_type, mat); if(argc > 0){ mat->dtype = sp_float64; @@ -23,19 +23,30 @@ VALUE dia_init(int argc, VALUE* argv, VALUE self) { VALUE dia_alloc(VALUE klass) { dia_matrix* mat = ALLOC(dia_matrix); - return Data_Wrap_Struct(klass, NULL, dia_free, mat); + return TypedData_Wrap_Struct(klass, &dia_data_type, mat); } -void dia_free(dia_matrix* mat) { +void dia_free(void* ptr) { + dia_matrix *mat = (dia_matrix*)ptr; + if (mat->shape) xfree(mat->shape); + if (mat->elements) xfree(mat->elements); xfree(mat); } +size_t dia_memsize(const void* ptr) { + dia_matrix *mat = (dia_matrix*)ptr; + size_t size = sizeof(mat); + if (mat->shape) size += mat->ndims; + if (mat->elements) size += mat->count; + return size; +} + VALUE dia_get_elements(VALUE self) { dia_matrix* input; - Data_Get_Struct(self, dia_matrix, input); + TypedData_Get_Struct(self, dia_matrix, &dia_data_type, input); VALUE* array = ALLOC_N(VALUE, input->count); for (size_t index = 0; index < input->count; index++) { @@ -49,7 +60,7 @@ VALUE dia_get_elements(VALUE self) { VALUE dia_get_count(VALUE self) { dia_matrix* input; - Data_Get_Struct(self, dia_matrix, input); + TypedData_Get_Struct(self, dia_matrix, &dia_data_type, input); return SIZET2NUM(input->count); } @@ -58,7 +69,7 @@ VALUE dia_get_count(VALUE self) { VALUE dia_get_ndims(VALUE self) { dia_matrix* input; - Data_Get_Struct(self, dia_matrix, input); + TypedData_Get_Struct(self, dia_matrix, &dia_data_type, input); return SIZET2NUM(input->ndims); } @@ -67,7 +78,7 @@ VALUE dia_get_ndims(VALUE self) { VALUE dia_get_dtype(VALUE self) { dia_matrix* input; - Data_Get_Struct(self, dia_matrix, input); + TypedData_Get_Struct(self, dia_matrix, &dia_data_type, input); return ID2SYM(rb_intern(DTYPE_NAMES[input->dtype])); } @@ -76,7 +87,7 @@ VALUE dia_get_dtype(VALUE self) { VALUE dia_get_shape(VALUE self) { dia_matrix* input; - Data_Get_Struct(self, dia_matrix, input); + TypedData_Get_Struct(self, dia_matrix, &dia_data_type, input); VALUE* array = ALLOC_N(VALUE, input->ndims); for (size_t index = 0; index < input->ndims; index++) { @@ -87,4 +98,5 @@ VALUE dia_get_shape(VALUE self) { } -#include "elementwise.c" \ No newline at end of file +#include "elementwise.c" +#include "iteration.c" \ No newline at end of file diff --git a/ext/dia/elementwise.c b/ext/dia/elementwise.c index 7e14398..4fb8839 100644 --- a/ext/dia/elementwise.c +++ b/ext/dia/elementwise.c @@ -1,3 +1,7 @@ +#include "string.h" +#include "math.h" + + double dia_perform_oper(double val_a, double val_b, char oper) { switch(oper) { case '+': @@ -18,8 +22,8 @@ double dia_perform_oper(double val_a, double val_b, char oper) { VALUE dia_elementwise_binary(VALUE self, VALUE another, char oper) { dia_matrix* left; dia_matrix* right; - Data_Get_Struct(self, dia_matrix, left); - Data_Get_Struct(another, dia_matrix, right); + TypedData_Get_Struct(self, dia_matrix, &dia_data_type, left); + TypedData_Get_Struct(another, dia_matrix, &dia_data_type, right); dia_matrix* result = ALLOC(dia_matrix); result->count = left->count; @@ -37,7 +41,7 @@ VALUE dia_elementwise_binary(VALUE self, VALUE another, char oper) { result->elements[result_index] = result_val; } - return Data_Wrap_Struct(DIA, NULL, dia_free, result); + return TypedData_Wrap_Struct(DIA, &dia_data_type, result); } VALUE dia_add(VALUE self, VALUE another) { @@ -50,4 +54,53 @@ VALUE dia_sub(VALUE self, VALUE another) { VALUE dia_mul(VALUE self, VALUE another) { return dia_elementwise_binary(self, another, '*'); -} \ No newline at end of file +} + +double dia_unary_oper(double val, const char* oper) { + if (strcmp(oper, "sin") == 0) + return sin(val); + else if (strcmp(oper, "cos") == 0) + return cos(val); + else if (strcmp(oper, "tan") == 0) + return tan(val); + else + return 0.00; +} + +/* + Takes the matrix and performs unary operator elementwise +*/ +VALUE dia_elementwise_unary(VALUE self, const char* oper) { + dia_matrix* left; + TypedData_Get_Struct(self, dia_matrix, &dia_data_type, left); + + dia_matrix* result = ALLOC(dia_matrix); + result->dtype = left->dtype; + result->count = left->count; + result->ndims = left->ndims; + result->shape = ALLOC_N(size_t, result->ndims); + + for (size_t index = 0; index < result->ndims; index++) { + result->shape[index] = left->shape[index]; + } + + result->elements = ALLOC_N(double, left->count); + + for(size_t index = 0; index < result->count; index++) { + result->elements[index] = dia_unary_oper(left->elements[index], oper); + } + + return TypedData_Wrap_Struct(DIA, &dia_data_type, result); +} + +VALUE dia_sin(VALUE self) { + return dia_elementwise_unary(self, "sin"); +} + +VALUE dia_cos(VALUE self) { + return dia_elementwise_unary(self, "cos"); +} + +VALUE dia_tan(VALUE self) { + return dia_elementwise_unary(self, "tan"); +} diff --git a/ext/dia/iteration.c b/ext/dia/iteration.c new file mode 100644 index 0000000..e27ff3a --- /dev/null +++ b/ext/dia/iteration.c @@ -0,0 +1,11 @@ +VALUE dia_each(VALUE self) { + dia_matrix* input; + + TypedData_Get_Struct(self, dia_matrix, &dia_data_type, input); + + for (size_t index = 0; index < input->count; index++) { + rb_yield(DBL2NUM(input->elements[index])); + } + + return self; +} \ No newline at end of file diff --git a/ext/ruby_sparse.c b/ext/ruby_sparse.c index 20eb943..98d2662 100644 --- a/ext/ruby_sparse.c +++ b/ext/ruby_sparse.c @@ -41,6 +41,67 @@ typedef struct COO_STRUCT size_t* ja; //col index }coo_matrix; +void coo_free(void* ptr); +size_t coo_memsize(const void* ptr); + +static const rb_data_type_t coo_data_type = { + "ruby-sparse/coo", + { + 0, + coo_free, + coo_memsize, + }, + 0, 0, RUBY_TYPED_FREE_IMMEDIATELY +}; + +typedef struct CSR_STRUCT +{ + sp_dtype dtype; + size_t ndims; + size_t count; //count of non-zero elements + size_t* shape; + double* elements; //elements array + size_t* ip; //row pointer vals + size_t* ja; //col index +}csr_matrix; + +void csr_free(void* ptr); +size_t csr_memsize(const void* ptr); + +static const rb_data_type_t csr_data_type = { + "ruby-sparse/csr", + { + 0, + csr_free, + csr_memsize, + }, + 0, 0, RUBY_TYPED_FREE_IMMEDIATELY +}; + +typedef struct CSC_STRUCT +{ + sp_dtype dtype; + size_t ndims; + size_t count; //count of non-zero elements + size_t* shape; + double* elements; //elements array + size_t* ia; //row index + size_t* jp; //col pointer vals +}csc_matrix; + +void csc_free(void* ptr); +size_t csc_memsize(const void* ptr); + +static const rb_data_type_t csc_data_type = { + "ruby-sparse/csc", + { + 0, + csc_free, + csc_memsize, + }, + 0, 0, RUBY_TYPED_FREE_IMMEDIATELY +}; + typedef struct DIA_STRUCT { sp_dtype dtype; @@ -50,6 +111,19 @@ typedef struct DIA_STRUCT double* elements; //elements array }dia_matrix; +void dia_free(void* ptr); +size_t dia_memsize(const void* ptr); + +static const rb_data_type_t dia_data_type = { + "ruby-sparse/dia", + { + 0, + dia_free, + dia_memsize, + }, + 0, 0, RUBY_TYPED_FREE_IMMEDIATELY +}; + // typedef struct GCXS_STRUCT // { // sp_dtype dtype; @@ -59,11 +133,13 @@ typedef struct DIA_STRUCT // double* elements; //elements array // size_t* ia; //row index // size_t* ja; //col index -// }coo_matrix; +// }gcxs_matrix; VALUE RubySparse = Qnil; VALUE SparseArray = Qnil; VALUE COO = Qnil; +VALUE CSR = Qnil; +VALUE CSC = Qnil; VALUE DIA = Qnil; // VALUE GCXS = Qnil; @@ -78,14 +154,77 @@ VALUE coo_get_coords(VALUE self); VALUE coo_get_count(VALUE self); VALUE coo_get_ndims(VALUE self); VALUE coo_alloc(VALUE klass); -void coo_free(coo_matrix* mat); + +VALUE coo_to_csr(VALUE self); +VALUE coo_to_csc(VALUE self); VALUE coo_add(VALUE self, VALUE another); VALUE coo_sub(VALUE self, VALUE another); VALUE coo_mul(VALUE self, VALUE another); +VALUE coo_sin(VALUE self); +VALUE coo_cos(VALUE self); +VALUE coo_tan(VALUE self); + +VALUE coo_each(VALUE self); +VALUE coo_each_with_indices(VALUE self); + VALUE coo_from_nmatrix(VALUE self, VALUE nmat); +// csr methods declaration +VALUE csr_init(int argc, VALUE* argv, VALUE self); +VALUE csr_get_dtype(VALUE self); +VALUE csr_get_shape(VALUE self); +VALUE csr_get_elements(VALUE self); +VALUE csr_get_indices(VALUE self); +VALUE csr_get_indptr(VALUE self); +VALUE csr_get_count(VALUE self); +VALUE csr_get_ndims(VALUE self); +VALUE csr_alloc(VALUE klass); + +VALUE csr_to_coo(VALUE self); + +VALUE csr_add(VALUE self, VALUE another); +VALUE csr_sub(VALUE self, VALUE another); +VALUE csr_mul(VALUE self, VALUE another); + +VALUE csr_sin(VALUE self); +VALUE csr_cos(VALUE self); +VALUE csr_tan(VALUE self); + +VALUE csr_each(VALUE self); +VALUE csr_each_row(VALUE self); +VALUE csr_each_with_indices(VALUE self); + +VALUE csr_from_nmatrix(VALUE self, VALUE nmat); + +// csc methods declaration +VALUE csc_init(int argc, VALUE* argv, VALUE self); +VALUE csc_get_dtype(VALUE self); +VALUE csc_get_shape(VALUE self); +VALUE csc_get_elements(VALUE self); +VALUE csc_get_indices(VALUE self); +VALUE csc_get_indptr(VALUE self); +VALUE csc_get_count(VALUE self); +VALUE csc_get_ndims(VALUE self); +VALUE csc_alloc(VALUE klass); + +VALUE csc_to_coo(VALUE self); + +VALUE csc_add(VALUE self, VALUE another); +VALUE csc_sub(VALUE self, VALUE another); +VALUE csc_mul(VALUE self, VALUE another); + +VALUE csc_sin(VALUE self); +VALUE csc_cos(VALUE self); +VALUE csc_tan(VALUE self); + +VALUE csc_each(VALUE self); +VALUE csc_each_column(VALUE self); +VALUE csc_each_with_indices(VALUE self); + +VALUE csc_from_nmatrix(VALUE self, VALUE nmat); + // dia methods declaration VALUE dia_init(int argc, VALUE* argv, VALUE self); VALUE dia_get_dtype(VALUE self); @@ -94,12 +233,17 @@ VALUE dia_get_elements(VALUE self); VALUE dia_get_count(VALUE self); VALUE dia_get_ndims(VALUE self); VALUE dia_alloc(VALUE klass); -void dia_free(dia_matrix* mat); VALUE dia_add(VALUE self, VALUE another); VALUE dia_sub(VALUE self, VALUE another); VALUE dia_mul(VALUE self, VALUE another); +VALUE dia_sin(VALUE self); +VALUE dia_cos(VALUE self); +VALUE dia_tan(VALUE self); + +VALUE dia_each(VALUE self); + VALUE dia_from_nmatrix(VALUE self, VALUE nmat); @@ -108,6 +252,8 @@ void Init_ruby_sparse() { SparseArray = rb_define_class_under(RubySparse, "SparseArray", rb_cObject); COO = rb_define_class_under(RubySparse, "COO", SparseArray); + CSR = rb_define_class_under(RubySparse, "CSR", SparseArray); + CSC = rb_define_class_under(RubySparse, "CSC", SparseArray); DIA = rb_define_class_under(RubySparse, "DIA", SparseArray); // GCXS = rb_define_class_under(RubySparse, "GCXS", SparseArray); @@ -120,13 +266,77 @@ void Init_ruby_sparse() { rb_define_method(COO, "nzcount", coo_get_count, 0); rb_define_method(COO, "dim", coo_get_ndims, 0); + rb_define_method(COO, "to_csr", coo_to_csr, 0); + rb_define_method(COO, "to_csc", coo_to_csc, 0); + rb_define_method(COO, "+", coo_add, 1); rb_define_method(COO, "-", coo_sub, 1); rb_define_method(COO, "*", coo_mul, 1); + rb_define_method(COO, "sin", coo_sin, 0); + rb_define_method(COO, "cos", coo_cos, 0); + rb_define_method(COO, "tan", coo_tan, 0); + + rb_define_method(COO, "each", coo_each, 0); + rb_define_method(COO, "each_with_indices", coo_each_with_indices, 0); + //rb_define_singleton_method(COO, "from_nmatrix", coo_from_nmatrix, 1); + rb_define_alloc_func(CSR, csr_alloc); + rb_define_method(CSR, "initialize", csr_init, -1); + rb_define_method(CSR, "dtype", csr_get_dtype, 0); + rb_define_method(CSR, "shape", csr_get_shape, 0); + rb_define_method(CSR, "elements", csr_get_elements, 0); + rb_define_method(CSR, "indices", csr_get_indices, 0); + rb_define_method(CSR, "indptr", csr_get_indptr, 0); + rb_define_method(CSR, "nzcount", csr_get_count, 0); + rb_define_method(CSR, "dim", csr_get_ndims, 0); + + rb_define_method(CSR, "to_coo", csr_to_coo, 0); + + rb_define_method(CSR, "+", csr_add, 1); + rb_define_method(CSR, "-", csr_sub, 1); + rb_define_method(CSR, "*", csr_mul, 1); + + rb_define_method(CSR, "sin", csr_sin, 0); + rb_define_method(CSR, "cos", csr_cos, 0); + rb_define_method(CSR, "tan", csr_tan, 0); + + rb_define_method(CSR, "each", csr_each, 0); + rb_define_method(CSR, "each_row", csr_each_row, 0); + rb_define_method(CSR, "each_with_indices", csr_each_with_indices, 0); + + //rb_define_singleton_method(CSR, "from_nmatrix", csr_from_nmatrix, 1); + + + rb_define_alloc_func(CSC, csc_alloc); + rb_define_method(CSC, "initialize", csc_init, -1); + rb_define_method(CSC, "dtype", csc_get_dtype, 0); + rb_define_method(CSC, "shape", csc_get_shape, 0); + rb_define_method(CSC, "elements", csc_get_elements, 0); + rb_define_method(CSC, "indices", csc_get_indices, 0); + rb_define_method(CSC, "indptr", csc_get_indptr, 0); + rb_define_method(CSC, "nzcount", csc_get_count, 0); + rb_define_method(CSC, "dim", csc_get_ndims, 0); + + rb_define_method(CSC, "to_coo", csc_to_coo, 0); + + rb_define_method(CSC, "+", csc_add, 1); + rb_define_method(CSC, "-", csc_sub, 1); + rb_define_method(CSC, "*", csc_mul, 1); + + rb_define_method(CSC, "sin", csc_sin, 0); + rb_define_method(CSC, "cos", csc_cos, 0); + rb_define_method(CSC, "tan", csc_tan, 0); + + rb_define_method(CSC, "each", csc_each, 0); + rb_define_method(CSC, "each_column", csc_each_column, 0); + rb_define_method(CSC, "each_with_indices", csc_each_with_indices, 0); + + //rb_define_singleton_method(CSC, "from_nmatrix", csc_from_nmatrix, 1); + + rb_define_alloc_func(DIA, dia_alloc); rb_define_method(DIA, "initialize", dia_init, -1); rb_define_method(DIA, "dtype", dia_get_dtype, 0); @@ -138,10 +348,18 @@ void Init_ruby_sparse() { rb_define_method(DIA, "-", dia_sub, 1); rb_define_method(DIA, "*", dia_mul, 1); + rb_define_method(DIA, "sin", dia_sin, 0); + rb_define_method(DIA, "cos", dia_cos, 0); + rb_define_method(DIA, "tan", dia_tan, 0); + + rb_define_method(DIA, "each", dia_each, 0); + //rb_define_singleton_method(DIA, "from_nmatrix", dia_from_nmatrix, 1); } #include "coo/coo_def.c" +#include "csr/csr_def.c" +#include "csc/csc_def.c" #include "dia/dia_def.c" #include "interfaces/nmatrix.c" diff --git a/lib/ruby_sparse/ruby_sparse.rb b/lib/ruby_sparse/ruby_sparse.rb index 8d06620..c6760f5 100644 --- a/lib/ruby_sparse/ruby_sparse.rb +++ b/lib/ruby_sparse/ruby_sparse.rb @@ -20,6 +20,16 @@ class COO < SparseArray # @elements = elements # end + def pretty_print(q) + self.inspect.pretty_print(q) + end + + def inspect #:nodoc: + original_inspect = super() + original_inspect = original_inspect[0...original_inspect.size-1] + original_inspect + "; " + inspect_helper.join("; ") + ">" + end + def self.from_nmatrix(nmat) nm_elements = nmat.elements nm_shape = nmat.shape @@ -62,6 +72,44 @@ def to_nmatrix m = NMatrix.new self.shape, nm_elements return m end + + def _dump data + [ + dim, + nzcount, + dtype, + shape, + elements, + coords[0], + coords[1], + ].join ":" + end + + def self._load args + values = args.split(":") + dim = values[0].to_i + nzcount = values[1].to_i + dtype = values[2].to_sym + shape = values[3...(3+dim)] + elements = values[(3+dim)...(3+dim+nzcount)] + ia = values[(3+dim+nzcount)...(3+dim+(2*nzcount))] + ja = values[(3+dim+(2*nzcount))...(3+dim+(3*nzcount))] + + (0...dim).each { |index| shape[index] = shape[index].to_i } + (0...nzcount).each { |index| elements[index] = elements[index].to_f } + (0...nzcount).each { |index| ia[index] = ia[index].to_i } + (0...nzcount).each { |index| ja[index] = ja[index].to_i } + + self.new(shape, elements, ia, ja) + end + + protected + + def inspect_helper #:nodoc: + ary = [] + ary << "shape: [#{shape.join(',')}]" << "dtype: #{dtype}" << "nnz: #{nzcount}" + ary + end end @@ -74,5 +122,296 @@ class DIA < SparseArray # @elements = elements # end + def pretty_print(q) + self.inspect.pretty_print(q) + end + + def inspect #:nodoc: + original_inspect = super() + original_inspect = original_inspect[0...original_inspect.size-1] + original_inspect + "; " + inspect_helper.join("; ") + ">" + end + + def self.from_nmatrix(nmat) + nm_elements = nmat.elements + nm_shape = nmat.shape + + if nmat.dim != 2 + raise StandardError.new "NMatrix must be of 2 dimensions." + end + + dia_elements = [] + + for i in (0...nmat.shape[0]) + nm_index = (nmat.shape[1] * i) + i + nm_value = nm_elements[nm_index] + dia_elements.append(nm_value) + end + + dia_mat = self.new nm_shape, dia_elements + return dia_mat + end + + def to_nmatrix + nm_elements = Array.new(self.shape[0]*self.shape[1], 0) + for i in (0...self.shape[0]) + nm_index = (self.shape[1] * i) + i + nm_elements[nm_index] = self.elements[i] + end + + m = NMatrix.new self.shape, nm_elements + return m + end + + def _dump data + [ + dim, + dtype, + shape, + elements, + ].join ":" + end + + def self._load args + values = args.split(":") + dim = values[0].to_i + dtype = values[1].to_sym + shape = values[2...(2+dim)] + elements = values[(2+dim)...(2+dim+shape[0].to_i)] + + (0...dim).each { |index| shape[index] = shape[index].to_i } + (0...shape[0]).each { |index| elements[index] = elements[index].to_f } + + self.new(shape, elements) + end + + protected + + def inspect_helper #:nodoc: + ary = [] + ary << "shape: [#{shape.join(',')}]" << "dtype: #{dtype}" << "nnz: #{nzcount}" + ary + end + + end + + class CSR < SparseArray + + # attr_reader :shape, :elements + + # def initialize shape, elements + # @shape = shape + # @elements = elements + # end + + def pretty_print(q) + self.inspect.pretty_print(q) + end + + def inspect #:nodoc: + original_inspect = super() + original_inspect = original_inspect[0...original_inspect.size-1] + original_inspect + "; " + inspect_helper.join("; ") + ">" + end + + def self.from_nmatrix(nmat) + nm_elements = nmat.elements + nm_shape = nmat.shape + + if nmat.dim != 2 + raise StandardError.new "NMatrix must be of 2 dimensions." + end + + csr_elements = [] + csr_ia = [] + csr_ja = [] + + for i in (0...nmat.shape[0]) + for j in (0...nmat.shape[1]) + nm_index = (nmat.shape[1] * i) + j + nm_value = nm_elements[nm_index] + if nm_value == 0 + next + end + csr_elements.append(nm_value) + csr_ia.append(i) + csr_ja.append(j) + end + end + + csr_mat = self.new nm_shape, csr_elements, csr_ia, csr_ja + return csr_mat + end + + def to_nmatrix + ip = self.indptr + ja = self.indices + + nm_elements = Array.new(self.shape[0]*self.shape[1], 0) + for i in (0...self.shape[0]) + for j in (ip[i]...ip[i + 1]) + nm_index = (self.shape[1] * i) + ja[j] + nm_elements[nm_index] = self.elements[j] + end + end + + m = NMatrix.new self.shape, nm_elements + return m + end + + def to_csc + return self.to_coo.to_csc + end + + def _dump data + [ + dim, + nzcount, + dtype, + shape, + elements, + indptr, + indices, + ].join ":" + end + + def self._load args + values = args.split(":") + dim = values[0].to_i + nzcount = values[1].to_i + dtype = values[2].to_sym + shape = values[3...(3+dim)] + elements = values[(3+dim)...(3+dim+nzcount)] + ip = values[(3+dim+nzcount)...(3+dim+nzcount+shape[0].to_i)] + ja = values[(3+dim+nzcount+shape[0].to_i)...(3+dim+(2*nzcount)+shape[0].to_i)] + + (0...dim).each { |index| shape[index] = shape[index].to_i } + (0...nzcount).each { |index| elements[index] = elements[index].to_f } + (0...shape[0]).each { |index| ip[index] = ip[index].to_i } + (0...nzcount).each { |index| ja[index] = ja[index].to_i } + + # convert ip to ia here and then pass ia below instead + + self.new(shape, elements, ip, ja) + end + + protected + + def inspect_helper #:nodoc: + ary = [] + ary << "shape: [#{shape.join(',')}]" << "dtype: #{dtype}" << "nnz: #{nzcount}" + ary + end + + end + + class CSC < SparseArray + + # attr_reader :shape, :elements + + # def initialize shape, elements + # @shape = shape + # @elements = elements + # end + + def pretty_print(q) + self.inspect.pretty_print(q) + end + + def inspect #:nodoc: + original_inspect = super() + original_inspect = original_inspect[0...original_inspect.size-1] + original_inspect + "; " + inspect_helper.join("; ") + ">" + end + + def self.from_nmatrix(nmat) + nm_elements = nmat.elements + nm_shape = nmat.shape + + if nmat.dim != 2 + raise StandardError.new "NMatrix must be of 2 dimensions." + end + + csc_elements = [] + csc_ia = [] + csc_ja = [] + + for i in (0...nmat.shape[0]) + for j in (0...nmat.shape[1]) + nm_index = (nmat.shape[1] * i) + j + nm_value = nm_elements[nm_index] + if nm_value == 0 + next + end + csc_elements.append(nm_value) + csc_ia.append(i) + csc_ja.append(j) + end + end + + csc_mat = self.new nm_shape, csc_elements, csc_ia, csc_ja + return csc_mat + end + + def to_nmatrix + ia = self.indices + jp = self.indptr + + nm_elements = Array.new(self.shape[0]*self.shape[1], 0) + for j in (0...self.shape[1]) + for i in (jp[j]...jp[j + 1]) + nm_index = (self.shape[0] * j) + ia[i] + nm_elements[nm_index] = self.elements[i] + end + end + + m = NMatrix.new self.shape, nm_elements + return m + end + + def to_csr + return self.to_coo.to_csr + end + + def _dump data + [ + dim, + nzcount, + dtype, + shape, + elements, + indices, + indptr, + ].join ":" + end + + def self._load args + values = args.split(":") + dim = values[0].to_i + nzcount = values[1].to_i + dtype = values[2].to_sym + shape = values[3...(3+dim)] + elements = values[(3+dim)...(3+dim+nzcount)] + ia = values[(3+dim+nzcount)...(3+dim+(2*nzcount))] + jp = values[(3+dim+(2*nzcount))...(3+dim+(2*nzcount)+shape[1].to_i)] + + (0...dim).each { |index| shape[index] = shape[index].to_i } + (0...nzcount).each { |index| elements[index] = elements[index].to_f } + (0...nzcount).each { |index| ia[index] = ia[index].to_i } + (0...shape[1]).each { |index| jp[index] = jp[index].to_i } + + # convert jp to ja here and then pass ja below instead + + self.new(shape, elements, ia, jp) + end + + protected + + def inspect_helper #:nodoc: + ary = [] + ary << "shape: [#{shape.join(',')}]" << "dtype: #{dtype}" << "nnz: #{nzcount}" + ary + end + end end \ No newline at end of file diff --git a/test/elementwise_test.rb b/test/elementwise_test.rb index 718ea4e..67e7e03 100644 --- a/test/elementwise_test.rb +++ b/test/elementwise_test.rb @@ -37,6 +37,16 @@ def test_mul assert_equal answer.coords, result.coords end + def test_sin + result = RubySparse::COO.new [3, 3], [0.8414709848078965, 0.9092974268256817, 0.1411200080598672], [0, 1, 2], [0, 1, 2] + answer = @left.sin + assert_equal answer.dim, result.dim + assert_equal answer.shape, result.shape + assert_equal answer.nzcount, result.nzcount + assert_equal answer.elements, result.elements + assert_equal answer.coords, result.coords + end + end class RubySparse::DIA::ElementWiseTest < Minitest::Test @@ -70,4 +80,118 @@ def test_mul assert_equal answer.elements, result.elements end -end \ No newline at end of file + def test_sin + result = RubySparse::DIA.new [3, 3], [0.8414709848078965, 0.0, 0.1411200080598672] + answer = @left.sin + assert_equal answer.dim, result.dim + assert_equal answer.shape, result.shape + assert_equal answer.elements, result.elements + end + +end + +class RubySparse::CSR::ElementWiseTest < Minitest::Test + + def setup + @left = RubySparse::CSR.new [3, 3], [1, 2, 3], [0, 1, 2], [0, 1, 2] + @right = RubySparse::CSR.new [3, 3], [3, 2, 3], [0, 1, 2], [2, 2, 2] + end + + def test_add + result = RubySparse::CSR.new [3, 3], [1, 3, 2, 2, 6], [0, 0, 1, 1, 2], [0, 2, 1, 2, 2] + answer = @left + @right + assert_equal answer.dim, result.dim + assert_equal answer.shape, result.shape + assert_equal answer.nzcount, result.nzcount + assert_equal answer.elements, result.elements + assert_equal answer.indices, result.indices + assert_equal answer.indptr, result.indptr + end + + def test_sub + result = RubySparse::CSR.new [3, 3], [1, -3, 2, -2], [0, 0, 1, 1], [0, 2, 1, 2] + answer = @left - @right + assert_equal answer.dim, result.dim + assert_equal answer.shape, result.shape + assert_equal answer.nzcount, result.nzcount + assert_equal answer.elements, result.elements + assert_equal answer.indices, result.indices + assert_equal answer.indptr, result.indptr + end + + def test_mul + result = RubySparse::CSR.new [3, 3], [9], [2], [2] + answer = @left * @right + assert_equal answer.dim, result.dim + assert_equal answer.shape, result.shape + assert_equal answer.nzcount, result.nzcount + assert_equal answer.elements, result.elements + assert_equal answer.indices, result.indices + assert_equal answer.indptr, result.indptr + end + + def test_sin + result = RubySparse::CSR.new [3, 3], [0.8414709848078965, 0.9092974268256817, 0.1411200080598672], [0, 1, 2], [0, 1, 2] + answer = @left.sin + assert_equal answer.dim, result.dim + assert_equal answer.shape, result.shape + assert_equal answer.nzcount, result.nzcount + assert_equal answer.elements, result.elements + assert_equal answer.indices, result.indices + assert_equal answer.indptr, result.indptr + end + +end + +class RubySparse::CSC::ElementWiseTest < Minitest::Test + + def setup + @left = RubySparse::CSC.new [3, 3], [1, 2, 3], [0, 1, 2], [0, 1, 2] + @right = RubySparse::CSC.new [3, 3], [3, 2, 3], [0, 1, 2], [2, 2, 2] + end + + def test_add + result = RubySparse::CSC.new [3, 3], [1, 3, 2, 2, 6], [0, 0, 1, 1, 2], [0, 2, 1, 2, 2] + answer = @left + @right + assert_equal answer.dim, result.dim + assert_equal answer.shape, result.shape + assert_equal answer.nzcount, result.nzcount + assert_equal answer.elements, result.elements + assert_equal answer.indices, result.indices + assert_equal answer.indptr, result.indptr + end + + def test_sub + result = RubySparse::CSC.new [3, 3], [1, -3, 2, -2], [0, 0, 1, 1], [0, 2, 1, 2] + answer = @left - @right + assert_equal answer.dim, result.dim + assert_equal answer.shape, result.shape + assert_equal answer.nzcount, result.nzcount + assert_equal answer.elements, result.elements + assert_equal answer.indices, result.indices + assert_equal answer.indptr, result.indptr + end + + def test_mul + result = RubySparse::CSC.new [3, 3], [9], [2], [2] + answer = @left * @right + assert_equal answer.dim, result.dim + assert_equal answer.shape, result.shape + assert_equal answer.nzcount, result.nzcount + assert_equal answer.elements, result.elements + assert_equal answer.indices, result.indices + assert_equal answer.indptr, result.indptr + end + + def test_sin + result = RubySparse::CSC.new [3, 3], [0.8414709848078965, 0.9092974268256817, 0.1411200080598672], [0, 1, 2], [0, 1, 2] + answer = @left.sin + assert_equal answer.dim, result.dim + assert_equal answer.shape, result.shape + assert_equal answer.nzcount, result.nzcount + assert_equal answer.elements, result.elements + assert_equal answer.indices, result.indices + assert_equal answer.indptr, result.indptr + end + +end diff --git a/test/sparse_test.rb b/test/sparse_test.rb index 0362a52..cdadb0d 100644 --- a/test/sparse_test.rb +++ b/test/sparse_test.rb @@ -46,4 +46,70 @@ def test_elements assert_equal [1, 0, 3], @n.elements end +end + +class RubySparse::CSR::CreationTest < Minitest::Test + + def setup + @n = RubySparse::CSR.new [3, 3], [1, 2, 3], [0, 1, 2], [0, 1, 2] + @m = RubySparse::CSR.new [4, 6], [1, 2, 3, 4, 5, 6, 7, 8], [0, 0, 1, 1, 2, 2, 2, 3], [0, 1, 1, 3, 2, 3, 4, 5] + end + + def test_ndims + assert_equal 2, @n.dim + end + + def test_shape + assert_equal [3, 3], @n.shape + end + + def test_elements + assert_equal [1, 2, 3], @n.elements + end + + def test_indices + assert_equal [0, 1, 1, 3, 2, 3, 4, 5], @m.indices + end + + def test_indptr + assert_equal [0, 2, 4, 7, 8], @m.indptr + end + + def test_count + assert_equal 3, @n.nzcount + end + +end + +class RubySparse::CSC::CreationTest < Minitest::Test + + def setup + @n = RubySparse::CSC.new [3, 3], [1, 2, 3], [0, 1, 2], [0, 1, 2] + @m = RubySparse::CSC.new [4, 6], [1, 2, 3, 4, 5, 6, 7, 8], [0, 0, 1, 1, 2, 2, 2, 3], [0, 1, 1, 3, 2, 3, 4, 5] + end + + def test_ndims + assert_equal 2, @n.dim + end + + def test_shape + assert_equal [3, 3], @n.shape + end + + def test_elements + assert_equal [1, 2, 3], @n.elements + end + + def test_indices + assert_equal [0, 0, 1, 2, 1, 2, 2, 3], @m.indices + end + + def test_indptr + assert_equal [0, 1, 3, 4, 6, 7, 8], @m.indptr + end + + def test_count + assert_equal 3, @n.nzcount + end + end \ No newline at end of file