text:
#include <php.h>
#include <Zend/zend_interfaces.h>
#include <Zend/zend_exceptions.h>
#include <ext/spl/spl_iterators.h>
#include <ext/spl/spl_exceptions.h>
#include <cblas.h>
#include <Interop/Polite/Math/Matrix.h>

#ifdef HAVE_CONFIG_H
#include "config.h"
#endif

#include "php_rindow_openblas.h"

static zend_object_handlers rindow_openblas_blas_object_handlers;

// destractor
static void php_rindow_openblas_blas_free_object(zend_object* object)
{
    zend_object_std_dtor(object);
}

// constructor
static zend_object* php_rindow_openblas_blas_create_object(zend_class_entry* class_type) /* {{{ */
{
    zend_object* intern = NULL;

    intern = (zend_object*)ecalloc(1, sizeof(zend_object) + zend_object_properties_size(class_type));

    zend_object_std_init(intern, class_type);
    object_properties_init(intern, class_type);

    intern->handlers = &rindow_openblas_blas_object_handlers;

    return intern;
} /* }}} */

/* Method Rindow\OpenBLAS\Blas::
    public function getNumThreads() : int
 {{{ */
static PHP_METHOD(Blas, getNumThreads)
{
    int n;
    n = openblas_get_num_threads();
    RETURN_LONG(n);
}
/* }}} */

/* Method Rindow\OpenBLAS\Blas::
    public function getNumProcs() : int
 {{{ */
static PHP_METHOD(Blas, getNumProcs)
{
    int n;
    n = openblas_get_num_procs();
    RETURN_LONG(n);
}
/* }}} */

/* Method Rindow\OpenBLAS\Blas::
    public function getConfig() : string
 {{{ */
static PHP_METHOD(Blas, getConfig)
{
    char *s;
    s = openblas_get_config();
    RETURN_STRING(s);
}
/* }}} */

/* Method Rindow\OpenBLAS\Blas::
    public function getCorename() : string
 {{{ */
static PHP_METHOD(Blas, getCorename)
{
    char *s;
    s = openblas_get_corename();
    RETURN_STRING(s);
}
/* }}} */

/* Method Rindow\OpenBLAS\Blas::
    public function scal(
        int $n,
        float $alpha,
        Buffer $X, int $offsetX, int $incX) : void
 {{{ */
static PHP_METHOD(Blas, scal)
{
    php_interop_polite_math_matrix_linear_buffer_t* buffer;
    zend_long n;
    double alpha;
    zval* x=NULL;
    zend_long offsetX;
    zend_long incX;

    ZEND_PARSE_PARAMETERS_START_EX(ZEND_PARSE_PARAMS_THROW, 5, 5)
        Z_PARAM_LONG(n)
        Z_PARAM_DOUBLE(alpha)
        Z_PARAM_OBJECT(x) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetX)
        Z_PARAM_LONG(incX)
    ZEND_PARSE_PARAMETERS_END();

    buffer = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(x);
    if(php_rindow_openblas_assert_buffer_type(buffer,"x")) {
        return;
    }
    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_N, n)) {
        return;
    }
    if(php_rindow_openblas_assert_vector_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_X, buffer,n,offsetX,incX)) {
        return;
    }
    switch (buffer->dtype) {
        case php_interop_polite_math_matrix_dtype_float32:
            cblas_sscal((blasint)n, (float)alpha, &(((float *)buffer->data)[offsetX]), (blasint)incX);
            break;
        case php_interop_polite_math_matrix_dtype_float64:
            cblas_dscal((blasint)n, (double)alpha, &(((double *)buffer->data)[offsetX]), (blasint)incX);
            break;
        default:
            zend_throw_exception(spl_ce_RuntimeException, "Unsupported data type.", 0);
            return;
    }
}
/* }}} */

/* Method Rindow\OpenBLAS\Blas::
    public function axpy(
        int $n,
        float $alpha,
        Buffer $X, int $offsetX, int $incX,
        Buffer $Y, int $offsetY, int $incY ) : void
 {{{ */
static PHP_METHOD(Blas, axpy)
{
    php_interop_polite_math_matrix_linear_buffer_t* bufferX;
    php_interop_polite_math_matrix_linear_buffer_t* bufferY;
    zend_long n;
    double alpha;
    zval* x=NULL;
    zend_long offsetX;
    zend_long incX;
    zval* y=NULL;
    zend_long offsetY;
    zend_long incY;

    ZEND_PARSE_PARAMETERS_START_EX(ZEND_PARSE_PARAMS_THROW, 8, 8)
        Z_PARAM_LONG(n)
        Z_PARAM_DOUBLE(alpha)
        Z_PARAM_OBJECT(x) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetX)
        Z_PARAM_LONG(incX)
        Z_PARAM_OBJECT(y) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetY)
        Z_PARAM_LONG(incY)
    ZEND_PARSE_PARAMETERS_END();

    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_N, n)) {
        return;
    }
    // Check Buffer X
    bufferX = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(x);
    if(php_rindow_openblas_assert_buffer_type(bufferX,"x")) {
        return;
    }
    if(php_rindow_openblas_assert_vector_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_X, bufferX,n,offsetX,incX)) {
        return;
    }

    // Check Buffer Y
    bufferY = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(y);
    if(php_rindow_openblas_assert_buffer_type(bufferY,"y")) {
        return;
    }
    if(php_rindow_openblas_assert_vector_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_Y, bufferY,n,offsetY,incY)) {
        return;
    }

    // Check Buffer X and Y
    if(bufferX->dtype!=bufferY->dtype) {
        zend_throw_exception(spl_ce_InvalidArgumentException, "Unmatch data type for X and Y", 0);
        return;
    }

    switch (bufferX->dtype) {
        case php_interop_polite_math_matrix_dtype_float32:
            cblas_saxpy((blasint)n, (float)alpha,
                &(((float *)bufferX->data)[offsetX]), (blasint)incX,
                &(((float *)bufferY->data)[offsetY]), (blasint)incY);
            break;
        case php_interop_polite_math_matrix_dtype_float64:
            cblas_daxpy((blasint)n, (double)alpha,
                &(((double *)bufferX->data)[offsetX]), (blasint)incX,
                &(((double *)bufferY->data)[offsetY]), (blasint)incY);
            break;
        default:
            zend_throw_exception(spl_ce_RuntimeException, "Unsupported data type.", 0);
            return;
    }
}
/* }}} */

/* Method Rindow\OpenBLAS\Blas::
    public function dot(
        int $n,
        Buffer $X, int $offsetX, int $incX,
        Buffer $Y, int $offsetY, int $incY ) : float
 {{{ */
static PHP_METHOD(Blas, dot)
{
    php_interop_polite_math_matrix_linear_buffer_t* bufferX;
    php_interop_polite_math_matrix_linear_buffer_t* bufferY;
    zend_long n;
    zval* x=NULL;
    zend_long offsetX;
    zend_long incX;
    zval* y=NULL;
    zend_long offsetY;
    zend_long incY;
    double result;

    ZEND_PARSE_PARAMETERS_START_EX(ZEND_PARSE_PARAMS_THROW, 7, 7)
        Z_PARAM_LONG(n)
        Z_PARAM_OBJECT(x) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetX)
        Z_PARAM_LONG(incX)
        Z_PARAM_OBJECT(y) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetY)
        Z_PARAM_LONG(incY)
    ZEND_PARSE_PARAMETERS_END();

    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_N, n)) {
        return;
    }
    // Check Buffer X
    bufferX = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(x);
    if(php_rindow_openblas_assert_buffer_type(bufferX,"x")) {
        return;
    }
    if(php_rindow_openblas_assert_vector_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_X, bufferX,n,offsetX,incX)) {
        return;
    }

    // Check Buffer Y
    bufferY = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(y);
    if(php_rindow_openblas_assert_buffer_type(bufferY,"y")) {
        return;
    }
    if(php_rindow_openblas_assert_vector_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_Y, bufferY,n,offsetY,incY)) {
        return;
    }

    // Check Buffer X and Y
    if(bufferX->dtype!=bufferY->dtype) {
        zend_throw_exception(spl_ce_InvalidArgumentException, "Unmatch data type for X and Y", 0);
        return;
    }

    switch (bufferX->dtype) {
        case php_interop_polite_math_matrix_dtype_float32:
            result = (double)cblas_sdot((blasint)n,
                &(((float *)bufferX->data)[offsetX]), (blasint)incX,
                &(((float *)bufferY->data)[offsetY]), (blasint)incY);
            break;
        case php_interop_polite_math_matrix_dtype_float64:
            result = (double)cblas_ddot((blasint)n,
                &(((double *)bufferX->data)[offsetX]), (blasint)incX,
                &(((double *)bufferY->data)[offsetY]), (blasint)incY);
            break;
        default:
            zend_throw_exception(spl_ce_RuntimeException, "Unsupported data type.", 0);
            return;
    }
    RETURN_DOUBLE(result);
}
/* }}} */

/* Method Rindow\OpenBLAS\Blas::
    public function asum(
        int $n,
        Buffer $X, int $offsetX, int $incX ) : float
 {{{ */
static PHP_METHOD(Blas, asum)
{
    php_interop_polite_math_matrix_linear_buffer_t* buffer;
    zend_long n;
    zval* x=NULL;
    zend_long offsetX;
    zend_long incX;
    double result;

    ZEND_PARSE_PARAMETERS_START_EX(ZEND_PARSE_PARAMS_THROW, 4, 4)
        Z_PARAM_LONG(n)
        Z_PARAM_OBJECT(x) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetX)
        Z_PARAM_LONG(incX)
    ZEND_PARSE_PARAMETERS_END();

    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_N, n)) {
        return;
    }
    buffer = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(x);
    if(php_rindow_openblas_assert_buffer_type(buffer,"x")) {
        return;
    }
    if(php_rindow_openblas_assert_vector_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_X, buffer,n,offsetX,incX)) {
        return;
    }
    switch (buffer->dtype) {
        case php_interop_polite_math_matrix_dtype_float32:
            result = (double)cblas_sasum((blasint)n, &(((float *)buffer->data)[offsetX]), (blasint)incX);
            break;
        case php_interop_polite_math_matrix_dtype_float64:
            result = (double)cblas_dasum((blasint)n, &(((double *)buffer->data)[offsetX]), (blasint)incX);
            break;
        default:
            zend_throw_exception(spl_ce_RuntimeException, "Unsupported data type.", 0);
            return;
    }
    RETURN_DOUBLE(result);
}
/* }}} */

/* Method Rindow\OpenBLAS\Blas::
    public function iamax(
        int $n,
        Buffer $X, int $offsetX, int $incX ) : int
 {{{ */
static PHP_METHOD(Blas, iamax)
{
    php_interop_polite_math_matrix_linear_buffer_t* buffer;
    zend_long n;
    zval* x=NULL;
    zend_long offsetX;
    zend_long incX;
    zend_long resultIdx;

    ZEND_PARSE_PARAMETERS_START_EX(ZEND_PARSE_PARAMS_THROW, 4, 4)
        Z_PARAM_LONG(n)
        Z_PARAM_OBJECT(x) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetX)
        Z_PARAM_LONG(incX)
    ZEND_PARSE_PARAMETERS_END();

    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_N, n)) {
        return;
    }
    buffer = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(x);
    if(php_rindow_openblas_assert_buffer_type(buffer,"x")) {
        return;
    }
    if(php_rindow_openblas_assert_vector_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_X, buffer,n,offsetX,incX)) {
        return;
    }
    switch (buffer->dtype) {
        case php_interop_polite_math_matrix_dtype_float32:
            resultIdx = (zend_long)cblas_isamax((blasint)n, &(((float *)buffer->data)[offsetX]), (blasint)incX);
            break;
        case php_interop_polite_math_matrix_dtype_float64:
            resultIdx = (zend_long)cblas_idamax((blasint)n, &(((double *)buffer->data)[offsetX]), (blasint)incX);
            break;
        default:
            zend_throw_exception(spl_ce_RuntimeException, "Unsupported data type.", 0);
            return;
    }
    RETURN_LONG(resultIdx);
}
/* }}} */

#ifdef OPENBLAS_HAVE_IAMIN
/* Method Rindow\OpenBLAS\Blas::
    public function iamin(
        int $n,
        Buffer $X, int $offsetX, int $incX ) : int
 {{{ */
static PHP_METHOD(Blas, iamin)
{
    php_interop_polite_math_matrix_linear_buffer_t* buffer;
    zend_long n;
    zval* x=NULL;
    zend_long offsetX;
    zend_long incX;
    zend_long resultIdx;

    ZEND_PARSE_PARAMETERS_START_EX(ZEND_PARSE_PARAMS_THROW, 4, 4)
        Z_PARAM_LONG(n)
        Z_PARAM_OBJECT(x) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetX)
        Z_PARAM_LONG(incX)
    ZEND_PARSE_PARAMETERS_END();

    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_N, n)) {
        return;
    }
    buffer = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(x);
    if(php_rindow_openblas_assert_buffer_type(buffer,"x")) {
        return;
    }
    if(php_rindow_openblas_assert_vector_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_X, buffer,n,offsetX,incX)) {
        return;
    }
    switch (buffer->dtype) {
        case php_interop_polite_math_matrix_dtype_float32:
            resultIdx = (zend_long)cblas_isamin((blasint)n, &(((float *)buffer->data)[offsetX]), (blasint)incX);
            break;
        case php_interop_polite_math_matrix_dtype_float64:
            resultIdx = (zend_long)cblas_idamin((blasint)n, &(((double *)buffer->data)[offsetX]), (blasint)incX);
            break;
        default:
            zend_throw_exception(spl_ce_RuntimeException, "Unsupported data type.", 0);
            return;
    }
    RETURN_LONG(resultIdx);
}
/* }}} */
#endif

/* Method Rindow\OpenBLAS\Blas::
    public function copy(
        int $n,
        Buffer $X, int $offsetX, int $incX,
        Buffer $Y, int $offsetY, int $incY ) : void
 {{{ */
static PHP_METHOD(Blas, copy)
{
    php_interop_polite_math_matrix_linear_buffer_t* bufferX;
    php_interop_polite_math_matrix_linear_buffer_t* bufferY;
    zend_long n;
    zval* x=NULL;
    zend_long offsetX;
    zend_long incX;
    zval* y=NULL;
    zend_long offsetY;
    zend_long incY;

    ZEND_PARSE_PARAMETERS_START_EX(ZEND_PARSE_PARAMS_THROW, 7, 7)
        Z_PARAM_LONG(n)
        Z_PARAM_OBJECT(x) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetX)
        Z_PARAM_LONG(incX)
        Z_PARAM_OBJECT(y) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetY)
        Z_PARAM_LONG(incY)
    ZEND_PARSE_PARAMETERS_END();

    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_N, n)) {
        return;
    }
    // Check Buffer X
    bufferX = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(x);
    if(php_rindow_openblas_assert_buffer_type(bufferX,"x")) {
        return;
    }
    if(php_rindow_openblas_assert_vector_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_X, bufferX,n,offsetX,incX)) {
        return;
    }

    // Check Buffer Y
    bufferY = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(y);
    if(php_rindow_openblas_assert_buffer_type(bufferY,"y")) {
        return;
    }
    if(php_rindow_openblas_assert_vector_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_Y, bufferY,n,offsetY,incY)) {
        return;
    }

    // Check Buffer X and Y
    if(bufferX->dtype!=bufferY->dtype) {
        zend_throw_exception(spl_ce_InvalidArgumentException, "Unmatch data type for X and Y", 0);
        return;
    }

    switch (bufferX->dtype) {
        case php_interop_polite_math_matrix_dtype_float32:
            cblas_scopy((blasint)n,
                &(((float *)bufferX->data)[offsetX]), (blasint)incX,
                &(((float *)bufferY->data)[offsetY]), (blasint)incY);
            break;
        case php_interop_polite_math_matrix_dtype_float64:
            cblas_dcopy((blasint)n,
                &(((double *)bufferX->data)[offsetX]), (blasint)incX,
                &(((double *)bufferY->data)[offsetY]), (blasint)incY);
            break;
        default:
            {
                zend_long i,idX,idY;
                int valueSize;
                uint8_t *x,*y;
                valueSize = php_rindow_openblas_common_dtype_to_valuesize(bufferX->dtype);
                x = php_rindow_openblas_get_address(bufferX,offsetX,valueSize);
                y = php_rindow_openblas_get_address(bufferY,offsetY,valueSize);
                if(incX==1 && incY==1) {
                    memcpy(y,x,valueSize*n);
                } else {
                    for(i=0,idX=0,idY=0; i<n; i++,idX+=incX,idY+=incY) {
                        memcpy(&y[valueSize*idY],&x[valueSize*idX],valueSize);
                    }
                }
            }
            break;
    }
}
/* }}} */

/* Method Rindow\OpenBLAS\Blas::
    public function nrm2(
        int $n,
        Buffer $X, int $offsetX, int $incX ) : float
 {{{ */
static PHP_METHOD(Blas, nrm2)
{
    php_interop_polite_math_matrix_linear_buffer_t* buffer;
    zend_long n;
    zval* x=NULL;
    zend_long offsetX;
    zend_long incX;
    double result;

    ZEND_PARSE_PARAMETERS_START_EX(ZEND_PARSE_PARAMS_THROW, 4, 4)
        Z_PARAM_LONG(n)
        Z_PARAM_OBJECT(x) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetX)
        Z_PARAM_LONG(incX)
    ZEND_PARSE_PARAMETERS_END();

    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_N, n)) {
        return;
    }
    buffer = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(x);
    if(php_rindow_openblas_assert_buffer_type(buffer,"x")) {
        return;
    }
    if(php_rindow_openblas_assert_vector_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_X, buffer,n,offsetX,incX)) {
        return;
    }
    switch (buffer->dtype) {
        case php_interop_polite_math_matrix_dtype_float32:
            result = (double)cblas_snrm2((blasint)n, &(((float *)buffer->data)[offsetX]), (blasint)incX);
            break;
        case php_interop_polite_math_matrix_dtype_float64:
            result = (double)cblas_dnrm2((blasint)n, &(((double *)buffer->data)[offsetX]), (blasint)incX);
            break;
        default:
            zend_throw_exception(spl_ce_RuntimeException, "Unsupported data type.", 0);
            return;
    }
    RETURN_DOUBLE(result);
}
/* }}} */

/* Method Rindow\OpenBLAS\Blas::
    public function rotg(
        Buffer $A, int $offsetA,
        Buffer $B, int $offsetB,
        Buffer $C, int $offsetC,
        Buffer $S, int $offsetS,
    ) : void
 {{{ */
static PHP_METHOD(Blas, rotg)
{
    zval* a;
    zend_long offsetA;
    zval* b;
    zend_long offsetB;
    zval* c;
    zend_long offsetC;
    zval* s;
    zend_long offsetS;
    php_interop_polite_math_matrix_linear_buffer_t* bufferA;
    php_interop_polite_math_matrix_linear_buffer_t* bufferB;
    php_interop_polite_math_matrix_linear_buffer_t* bufferC;
    php_interop_polite_math_matrix_linear_buffer_t* bufferS;

    ZEND_PARSE_PARAMETERS_START_EX(ZEND_PARSE_PARAMS_THROW, 8, 8)
        Z_PARAM_OBJECT(a) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetA)
        Z_PARAM_OBJECT(b) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetB)
        Z_PARAM_OBJECT(c) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetC)
        Z_PARAM_OBJECT(s) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetS)
    ZEND_PARSE_PARAMETERS_END();

    // Check Buffer A
    bufferA = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(a);
    if(php_rindow_openblas_assert_buffer_type(bufferA,"a")) {
        return;
    }
    if(php_rindow_openblas_assert_vector_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_X, bufferA,1,offsetA,1)) {
        return;
    }

    // Check Buffer B
    bufferB = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(b);
    if(php_rindow_openblas_assert_buffer_type(bufferB,"b")) {
        return;
    }
    if(php_rindow_openblas_assert_vector_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_Y, bufferB,1,offsetB,1)) {
        return;
    }

    // Check Buffer C
    bufferC = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(c);
    if(php_rindow_openblas_assert_buffer_type(bufferC,"c")) {
        return;
    }
    if(php_rindow_openblas_assert_vector_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_X, bufferC,1,offsetC,1)) {
        return;
    }

    // Check Buffer S
    bufferS = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(s);
    if(php_rindow_openblas_assert_buffer_type(bufferS,"s")) {
        return;
    }
    if(php_rindow_openblas_assert_vector_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_Y, bufferS,1,offsetS,1)) {
        return;
    }

    // Check Buffer A and B and C and S
    if(bufferA->dtype!=bufferB->dtype || bufferB->dtype!=bufferC->dtype ||
        bufferC->dtype!=bufferS->dtype || bufferS->dtype!=bufferA->dtype) {
        zend_throw_exception(spl_ce_InvalidArgumentException, "Unmatch data type for A,B,C and S", 0);
        return;
    }

    switch (bufferA->dtype) {
        case php_interop_polite_math_matrix_dtype_float32:
            cblas_srotg(
                &(((float *)bufferA->data)[offsetA]),
                &(((float *)bufferB->data)[offsetB]),
                &(((float *)bufferC->data)[offsetC]),
                &(((float *)bufferS->data)[offsetS])
            );
            break;
        case php_interop_polite_math_matrix_dtype_float64:
            cblas_drotg(
                &(((double *)bufferA->data)[offsetA]),
                &(((double *)bufferB->data)[offsetB]),
                &(((double *)bufferC->data)[offsetC]),
                &(((double *)bufferS->data)[offsetS])
            );
            break;
        default:
            zend_throw_exception(spl_ce_RuntimeException, "Unsupported data type.", 0);
            return;
    }
}
/* }}} */

/* Method Rindow\OpenBLAS\Blas::
    public function rot(
        int $n,
        Buffer $X, int $offsetX, int $incX,
        Buffer $Y, int $offsetY, int $incY,
        Buffer $C, int $offsetC,
        Buffer $S, int $offsetS,
    ) : void
 {{{ */
static PHP_METHOD(Blas, rot)
{
    zend_long n;
    zval* x;
    zend_long offsetX;
    zend_long incX;
    zval* y;
    zend_long offsetY;
    zend_long incY;
    zval* c;
    zend_long offsetC;
    zval* s;
    zend_long offsetS;
    php_interop_polite_math_matrix_linear_buffer_t* bufferX;
    php_interop_polite_math_matrix_linear_buffer_t* bufferY;
    php_interop_polite_math_matrix_linear_buffer_t* bufferC;
    php_interop_polite_math_matrix_linear_buffer_t* bufferS;

    ZEND_PARSE_PARAMETERS_START_EX(ZEND_PARSE_PARAMS_THROW, 11, 11)
        Z_PARAM_LONG(n)
        Z_PARAM_OBJECT(x) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetX)
        Z_PARAM_LONG(incX)
        Z_PARAM_OBJECT(y) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetY)
        Z_PARAM_LONG(incY)
        Z_PARAM_OBJECT(c) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetC)
        Z_PARAM_OBJECT(s) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetS)
    ZEND_PARSE_PARAMETERS_END();

    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_N, n)) {
        return;
    }
    // Check Buffer X
    bufferX = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(x);
    if(php_rindow_openblas_assert_buffer_type(bufferX,"x")) {
        return;
    }
    if(php_rindow_openblas_assert_vector_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_X, bufferX,n,offsetX,incX)) {
        return;
    }

    // Check Buffer Y
    bufferY = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(y);
    if(php_rindow_openblas_assert_buffer_type(bufferY,"y")) {
        return;
    }
    if(php_rindow_openblas_assert_vector_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_Y, bufferY,n,offsetY,incY)) {
        return;
    }

    // Check Buffer C
    bufferC = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(c);
    if(php_rindow_openblas_assert_buffer_type(bufferC,"c")) {
        return;
    }
    if(php_rindow_openblas_assert_vector_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_X, bufferC,1,offsetC,1)) {
        return;
    }

    // Check Buffer S
    bufferS = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(s);
    if(php_rindow_openblas_assert_buffer_type(bufferS,"s")) {
        return;
    }
    if(php_rindow_openblas_assert_vector_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_Y, bufferS,1,offsetS,1)) {
        return;
    }

    // Check Buffer X and Y
    if(bufferX->dtype!=bufferY->dtype||bufferY->dtype!=bufferC->dtype||
        bufferC->dtype!=bufferS->dtype) {
        zend_throw_exception(spl_ce_InvalidArgumentException, "Unmatch data type for X,Y,C and S", 0);
        return;
    }

    switch (bufferX->dtype) {
        case php_interop_polite_math_matrix_dtype_float32:
            cblas_srot((blasint)n,
                &(((float *)bufferX->data)[offsetX]), (blasint)incX,
                &(((float *)bufferY->data)[offsetY]), (blasint)incY,
                ((float *)bufferC->data)[offsetC],
                ((float *)bufferS->data)[offsetS]
            );
            break;
        case php_interop_polite_math_matrix_dtype_float64:
            cblas_drot((blasint)n,
                &(((double *)bufferX->data)[offsetX]), (blasint)incX,
                &(((double *)bufferY->data)[offsetY]), (blasint)incY,
                ((double *)bufferC->data)[offsetC],
                ((double *)bufferS->data)[offsetS]
            );
            break;
        default:
            zend_throw_exception(spl_ce_RuntimeException, "Unsupported data type.", 0);
            return;
    }
}
/* }}} */

/* Method Rindow\OpenBLAS\Blas::
    public function swap(
        int $n,
        Buffer $X, int $offsetX, int $incX,
        Buffer $Y, int $offsetY, int $incY ) : void
 {{{ */
static PHP_METHOD(Blas, swap)
{
    php_interop_polite_math_matrix_linear_buffer_t* bufferX;
    php_interop_polite_math_matrix_linear_buffer_t* bufferY;
    zend_long n;
    zval* x=NULL;
    zend_long offsetX;
    zend_long incX;
    zval* y=NULL;
    zend_long offsetY;
    zend_long incY;

    ZEND_PARSE_PARAMETERS_START_EX(ZEND_PARSE_PARAMS_THROW, 7, 7)
        Z_PARAM_LONG(n)
        Z_PARAM_OBJECT(x) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetX)
        Z_PARAM_LONG(incX)
        Z_PARAM_OBJECT(y) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetY)
        Z_PARAM_LONG(incY)
    ZEND_PARSE_PARAMETERS_END();

    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_N, n)) {
        return;
    }
    // Check Buffer X
    bufferX = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(x);
    if(php_rindow_openblas_assert_buffer_type(bufferX,"x")) {
        return;
    }
    if(php_rindow_openblas_assert_vector_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_X, bufferX,n,offsetX,incX)) {
        return;
    }

    // Check Buffer Y
    bufferY = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(y);
    if(php_rindow_openblas_assert_buffer_type(bufferY,"y")) {
        return;
    }
    if(php_rindow_openblas_assert_vector_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_Y, bufferY,n,offsetY,incY)) {
        return;
    }

    // Check Buffer X and Y
    if(bufferX->dtype!=bufferY->dtype) {
        zend_throw_exception(spl_ce_InvalidArgumentException, "Unmatch data type for X and Y", 0);
        return;
    }

    switch (bufferX->dtype) {
        case php_interop_polite_math_matrix_dtype_float32:
            cblas_sswap((blasint)n,
                &(((float *)bufferX->data)[offsetX]), (blasint)incX,
                &(((float *)bufferY->data)[offsetY]), (blasint)incY);
            break;
        case php_interop_polite_math_matrix_dtype_float64:
            cblas_dswap((blasint)n,
                &(((double *)bufferX->data)[offsetX]), (blasint)incX,
                &(((double *)bufferY->data)[offsetY]), (blasint)incY);
            break;
        default:
            zend_throw_exception(spl_ce_RuntimeException, "Unsupported data type.", 0);
            return;
    }
}
/* }}} */

/* Method Rindow\OpenBLAS\Blas::
    public function gemv(
        int $order,
        int $trans,
        int $m,
        int $n,
        float $alpha,
        Buffer $A, int $offsetA, int $ldA,
        Buffer $X, int $offsetX, int $incX,
        float $beta,
        Buffer $Y, int $offsetY, int $incY ) : void
 {{{ */
static PHP_METHOD(Blas, gemv)
{
    php_interop_polite_math_matrix_linear_buffer_t* bufferA;
    php_interop_polite_math_matrix_linear_buffer_t* bufferX;
    php_interop_polite_math_matrix_linear_buffer_t* bufferY;
    zend_long order;
    zend_long trans;
    zend_long m;
    zend_long n;
    double alpha;
    zval* a=NULL;
    zend_long offsetA;
    zend_long ldA;
    zval* x=NULL;
    zend_long offsetX;
    zend_long incX;
    double beta;
    zval* y=NULL;
    zend_long offsetY;
    zend_long incY;

    ZEND_PARSE_PARAMETERS_START_EX(ZEND_PARSE_PARAMS_THROW, 15, 15)
        Z_PARAM_LONG(order)
        Z_PARAM_LONG(trans)
        Z_PARAM_LONG(m)
        Z_PARAM_LONG(n)
        Z_PARAM_DOUBLE(alpha)
        Z_PARAM_OBJECT(a) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetA)
        Z_PARAM_LONG(ldA)
        Z_PARAM_OBJECT(x) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetX)
        Z_PARAM_LONG(incX)
        Z_PARAM_DOUBLE(beta)
        Z_PARAM_OBJECT(y) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetY)
        Z_PARAM_LONG(incY)
    ZEND_PARSE_PARAMETERS_END();

    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_M, m)) {
        return;
    }
    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_N, n)) {
        return;
    }
    // Check Buffer A
    bufferA = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(a);
    if(php_rindow_openblas_assert_buffer_type(bufferA,"a")) {
        return;
    }
    if(php_rindow_openblas_assert_matrix_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_A, bufferA,m,n,offsetA,ldA)) {
        return;
    }

    // Check Buffer X
    bufferX = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(x);
    if(php_rindow_openblas_assert_buffer_type(bufferX,"x")) {
        return;
    }

    // Check Buffer Y
    bufferY = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(y);
    if(php_rindow_openblas_assert_buffer_type(bufferY,"y")) {
        return;
    }

    // Check Buffer size X and Y
    {
        zend_long rows,cols;
        if(trans==CblasNoTrans || trans==CblasConjNoTrans ) {
            rows = m; cols = n;
        } else if(trans==CblasTrans || trans==CblasConjTrans) {
            rows = n; cols = m;
        } else {
            zend_throw_exception(spl_ce_RuntimeException, "unknown transpose mode for bufferA.", 0);
            return;
        }
        if(php_rindow_openblas_assert_vector_buffer_spec(
            PHP_RINDOW_OPENBLAS_ASSERT_X, bufferX,cols,offsetX,incX)) {
            return;
        }
        if(php_rindow_openblas_assert_vector_buffer_spec(
            PHP_RINDOW_OPENBLAS_ASSERT_Y, bufferY,rows,offsetY,incY)) {
            return;
        }
    }

    // Check Buffer A and X and Y
    if(bufferA->dtype!=bufferX->dtype || bufferX->dtype!=bufferY->dtype) {
        zend_throw_exception(spl_ce_InvalidArgumentException, "Unmatch data type for A and X and Y", 0);
        return;
    }

    switch (bufferX->dtype) {
        case php_interop_polite_math_matrix_dtype_float32:
            cblas_sgemv(
                (OPENBLAS_CONST enum CBLAS_ORDER)order,
                (OPENBLAS_CONST enum CBLAS_TRANSPOSE)trans,
                (blasint)m,(blasint)n,
                (float)alpha,
                &(((float *)bufferA->data)[offsetA]), (blasint)ldA,
                &(((float *)bufferX->data)[offsetX]), (blasint)incX,
                (float)beta,
                &(((float *)bufferY->data)[offsetY]), (blasint)incY);
            break;
        case php_interop_polite_math_matrix_dtype_float64:
            cblas_dgemv(
                (OPENBLAS_CONST enum CBLAS_ORDER)order,
                (OPENBLAS_CONST enum CBLAS_TRANSPOSE)trans,
                (blasint)m,(blasint)n,
                (double)alpha,
                &(((double *)bufferA->data)[offsetA]), (blasint)ldA,
                &(((double *)bufferX->data)[offsetX]), (blasint)incX,
                (double)beta,
                &(((double *)bufferY->data)[offsetY]), (blasint)incY);
            break;
        default:
            zend_throw_exception(spl_ce_RuntimeException, "Unsupported data type.", 0);
            return;
    }
}
/* }}} */

/* Method Rindow\OpenBLAS\Blas::
    public function gemm(
        int $order,
        int $transA,
        int $transB,
        int $m,
        int $n,
        int $k,
        float $alpha,
        Buffer $A, int $offsetA, int $ldA,
        Buffer $B, int $offsetB, int $ldB,
        float $beta,
        Buffer $C, int $offsetC, int $ldC ) : void
 {{{ */
static PHP_METHOD(Blas, gemm)
{
    php_interop_polite_math_matrix_linear_buffer_t* bufferA;
    php_interop_polite_math_matrix_linear_buffer_t* bufferB;
    php_interop_polite_math_matrix_linear_buffer_t* bufferC;
    zend_long order;
    zend_long transA;
    zend_long transB;
    zend_long m;
    zend_long n;
    zend_long k;
    double alpha;
    zval* a=NULL;
    zend_long offsetA;
    zend_long ldA;
    zval* b=NULL;
    zend_long offsetB;
    zend_long ldB;
    double beta;
    zval* c=NULL;
    zend_long offsetC;
    zend_long ldC;

    ZEND_PARSE_PARAMETERS_START_EX(ZEND_PARSE_PARAMS_THROW, 17, 17)
        Z_PARAM_LONG(order)
        Z_PARAM_LONG(transA)
        Z_PARAM_LONG(transB)
        Z_PARAM_LONG(m)
        Z_PARAM_LONG(n)
        Z_PARAM_LONG(k)
        Z_PARAM_DOUBLE(alpha)
        Z_PARAM_OBJECT(a) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetA)
        Z_PARAM_LONG(ldA)
        Z_PARAM_OBJECT(b) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetB)
        Z_PARAM_LONG(ldB)
        Z_PARAM_DOUBLE(beta)
        Z_PARAM_OBJECT(c) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetC)
        Z_PARAM_LONG(ldC)
    ZEND_PARSE_PARAMETERS_END();

    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_M, m)) {
        return;
    }
    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_N, n)) {
        return;
    }
    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_K, k)) {
        return;
    }
    // Check Buffer A
    bufferA = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(a);
    if(php_rindow_openblas_assert_buffer_type(bufferA,"a")) {
        return;
    }
    {
        zend_long rows,cols;
        if(transA==CblasNoTrans || transA==CblasConjNoTrans) {
            rows = m; cols = k;
        } else if(transA==CblasTrans || transA==CblasConjTrans) {
            rows = k; cols = m;
        } else {
            zend_throw_exception(spl_ce_RuntimeException, "unknown transpose mode for bufferA.", 0);
            return;
        }
        if(php_rindow_openblas_assert_matrix_buffer_spec(
            PHP_RINDOW_OPENBLAS_ASSERT_A, bufferA,rows,cols,offsetA,ldA)) {
            return;
        }
    }

    // Check Buffer B
    bufferB = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(b);
    if(php_rindow_openblas_assert_buffer_type(bufferB,"b")) {
        return;
    }
    {
        zend_long rows,cols;
        if(transB==CblasNoTrans || transB==CblasConjNoTrans) {
            rows = k; cols = n;
        } else if(transB==CblasTrans || transB==CblasConjTrans) {
            rows = n; cols = k;
        } else {
            zend_throw_exception(spl_ce_RuntimeException, "unknown transpose mode for bufferB.", 0);
            return;
        }
        if(php_rindow_openblas_assert_matrix_buffer_spec(
            PHP_RINDOW_OPENBLAS_ASSERT_B, bufferB,rows,cols,offsetB,ldB)) {
            return;
        }
    }

    // Check Buffer C
    bufferC = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(c);
    if(php_rindow_openblas_assert_buffer_type(bufferC,"c")) {
        return;
    }
    if(php_rindow_openblas_assert_matrix_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_C, bufferC,m,n,offsetC,ldC)) {
        return;
    }

    // Check Buffer A and B and C
    if(bufferA->dtype!=bufferB->dtype || bufferB->dtype!=bufferC->dtype) {
        zend_throw_exception(spl_ce_InvalidArgumentException, "Unmatch data type for A and B and C", 0);
        return;
    }

    switch (bufferA->dtype) {
        case php_interop_polite_math_matrix_dtype_float32:
            cblas_sgemm(
                (OPENBLAS_CONST enum CBLAS_ORDER)order,
                (OPENBLAS_CONST enum CBLAS_TRANSPOSE)transA,
                (OPENBLAS_CONST enum CBLAS_TRANSPOSE)transB,
                (blasint)m,(blasint)n,(blasint)k,
                (float)alpha,
                &(((float *)bufferA->data)[offsetA]), (blasint)ldA,
                &(((float *)bufferB->data)[offsetB]), (blasint)ldB,
                (float)beta,
                &(((float *)bufferC->data)[offsetC]), (blasint)ldC);
            break;
        case php_interop_polite_math_matrix_dtype_float64:
            cblas_dgemm(
                (OPENBLAS_CONST enum CBLAS_ORDER)order,
                (OPENBLAS_CONST enum CBLAS_TRANSPOSE)transA,
                (OPENBLAS_CONST enum CBLAS_TRANSPOSE)transB,
                (blasint)m,(blasint)n,(blasint)k,
                (double)alpha,
                &(((double *)bufferA->data)[offsetA]), (blasint)ldA,
                &(((double *)bufferB->data)[offsetB]), (blasint)ldB,
                (double)beta,
                &(((double *)bufferC->data)[offsetC]), (blasint)ldC);
            break;
        default:
            zend_throw_exception(spl_ce_RuntimeException, "Unsupported data type.", 0);
            return;
    }
}
/* }}} */

/* Method Rindow\OpenBLAS\Blas::
    public function symm(
        int $order,
        int $side,
        int $uplo,
        int $m,
        int $n,
        float $alpha,
        Buffer $A, int $offsetA, int $ldA,
        Buffer $B, int $offsetB, int $ldB,
        float $beta,
        Buffer $C, int $offsetC, int $ldC ) : void
 {{{ */
static PHP_METHOD(Blas, symm)
{
    php_interop_polite_math_matrix_linear_buffer_t* bufferA;
    php_interop_polite_math_matrix_linear_buffer_t* bufferB;
    php_interop_polite_math_matrix_linear_buffer_t* bufferC;
    zend_long order;
    zend_long side;
    zend_long uplo;
    zend_long m;
    zend_long n;
    double alpha;
    zval* a=NULL;
    zend_long offsetA;
    zend_long ldA;
    zval* b=NULL;
    zend_long offsetB;
    zend_long ldB;
    double beta;
    zval* c=NULL;
    zend_long offsetC;
    zend_long ldC;

    ZEND_PARSE_PARAMETERS_START_EX(ZEND_PARSE_PARAMS_THROW, 16, 16)
        Z_PARAM_LONG(order)
        Z_PARAM_LONG(side)
        Z_PARAM_LONG(uplo)
        Z_PARAM_LONG(m)
        Z_PARAM_LONG(n)
        Z_PARAM_DOUBLE(alpha)
        Z_PARAM_OBJECT(a) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetA)
        Z_PARAM_LONG(ldA)
        Z_PARAM_OBJECT(b) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetB)
        Z_PARAM_LONG(ldB)
        Z_PARAM_DOUBLE(beta)
        Z_PARAM_OBJECT(c) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetC)
        Z_PARAM_LONG(ldC)
    ZEND_PARSE_PARAMETERS_END();

    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_M, m)) {
        return;
    }
    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_N, n)) {
        return;
    }
    // Check Buffer A
    bufferA = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(a);
    if(php_rindow_openblas_assert_buffer_type(bufferA,"a")) {
        return;
    }
    if(!php_interop_polite_math_matrix_is_linear_buffer(bufferA)) {
        zend_throw_exception(zend_ce_type_error, "invalid type.", 0);
        return;
    }
    {
        zend_long rows;
        if(side==CblasLeft) {
            rows = m;
        } else if(side==CblasRight) {
            rows = n;
        } else {
            zend_throw_exception(spl_ce_RuntimeException, "unknown side mode for bufferA.", 0);
            return;
        }
        if(php_rindow_openblas_assert_matrix_buffer_spec(
            PHP_RINDOW_OPENBLAS_ASSERT_A, bufferA,rows,rows,offsetA,ldA)) {
            return;
        }
    }

    // Check Buffer B
    bufferB = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(b);
    if(php_rindow_openblas_assert_buffer_type(bufferB,"b")) {
        return;
    }
    {
        if(php_rindow_openblas_assert_matrix_buffer_spec(
            PHP_RINDOW_OPENBLAS_ASSERT_B, bufferB,m,n,offsetB,ldB)) {
            return;
        }
    }

    // Check Buffer C
    bufferC = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(c);
    if(php_rindow_openblas_assert_buffer_type(bufferC,"c")) {
        return;
    }
    if(php_rindow_openblas_assert_matrix_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_C, bufferC,m,n,offsetC,ldC)) {
        return;
    }

    // Check Buffer A and B and C
    if(bufferA->dtype!=bufferB->dtype || bufferB->dtype!=bufferC->dtype) {
        zend_throw_exception(spl_ce_InvalidArgumentException, "Unmatch data type for A and B and C", 0);
        return;
    }

    switch (bufferA->dtype) {
        case php_interop_polite_math_matrix_dtype_float32:
            cblas_ssymm(
                (OPENBLAS_CONST enum CBLAS_ORDER)order,
                (OPENBLAS_CONST enum CBLAS_SIDE)side,
                (OPENBLAS_CONST enum CBLAS_UPLO)uplo,
                (blasint)m,(blasint)n,
                (float)alpha,
                &(((float *)bufferA->data)[offsetA]), (blasint)ldA,
                &(((float *)bufferB->data)[offsetB]), (blasint)ldB,
                (float)beta,
                &(((float *)bufferC->data)[offsetC]), (blasint)ldC);
            break;
        case php_interop_polite_math_matrix_dtype_float64:
            cblas_dsymm(
                (OPENBLAS_CONST enum CBLAS_ORDER)order,
                (OPENBLAS_CONST enum CBLAS_SIDE)side,
                (OPENBLAS_CONST enum CBLAS_UPLO)uplo,
                (blasint)m,(blasint)n,
                (double)alpha,
                &(((double *)bufferA->data)[offsetA]), (blasint)ldA,
                &(((double *)bufferB->data)[offsetB]), (blasint)ldB,
                (double)beta,
                &(((double *)bufferC->data)[offsetC]), (blasint)ldC);
            break;
        default:
            zend_throw_exception(spl_ce_RuntimeException, "Unsupported data type.", 0);
            return;
    }
}
/* }}} */

/* Method Rindow\OpenBLAS\Blas::
    public function syrk(
        int $order,
        int $uplo,
        int $trans,
        int $n,
        int $k,
        float $alpha,
        Buffer $A, int $offsetA, int $ldA,
        float $beta,
        Buffer $C, int $offsetC, int $ldC ) : void
 {{{ */
static PHP_METHOD(Blas, syrk)
{
    php_interop_polite_math_matrix_linear_buffer_t* bufferA;
    php_interop_polite_math_matrix_linear_buffer_t* bufferC;
    zend_long order;
    zend_long uplo;
    zend_long trans;
    zend_long n;
    zend_long k;
    double alpha;
    zval* a=NULL;
    zend_long offsetA;
    zend_long ldA;
    double beta;
    zval* c=NULL;
    zend_long offsetC;
    zend_long ldC;

    ZEND_PARSE_PARAMETERS_START_EX(ZEND_PARSE_PARAMS_THROW, 13, 13)
        Z_PARAM_LONG(order)
        Z_PARAM_LONG(uplo)
        Z_PARAM_LONG(trans)
        Z_PARAM_LONG(n)
        Z_PARAM_LONG(k)
        Z_PARAM_DOUBLE(alpha)
        Z_PARAM_OBJECT(a) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetA)
        Z_PARAM_LONG(ldA)
        Z_PARAM_DOUBLE(beta)
        Z_PARAM_OBJECT(c) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetC)
        Z_PARAM_LONG(ldC)
    ZEND_PARSE_PARAMETERS_END();

    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_N, n)) {
        return;
    }
    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_K, k)) {
        return;
    }
    // Check Buffer A
    bufferA = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(a);
    if(php_rindow_openblas_assert_buffer_type(bufferA,"a")) {
        return;
    }
    {
        zend_long rows,cols;
        if(trans==CblasNoTrans || trans==CblasConjNoTrans) {
            rows = n; cols = k;
        } else if(trans==CblasTrans || trans==CblasConjTrans) {
            rows = k; cols = n;
        } else {
            zend_throw_exception(spl_ce_RuntimeException, "unknown transpose mode for bufferA.", 0);
            return;
        }
        if(php_rindow_openblas_assert_matrix_buffer_spec(
            PHP_RINDOW_OPENBLAS_ASSERT_A, bufferA,rows,cols,offsetA,ldA)) {
            return;
        }
    }

    // Check Buffer C
    bufferC = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(c);
    if(php_rindow_openblas_assert_buffer_type(bufferC,"c")) {
        return;
    }
    if(php_rindow_openblas_assert_matrix_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_C, bufferC,n,n,offsetC,ldC)) {
        return;
    }

    // Check Buffer A and C
    if(bufferA->dtype!=bufferC->dtype) {
        zend_throw_exception(spl_ce_InvalidArgumentException, "Unmatch data type for A and C", 0);
        return;
    }

    switch (bufferA->dtype) {
        case php_interop_polite_math_matrix_dtype_float32:
            cblas_ssyrk(
                (OPENBLAS_CONST enum CBLAS_ORDER)order,
                (OPENBLAS_CONST enum CBLAS_UPLO)uplo,
                (OPENBLAS_CONST enum CBLAS_TRANSPOSE)trans,
                (blasint)n,(blasint)k,
                (float)alpha,
                &(((float *)bufferA->data)[offsetA]), (blasint)ldA,
                (float)beta,
                &(((float *)bufferC->data)[offsetC]), (blasint)ldC);
            break;
        case php_interop_polite_math_matrix_dtype_float64:
            cblas_dsyrk(
                (OPENBLAS_CONST enum CBLAS_ORDER)order,
                (OPENBLAS_CONST enum CBLAS_UPLO)uplo,
                (OPENBLAS_CONST enum CBLAS_TRANSPOSE)trans,
                (blasint)n,(blasint)k,
                (double)alpha,
                &(((double *)bufferA->data)[offsetA]), (blasint)ldA,
                (double)beta,
                &(((double *)bufferC->data)[offsetC]), (blasint)ldC);
            break;
        default:
            zend_throw_exception(spl_ce_RuntimeException, "Unsupported data type.", 0);
            return;
    }
}
/* }}} */

/* Method Rindow\OpenBLAS\Blas::
    public function syr2k(
        int $order,
        int $uplo,
        int $trans,
        int $n,
        int $k,
        float $alpha,
        Buffer $A, int $offsetA, int $ldA,
        Buffer $B, int $offsetB, int $ldB,
        float $beta,
        Buffer $C, int $offsetC, int $ldC ) : void
 {{{ */
static PHP_METHOD(Blas, syr2k)
{
    php_interop_polite_math_matrix_linear_buffer_t* bufferA;
    php_interop_polite_math_matrix_linear_buffer_t* bufferB;
    php_interop_polite_math_matrix_linear_buffer_t* bufferC;
    zend_long order;
    zend_long uplo;
    zend_long trans;
    zend_long n;
    zend_long k;
    double alpha;
    zval* a=NULL;
    zend_long offsetA;
    zend_long ldA;
    zval* b=NULL;
    zend_long offsetB;
    zend_long ldB;
    double beta;
    zval* c=NULL;
    zend_long offsetC;
    zend_long ldC;

    ZEND_PARSE_PARAMETERS_START_EX(ZEND_PARSE_PARAMS_THROW, 16, 16)
        Z_PARAM_LONG(order)
        Z_PARAM_LONG(uplo)
        Z_PARAM_LONG(trans)
        Z_PARAM_LONG(n)
        Z_PARAM_LONG(k)
        Z_PARAM_DOUBLE(alpha)
        Z_PARAM_OBJECT(a) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetA)
        Z_PARAM_LONG(ldA)
        Z_PARAM_OBJECT(b) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetB)
        Z_PARAM_LONG(ldB)
        Z_PARAM_DOUBLE(beta)
        Z_PARAM_OBJECT(c) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetC)
        Z_PARAM_LONG(ldC)
    ZEND_PARSE_PARAMETERS_END();

    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_N, n)) {
        return;
    }
    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_M, k)) {
        return;
    }
    // Check Buffer A and B
    bufferA = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(a);
    if(php_rindow_openblas_assert_buffer_type(bufferA,"a")) {
        return;
    }
    bufferB = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(b);
    if(php_rindow_openblas_assert_buffer_type(bufferB,"b")) {
        return;
    }
    {
        zend_long rows,cols;
        if(trans==CblasNoTrans || trans==CblasConjNoTrans) {
            rows = n; cols = k;
        } else if(trans==CblasTrans || trans==CblasConjTrans) {
            rows = k; cols = n;
        } else {
            zend_throw_exception(spl_ce_RuntimeException, "unknown transpose mode for bufferA.", 0);
            return;
        }
        if(php_rindow_openblas_assert_matrix_buffer_spec(
            PHP_RINDOW_OPENBLAS_ASSERT_A, bufferA,rows,cols,offsetA,ldA)) {
            return;
        }
        if(php_rindow_openblas_assert_matrix_buffer_spec(
            PHP_RINDOW_OPENBLAS_ASSERT_B, bufferB,rows,cols,offsetB,ldB)) {
            return;
        }
    }

    // Check Buffer C
    bufferC = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(c);
    if(php_rindow_openblas_assert_buffer_type(bufferC,"c")) {
        return;
    }
    if(php_rindow_openblas_assert_matrix_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_C, bufferC,n,n,offsetC,ldC)) {
        return;
    }

    // Check Buffer A and B and C
    if(bufferA->dtype!=bufferB->dtype || bufferB->dtype!=bufferC->dtype) {
        zend_throw_exception(spl_ce_InvalidArgumentException, "Unmatch data type for A and B and C", 0);
        return;
    }

    switch (bufferA->dtype) {
        case php_interop_polite_math_matrix_dtype_float32:
            cblas_ssyr2k(
                (OPENBLAS_CONST enum CBLAS_ORDER)order,
                (OPENBLAS_CONST enum CBLAS_UPLO)uplo,
                (OPENBLAS_CONST enum CBLAS_TRANSPOSE)trans,
                (blasint)n,(blasint)k,
                (float)alpha,
                &(((float *)bufferA->data)[offsetA]), (blasint)ldA,
                &(((float *)bufferB->data)[offsetB]), (blasint)ldB,
                (float)beta,
                &(((float *)bufferC->data)[offsetC]), (blasint)ldC);
            break;
        case php_interop_polite_math_matrix_dtype_float64:
            cblas_dsyr2k(
                (OPENBLAS_CONST enum CBLAS_ORDER)order,
                (OPENBLAS_CONST enum CBLAS_UPLO)uplo,
                (OPENBLAS_CONST enum CBLAS_TRANSPOSE)trans,
                (blasint)n,(blasint)k,
                (double)alpha,
                &(((double *)bufferA->data)[offsetA]), (blasint)ldA,
                &(((double *)bufferB->data)[offsetB]), (blasint)ldB,
                (double)beta,
                &(((double *)bufferC->data)[offsetC]), (blasint)ldC);
            break;
        default:
            zend_throw_exception(spl_ce_RuntimeException, "Unsupported data type.", 0);
            return;
    }
}
/* }}} */

/* Method Rindow\OpenBLAS\Blas::
    public function trmm(
        int $order,
        int $side,
        int $uplo,
        int $trans,
        int $diag,
        int $m,
        int $n,
        float $alpha,
        Buffer $A, int $offsetA, int $ldA,
        Buffer $B, int $offsetB, int $ldB) : void
 {{{ */
static PHP_METHOD(Blas, trmm)
{
    php_interop_polite_math_matrix_linear_buffer_t* bufferA;
    php_interop_polite_math_matrix_linear_buffer_t* bufferB;
    zend_long order;
    zend_long side;
    zend_long uplo;
    zend_long trans;
    zend_long diag;
    zend_long m;
    zend_long n;
    double alpha;
    zval* a=NULL;
    zend_long offsetA;
    zend_long ldA;
    zval* b=NULL;
    zend_long offsetB;
    zend_long ldB;

    ZEND_PARSE_PARAMETERS_START_EX(ZEND_PARSE_PARAMS_THROW, 14, 14)
        Z_PARAM_LONG(order)
        Z_PARAM_LONG(side)
        Z_PARAM_LONG(uplo)
        Z_PARAM_LONG(trans)
        Z_PARAM_LONG(diag)
        Z_PARAM_LONG(m)
        Z_PARAM_LONG(n)
        Z_PARAM_DOUBLE(alpha)
        Z_PARAM_OBJECT(a) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetA)
        Z_PARAM_LONG(ldA)
        Z_PARAM_OBJECT(b) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetB)
        Z_PARAM_LONG(ldB)
    ZEND_PARSE_PARAMETERS_END();

    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_M, m)) {
        return;
    }
    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_N, n)) {
        return;
    }
    // Check Buffer A
    bufferA = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(a);
    if(php_rindow_openblas_assert_buffer_type(bufferA,"a")) {
        return;
    }
    {
        zend_long sizeA;
        if(side==CblasLeft) {
            sizeA = m;
        } else if(side==CblasRight) {
            sizeA = n;
        } else {
            zend_throw_exception(spl_ce_RuntimeException, "unknown transpose mode for bufferA.", 0);
            return;
        }
        if(php_rindow_openblas_assert_matrix_buffer_spec(
            PHP_RINDOW_OPENBLAS_ASSERT_A, bufferA,sizeA,sizeA,offsetA,ldA)) {
            return;
        }
    }
    bufferB = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(b);
    if(php_rindow_openblas_assert_buffer_type(bufferB,"b")) {
        return;
    }
    if(php_rindow_openblas_assert_matrix_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_B, bufferB,m,n,offsetB,ldB)) {
        return;
    }

    // Check Buffer A and B
    if(bufferA->dtype!=bufferB->dtype) {
        zend_throw_exception(spl_ce_InvalidArgumentException, "Unmatch data type for A and B", 0);
        return;
    }

    switch (bufferA->dtype) {
        case php_interop_polite_math_matrix_dtype_float32:
            cblas_strmm(
                (OPENBLAS_CONST enum CBLAS_ORDER)order,
                (OPENBLAS_CONST enum CBLAS_SIDE)side,
                (OPENBLAS_CONST enum CBLAS_UPLO)uplo,
                (OPENBLAS_CONST enum CBLAS_TRANSPOSE)trans,
                (OPENBLAS_CONST enum CBLAS_DIAG)diag,
                (blasint)m,(blasint)n,
                (float)alpha,
                &(((float *)bufferA->data)[offsetA]), (blasint)ldA,
                &(((float *)bufferB->data)[offsetB]), (blasint)ldB);
            break;
        case php_interop_polite_math_matrix_dtype_float64:
            cblas_dtrmm(
                (OPENBLAS_CONST enum CBLAS_ORDER)order,
                (OPENBLAS_CONST enum CBLAS_SIDE)side,
                (OPENBLAS_CONST enum CBLAS_UPLO)uplo,
                (OPENBLAS_CONST enum CBLAS_TRANSPOSE)trans,
                (OPENBLAS_CONST enum CBLAS_DIAG)diag,
                (blasint)m,(blasint)n,
                (double)alpha,
                &(((double *)bufferA->data)[offsetA]), (blasint)ldA,
                &(((double *)bufferB->data)[offsetB]), (blasint)ldB);
            break;
        default:
            zend_throw_exception(spl_ce_RuntimeException, "Unsupported data type.", 0);
            return;
    }
}
/* }}} */

/* Method Rindow\OpenBLAS\Blas::
    public function trsm(
        int $order,
        int $side,
        int $uplo,
        int $trans,
        int $diag,
        int $m,
        int $n,
        float $alpha,
        Buffer $A, int $offsetA, int $ldA,
        Buffer $B, int $offsetB, int $ldB) : void
 {{{ */
static PHP_METHOD(Blas, trsm)
{
    php_interop_polite_math_matrix_linear_buffer_t* bufferA;
    php_interop_polite_math_matrix_linear_buffer_t* bufferB;
    zend_long order;
    zend_long side;
    zend_long uplo;
    zend_long trans;
    zend_long diag;
    zend_long m;
    zend_long n;
    double alpha;
    zval* a=NULL;
    zend_long offsetA;
    zend_long ldA;
    zval* b=NULL;
    zend_long offsetB;
    zend_long ldB;

    ZEND_PARSE_PARAMETERS_START_EX(ZEND_PARSE_PARAMS_THROW, 14, 14)
        Z_PARAM_LONG(order)
        Z_PARAM_LONG(side)
        Z_PARAM_LONG(uplo)
        Z_PARAM_LONG(trans)
        Z_PARAM_LONG(diag)
        Z_PARAM_LONG(m)
        Z_PARAM_LONG(n)
        Z_PARAM_DOUBLE(alpha)
        Z_PARAM_OBJECT(a) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetA)
        Z_PARAM_LONG(ldA)
        Z_PARAM_OBJECT(b) // Interop\Polite\Math\Matrix\LinearBuffer
        Z_PARAM_LONG(offsetB)
        Z_PARAM_LONG(ldB)
    ZEND_PARSE_PARAMETERS_END();

    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_M, m)) {
        return;
    }
    if(php_rindow_openblas_assert_shape_parameter(
        PHP_RINDOW_OPENBLAS_ASSERT_N, n)) {
        return;
    }
    // Check Buffer A
    bufferA = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(a);
    if(php_rindow_openblas_assert_buffer_type(bufferA,"a")) {
        return;
    }
    {
        zend_long sizeA;
        if(side==CblasLeft) {
            sizeA = m;
        } else if(side==CblasRight) {
            sizeA = n;
        } else {
            zend_throw_exception(spl_ce_RuntimeException, "unknown transpose mode for bufferA.", 0);
            return;
        }
        if(php_rindow_openblas_assert_matrix_buffer_spec(
            PHP_RINDOW_OPENBLAS_ASSERT_A, bufferA,sizeA,sizeA,offsetA,ldA)) {
            return;
        }
    }
    bufferB = Z_INTEROP_POLITE_MATH_MATRIX_LINEAR_BUFFER_OBJ_P(b);
    if(php_rindow_openblas_assert_buffer_type(bufferB,"b")) {
        return;
    }
    if(php_rindow_openblas_assert_matrix_buffer_spec(
        PHP_RINDOW_OPENBLAS_ASSERT_B, bufferB,m,n,offsetB,ldB)) {
        return;
    }

    // Check Buffer A and B
    if(bufferA->dtype!=bufferB->dtype) {
        zend_throw_exception(spl_ce_InvalidArgumentException, "Unmatch data type for A and B", 0);
        return;
    }

    switch (bufferA->dtype) {
        case php_interop_polite_math_matrix_dtype_float32:
            cblas_strsm(
                (OPENBLAS_CONST enum CBLAS_ORDER)order,
                (OPENBLAS_CONST enum CBLAS_SIDE)side,
                (OPENBLAS_CONST enum CBLAS_UPLO)uplo,
                (OPENBLAS_CONST enum CBLAS_TRANSPOSE)trans,
                (OPENBLAS_CONST enum CBLAS_DIAG)diag,
                (blasint)m,(blasint)n,
                (float)alpha,
                &(((float *)bufferA->data)[offsetA]), (blasint)ldA,
                &(((float *)bufferB->data)[offsetB]), (blasint)ldB);
            break;
        case php_interop_polite_math_matrix_dtype_float64:
            cblas_dtrsm(
                (OPENBLAS_CONST enum CBLAS_ORDER)order,
                (OPENBLAS_CONST enum CBLAS_SIDE)side,
                (OPENBLAS_CONST enum CBLAS_UPLO)uplo,
                (OPENBLAS_CONST enum CBLAS_TRANSPOSE)trans,
                (OPENBLAS_CONST enum CBLAS_DIAG)diag,
                (blasint)m,(blasint)n,
                (double)alpha,
                &(((double *)bufferA->data)[offsetA]), (blasint)ldA,
                &(((double *)bufferB->data)[offsetB]), (blasint)ldB);
            break;
        default:
            zend_throw_exception(spl_ce_RuntimeException, "Unsupported data type.", 0);
            return;
    }
}
/* }}} */

ZEND_BEGIN_ARG_INFO_EX(ai_Blas_scal, 0, 0, 5)
    ZEND_ARG_INFO(0, n)
    ZEND_ARG_INFO(0, alpha)
    ZEND_ARG_OBJ_INFO(0, x, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetX)
    ZEND_ARG_INFO(0, incX)
ZEND_END_ARG_INFO()

ZEND_BEGIN_ARG_INFO_EX(ai_Blas_axpy, 0, 0, 8)
    ZEND_ARG_INFO(0, n)
    ZEND_ARG_INFO(0, alpha)
    ZEND_ARG_OBJ_INFO(0, x, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetX)
    ZEND_ARG_INFO(0, incX)
    ZEND_ARG_OBJ_INFO(0, y, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetY)
    ZEND_ARG_INFO(0, incY)
ZEND_END_ARG_INFO()

ZEND_BEGIN_ARG_INFO_EX(ai_Blas_dot, 0, 0, 7)
    ZEND_ARG_INFO(0, n)
    ZEND_ARG_OBJ_INFO(0, x, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetX)
    ZEND_ARG_INFO(0, incX)
    ZEND_ARG_OBJ_INFO(0, y, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetY)
    ZEND_ARG_INFO(0, incY)
ZEND_END_ARG_INFO()

ZEND_BEGIN_ARG_INFO_EX(ai_Blas_asum, 0, 0, 4)
    ZEND_ARG_INFO(0, n)
    ZEND_ARG_OBJ_INFO(0, x, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetX)
    ZEND_ARG_INFO(0, incX)
ZEND_END_ARG_INFO()

ZEND_BEGIN_ARG_INFO_EX(ai_Blas_iamax, 0, 0, 4)
    ZEND_ARG_INFO(0, n)
    ZEND_ARG_OBJ_INFO(0, x, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetX)
    ZEND_ARG_INFO(0, incX)
ZEND_END_ARG_INFO()

#ifdef OPENBLAS_HAVE_IAMIN
ZEND_BEGIN_ARG_INFO_EX(ai_Blas_iamin, 0, 0, 4)
    ZEND_ARG_INFO(0, n)
    ZEND_ARG_OBJ_INFO(0, x, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetX)
    ZEND_ARG_INFO(0, incX)
ZEND_END_ARG_INFO()
#endif

ZEND_BEGIN_ARG_INFO_EX(ai_Blas_copy, 0, 0, 7)
    ZEND_ARG_INFO(0, n)
    ZEND_ARG_OBJ_INFO(0, x, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetX)
    ZEND_ARG_INFO(0, incX)
    ZEND_ARG_OBJ_INFO(0, y, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetY)
    ZEND_ARG_INFO(0, incY)
ZEND_END_ARG_INFO()

ZEND_BEGIN_ARG_INFO_EX(ai_Blas_nrm2, 0, 0, 4)
    ZEND_ARG_INFO(0, n)
    ZEND_ARG_OBJ_INFO(0, x, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetX)
    ZEND_ARG_INFO(0, incX)
ZEND_END_ARG_INFO()

ZEND_BEGIN_ARG_INFO_EX(ai_Blas_rotg, 0, 0, 8)
    ZEND_ARG_OBJ_INFO(0, a, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetA)
    ZEND_ARG_OBJ_INFO(0, b, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetB)
    ZEND_ARG_OBJ_INFO(0, c, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetC)
    ZEND_ARG_OBJ_INFO(0, s, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetS)
ZEND_END_ARG_INFO()

ZEND_BEGIN_ARG_INFO_EX(ai_Blas_rot, 0, 0, 11)
    ZEND_ARG_INFO(0, n)
    ZEND_ARG_OBJ_INFO(0, x, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetX)
    ZEND_ARG_INFO(0, incX)
    ZEND_ARG_OBJ_INFO(0, y, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetY)
    ZEND_ARG_INFO(0, incY)
    ZEND_ARG_OBJ_INFO(0, c, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetC)
    ZEND_ARG_OBJ_INFO(0, s, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetS)
ZEND_END_ARG_INFO()

ZEND_BEGIN_ARG_INFO_EX(ai_Blas_swap, 0, 0, 7)
    ZEND_ARG_INFO(0, n)
    ZEND_ARG_OBJ_INFO(0, x, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetX)
    ZEND_ARG_INFO(0, incX)
    ZEND_ARG_OBJ_INFO(0, y, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetY)
    ZEND_ARG_INFO(0, incY)
ZEND_END_ARG_INFO()

ZEND_BEGIN_ARG_INFO_EX(ai_Blas_gemv, 0, 0, 15)
    ZEND_ARG_INFO(0, order)
    ZEND_ARG_INFO(0, trans)
    ZEND_ARG_INFO(0, m)
    ZEND_ARG_INFO(0, n)
    ZEND_ARG_INFO(0, alpha)
    ZEND_ARG_OBJ_INFO(0, a, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetA)
    ZEND_ARG_INFO(0, ldA)
    ZEND_ARG_OBJ_INFO(0, x, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetX)
    ZEND_ARG_INFO(0, incX)
    ZEND_ARG_INFO(0, beta)
    ZEND_ARG_OBJ_INFO(0, x, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetY)
    ZEND_ARG_INFO(0, incY)
ZEND_END_ARG_INFO()

ZEND_BEGIN_ARG_INFO_EX(ai_Blas_gemm, 0, 0, 17)
    ZEND_ARG_INFO(0, order)
    ZEND_ARG_INFO(0, transA)
    ZEND_ARG_INFO(0, transB)
    ZEND_ARG_INFO(0, m)
    ZEND_ARG_INFO(0, n)
    ZEND_ARG_INFO(0, k)
    ZEND_ARG_INFO(0, alpha)
    ZEND_ARG_OBJ_INFO(0, a, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetA)
    ZEND_ARG_INFO(0, ldA)
    ZEND_ARG_OBJ_INFO(0, b, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetB)
    ZEND_ARG_INFO(0, ldB)
    ZEND_ARG_INFO(0, beta)
    ZEND_ARG_OBJ_INFO(0, c, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetC)
    ZEND_ARG_INFO(0, ldC)
ZEND_END_ARG_INFO()

ZEND_BEGIN_ARG_INFO_EX(ai_Blas_symm, 0, 0, 16)
    ZEND_ARG_INFO(0, order)
    ZEND_ARG_INFO(0, side)
    ZEND_ARG_INFO(0, uplo)
    ZEND_ARG_INFO(0, m)
    ZEND_ARG_INFO(0, n)
    ZEND_ARG_INFO(0, alpha)
    ZEND_ARG_OBJ_INFO(0, a, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetA)
    ZEND_ARG_INFO(0, ldA)
    ZEND_ARG_OBJ_INFO(0, b, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetB)
    ZEND_ARG_INFO(0, ldB)
    ZEND_ARG_INFO(0, beta)
    ZEND_ARG_OBJ_INFO(0, c, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetC)
    ZEND_ARG_INFO(0, ldC)
ZEND_END_ARG_INFO()

ZEND_BEGIN_ARG_INFO_EX(ai_Blas_syrk, 0, 0, 13)
    ZEND_ARG_INFO(0, order)
    ZEND_ARG_INFO(0, uplo)
    ZEND_ARG_INFO(0, trans)
    ZEND_ARG_INFO(0, n)
    ZEND_ARG_INFO(0, k)
    ZEND_ARG_INFO(0, alpha)
    ZEND_ARG_OBJ_INFO(0, a, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetA)
    ZEND_ARG_INFO(0, ldA)
    ZEND_ARG_INFO(0, beta)
    ZEND_ARG_OBJ_INFO(0, c, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetC)
    ZEND_ARG_INFO(0, ldC)
ZEND_END_ARG_INFO()

ZEND_BEGIN_ARG_INFO_EX(ai_Blas_syr2k, 0, 0, 16)
    ZEND_ARG_INFO(0, order)
    ZEND_ARG_INFO(0, uplo)
    ZEND_ARG_INFO(0, trans)
    ZEND_ARG_INFO(0, n)
    ZEND_ARG_INFO(0, k)
    ZEND_ARG_INFO(0, alpha)
    ZEND_ARG_OBJ_INFO(0, a, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetA)
    ZEND_ARG_INFO(0, ldA)
    ZEND_ARG_OBJ_INFO(0, b, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetB)
    ZEND_ARG_INFO(0, ldB)
    ZEND_ARG_INFO(0, beta)
    ZEND_ARG_OBJ_INFO(0, c, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetC)
    ZEND_ARG_INFO(0, ldC)
ZEND_END_ARG_INFO()

ZEND_BEGIN_ARG_INFO_EX(ai_Blas_trmm, 0, 0, 14)
    ZEND_ARG_INFO(0, order)
    ZEND_ARG_INFO(0, side)
    ZEND_ARG_INFO(0, uplo)
    ZEND_ARG_INFO(0, trans)
    ZEND_ARG_INFO(0, diag)
    ZEND_ARG_INFO(0, m)
    ZEND_ARG_INFO(0, n)
    ZEND_ARG_INFO(0, alpha)
    ZEND_ARG_OBJ_INFO(0, a, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetA)
    ZEND_ARG_INFO(0, ldA)
    ZEND_ARG_OBJ_INFO(0, b, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetB)
    ZEND_ARG_INFO(0, ldB)
ZEND_END_ARG_INFO()

ZEND_BEGIN_ARG_INFO_EX(ai_Blas_trsm, 0, 0, 14)
    ZEND_ARG_INFO(0, order)
    ZEND_ARG_INFO(0, side)
    ZEND_ARG_INFO(0, uplo)
    ZEND_ARG_INFO(0, trans)
    ZEND_ARG_INFO(0, diag)
    ZEND_ARG_INFO(0, m)
    ZEND_ARG_INFO(0, n)
    ZEND_ARG_INFO(0, alpha)
    ZEND_ARG_OBJ_INFO(0, a, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetA)
    ZEND_ARG_INFO(0, ldA)
    ZEND_ARG_OBJ_INFO(0, b, Interop\\Polite\\Math\\Matrix\\LinearBuffer, 0)
    ZEND_ARG_INFO(0, offsetB)
    ZEND_ARG_INFO(0, ldB)
ZEND_END_ARG_INFO()

ZEND_BEGIN_ARG_INFO_EX(ai_Blas_void, 0, 0, 0)
ZEND_END_ARG_INFO()

/* {{{ Rindow\OpenBLAS\Blas function entries */
static zend_function_entry php_rindow_openblas_blas_me[] = {
    /* clang-format off */
    PHP_ME(Blas, getNumThreads, ai_Blas_void, ZEND_ACC_PUBLIC)
    PHP_ME(Blas, getNumProcs,   ai_Blas_void, ZEND_ACC_PUBLIC)
    PHP_ME(Blas, getConfig,     ai_Blas_void, ZEND_ACC_PUBLIC)
    PHP_ME(Blas, getCorename,   ai_Blas_void, ZEND_ACC_PUBLIC)
    PHP_ME(Blas, scal,  ai_Blas_scal,  ZEND_ACC_PUBLIC)
    PHP_ME(Blas, axpy,  ai_Blas_axpy,  ZEND_ACC_PUBLIC)
    PHP_ME(Blas, dot,   ai_Blas_dot,   ZEND_ACC_PUBLIC)
    PHP_ME(Blas, asum,  ai_Blas_asum,  ZEND_ACC_PUBLIC)
    PHP_ME(Blas, iamax, ai_Blas_iamax, ZEND_ACC_PUBLIC)
#ifdef OPENBLAS_HAVE_IAMIN
    PHP_ME(Blas, iamin, ai_Blas_iamin, ZEND_ACC_PUBLIC)
#endif
    PHP_ME(Blas, copy,  ai_Blas_copy,  ZEND_ACC_PUBLIC)
    PHP_ME(Blas, nrm2,  ai_Blas_nrm2,  ZEND_ACC_PUBLIC)
    PHP_ME(Blas, rotg,  ai_Blas_rotg,  ZEND_ACC_PUBLIC)
    PHP_ME(Blas, rot,   ai_Blas_rot,   ZEND_ACC_PUBLIC)
    PHP_ME(Blas, swap,  ai_Blas_swap,  ZEND_ACC_PUBLIC)
    PHP_ME(Blas, gemv,  ai_Blas_gemv,  ZEND_ACC_PUBLIC)
    PHP_ME(Blas, gemm,  ai_Blas_gemm,  ZEND_ACC_PUBLIC)
    PHP_ME(Blas, symm,  ai_Blas_symm,  ZEND_ACC_PUBLIC)
    PHP_ME(Blas, syrk,  ai_Blas_syrk,  ZEND_ACC_PUBLIC)
    PHP_ME(Blas, syr2k, ai_Blas_syr2k, ZEND_ACC_PUBLIC)
    PHP_ME(Blas, trmm,  ai_Blas_trmm,  ZEND_ACC_PUBLIC)
    PHP_ME(Blas, trsm,  ai_Blas_trsm,  ZEND_ACC_PUBLIC)
    PHP_FE_END
    /* clang-format on */
};
/* }}} */

/* Class Rindow\OpenBLAS\Blas {{{ */
static zend_class_entry* rindow_openblas_blas_ce;

void php_rindow_openblas_blas_init_ce(INIT_FUNC_ARGS)
{
    zend_class_entry ce;

    INIT_NS_CLASS_ENTRY(ce, "Rindow\\OpenBLAS", "Blas", php_rindow_openblas_blas_me);
    rindow_openblas_blas_ce = zend_register_internal_class(&ce);
    rindow_openblas_blas_ce->create_object = php_rindow_openblas_blas_create_object;

    memcpy(&rindow_openblas_blas_object_handlers, zend_get_std_object_handlers(), sizeof(zend_object_handlers));
    rindow_openblas_blas_object_handlers.offset    = 0;
    rindow_openblas_blas_object_handlers.free_obj  = php_rindow_openblas_blas_free_object;
    rindow_openblas_blas_object_handlers.clone_obj = NULL;

    //zend_class_implements(rindow_openblas_blas_ce, 2, spl_ce_ArrayAccess, spl_ce_Countable);
}
/* }}} */

max_stars_repo_name: yuichiis/rindow-openblas

max_issues_repo_name: yuichiis/rindow-openblas

max_forks_repo_name: yuichiis/rindow-openblas


text:
#ifndef MATHLAB_H
#define MATHLAB_H

#define USE_MKL 0
#define USE_OPENBLAS 1

#if USE_MKL
#include <mkl.h>

#elif USE_OPENBLAS
#include <cblas.h>

#endif

#endif

ext: h


text:
/**
 * Finite element method solver. This solver solves laplace(u) = 6 on square
 * [0, 1]x[0, 1]. Boundary conditions are set to be 1 + x^2 + 2y^2.
 *
 * Exact solution of this problem is 1 + x^2 + 2y^2. We will test against this
 * solution in the end.
 *
 */

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <assert.h>
#include <stdbool.h>
#include <lapacke.h>
#include <triangle.h>

#define F -6.0

#define likely(x) __builtin_expect(!!(x), 1)
#define unlikely(x) __builtin_expect(!!(x), 0)

static void compute_local_a(uint *point_id, REAL *pointlist, REAL *matrix);
static void compute_local_f(uint *point_id, REAL *pointlist, REAL *matrix);
static void generate_mesh(struct triangulateio *io, struct triangulateio *out);
static REAL boundary_condition(uint point_id, struct triangulateio *out);
static REAL exact_solution(uint point_id, struct triangulateio *out);

int main(int argc, char **argv)
{
    struct triangulateio in, out;
    uint num_of_points;
    uint num_of_triangles;
    REAL *local_A, *local_F, *global_A, *global_F;

    /* generate geometry */
    generate_mesh(&in, &out);

    num_of_points = out.numberofpoints;
    num_of_triangles = out.numberoftriangles;

    /**
     * Allocate all matrix data structures
     */
    global_A = malloc(num_of_points * num_of_points * sizeof(REAL));
    global_F = malloc(num_of_points * sizeof(REAL));
    local_A = malloc(3 * 3 * sizeof(REAL));
    local_F = malloc(3 * 3 * sizeof(REAL));

    assert(global_A != NULL);
    assert(global_F != NULL);
    assert(local_A != NULL);
    assert(local_F != NULL);

    memset(global_A, 0, num_of_points * num_of_points * sizeof(REAL));
    memset(global_F, 0, num_of_points * sizeof(REAL));

    /**
     * Matrix assembly without boundary conditions
     */
    for (size_t n = 0; n < num_of_triangles; n++) {
        // three points that describe triangle
        uint point_id[3];
        memcpy(point_id, &out.trianglelist[n * 3], 3 * sizeof(int));
        compute_local_a(point_id, out.pointlist, local_A);
        compute_local_f(point_id, out.pointlist, local_F);

        for (size_t i = 0; i < 3; i++) {
            int i1 = point_id[i];
            for (int j = 0; j < 3; j++) {
                int j1 = point_id[j];
                global_A[i1 * num_of_points + j1] += local_A[i * 3 + j];
                global_F[i1] += local_F[i * 3 + j] * F;
            }
        }
    }

    for (size_t i = 0; i < num_of_points; i++) {
        if (out.pointmarkerlist[i] == 1) {
            for (size_t j = 0; j < num_of_points; j++) {
                global_F[j] -= global_A[j * num_of_points + i]
                    * boundary_condition(i, &out);
                global_A[i * num_of_points + j] = 0;
                global_A[j * num_of_points + i] = 0;
            }
            global_A[i * num_of_points + i] = 1.0;
            global_F[i] = boundary_condition(i, &out);
        }
    }

    /*************************************************************************
     ************************ Solve Ax = b ***********************************
     *************************************************************************/
    int ipiv[num_of_points];
    int info = LAPACKE_dgesv(LAPACK_ROW_MAJOR,
                             num_of_points,
                             1,
                             global_A,
                             num_of_points,
                             ipiv,
                             global_F,
                             1);
    if (info > 0) {
        printf("The factorization has been completed, but the factor U is\n");
        printf("exactly singular, so the solution could not be computed.");
        goto mem_free;
    } else if (info < 0) {
        printf("%d-th argument had an illegal value\n", -info);
        goto mem_free;
    }

    //TODO: implement max error
    printf("\nSolution in point vs exact solution\n");
    for (size_t i = 0; i < num_of_points; i++) {
        double temp = exact_solution(i, &out);
        printf(
            "%lf  \t  %lf  \t  %lf\n", temp, global_F[i], temp - global_F[i]);
    }

mem_free:
    /* free all Triangle data structures */
    free(out.edgelist);
    free(out.edgemarkerlist);
    free(out.neighborlist);
    free(out.pointlist);
    free(out.trianglelist);
    free(out.segmentlist);
    free(out.segmentmarkerlist);
    free(out.pointmarkerlist);

    /* free all input data structures */
    free(in.regionlist);
    free(in.pointlist);

    /* free all matrix arrays */
    free(global_A);
    free(local_A);
    free(global_F);
    free(local_F);

    return info;
}

static void generate_mesh(struct triangulateio *in, struct triangulateio *out)
{
    in->numberofpoints = 4;
    in->numberofpointattributes = 0;
    in->pointattributelist = (REAL *) NULL;
    in->pointlist = malloc(in->numberofpoints * 2 * sizeof(REAL));

    in->pointlist[0] = 0.0;
    in->pointlist[1] = 0.0;

    in->pointlist[2] = 0.0;
    in->pointlist[3] = 1.0;

    in->pointlist[4] = 1.0;
    in->pointlist[5] = 1.0;

    in->pointlist[6] = 1.0;
    in->pointlist[7] = 0.0;

    in->pointmarkerlist = (int *) NULL;

    in->numberofsegments = 0;
    in->numberofholes = 0;
    in->numberofregions = 1;
    in->regionlist = malloc(in->numberofregions * 4 * sizeof(REAL));
    in->regionlist[0] = 0.5;
    in->regionlist[1] = 0.5;
    in->regionlist[2] = 7.0; /* Regional attribute (for whole mesh). */
    in->regionlist[3] = 0.1; /* Area constraint. */

    out->pointlist = (REAL *) NULL;
    out->pointmarkerlist = (int *) NULL;
    out->trianglelist = (int *) NULL;
    out->triangleattributelist = (REAL *) NULL;
    out->neighborlist = (int *) NULL;
    out->segmentlist = (int *) NULL;
    out->segmentmarkerlist = (int *) NULL;
    out->edgelist = (int *) NULL;
    out->edgemarkerlist = (int *) NULL;

    /* Triangulate the points. read and write a PSLG (p), preserve the
     * convex hull (c), zero index (z), produce edge list (e), produce
     * neighbor list (n), be quiet (Q), area should be less than 0.1 (a.1)
     * and generate triangle for FEM (q).
     */
    triangulate("pzceQna.1q", in, out, (struct triangulateio *) NULL);
}

static void compute_local_a(uint *point_id, REAL *pointlist, REAL *matrix)
{
        REAL dx23 = pointlist[2 * point_id[1]] - pointlist[2 * point_id[2]];
        REAL dx31 = pointlist[2 * point_id[2]] - pointlist[2 * point_id[0]];
        REAL dx12 = pointlist[2 * point_id[0]] - pointlist[2 * point_id[1]];
        REAL dy23 = pointlist[2 * point_id[1] + 1] - pointlist[2 * point_id[2] + 1];
        REAL dy31 = pointlist[2 * point_id[2] + 1] - pointlist[2 * point_id[0] + 1];
        REAL dy12 = pointlist[2 * point_id[0] + 1] - pointlist[2 * point_id[1] + 1];;

        REAL area = 0.5 * (dx31 * dy12 - dy31 * dx12);
        REAL _tmp_mult = 0.25 / area;
        matrix[0] = _tmp_mult * (dx23 * dx23 + dy23 * dy23);
        matrix[1] = _tmp_mult * (dx23 * dx31 + dy23 * dy31);
        matrix[2] = _tmp_mult * (dx23 * dx12 + dy23 * dy12);
        matrix[3] = _tmp_mult * (dx31 * dx23 + dy31 * dy23);
        matrix[4] = _tmp_mult * (dx31 * dx31 + dy31 * dy31);
        matrix[5] = _tmp_mult * (dx31 * dx12 + dy31 * dy12);
        matrix[6] = _tmp_mult * (dx12 * dx23 + dy12 * dy23);
        matrix[7] = _tmp_mult * (dx12 * dx31 + dy12 * dy31);
        matrix[8] = _tmp_mult * (dx12 * dx12 + dy12 * dy12);
}

static inline void compute_local_f(uint *point_id, REAL *pointlist, REAL *matrix)
{
        REAL dx31 = pointlist[2 * point_id[2]] - pointlist[2 * point_id[0]];
        REAL dx12 = pointlist[2 * point_id[0]] - pointlist[2 * point_id[1]];
        REAL dy31 = pointlist[2 * point_id[2] + 1] - pointlist[2 * point_id[0] + 1];
        REAL dy12 = pointlist[2 * point_id[0] + 1] - pointlist[2 * point_id[1] + 1];

        REAL area = 0.5 * (dx31 * dy12 - dy31 * dx12);
        REAL c_diag = area / 6.0;
        REAL c_off = area / 12.0;

        matrix[0] = c_diag; matrix[1] = c_off;  matrix[2] = c_off;
        matrix[3] = c_off;  matrix[4] = c_diag; matrix[5] = c_off;
        matrix[6] = c_off;  matrix[7] = c_off;  matrix[8] = c_diag;
}

static inline REAL exact_solution(uint point_id, struct triangulateio *out)
{
    REAL x = out->pointlist[2 * point_id];
    REAL y = out->pointlist[2 * point_id + 1];
    return (1.0 + x * x + 2.0 * y * y);
}

static inline REAL boundary_condition(uint point_id, struct triangulateio *out)
{
    REAL x = out->pointlist[2 * point_id];
    REAL y = out->pointlist[2 * point_id + 1];
    return (1.0 + x * x + 2.0 * y * y);
}

ext: c


text:
#include <stdio.h>
#include "sweeny_dc.h"
#include <gsl/gsl_rng.h>
#include <math.h>
#include "../src/dyncon.h"
#include <ctype.h>

#define  MIN(a,b) a<=b ? a: b
#define MAX(a,b) a>=b? a:b
static char verbose = 0;
static char equilibration=1;
static char setup=0;

static __u32 DX;
static __u32 seed;
static __u32 cutoff;
static __u32 steps;
static double beta;
static double coupling;
static double q;
static double rnd_num;

static double v;
static __u32 j=0;
static __s8 dB = 0,dN=0;
static gsl_rng *r;
static Edge *cEdge;
static double rcWeights[4]; // array with precalculated mc weights
static double p_min_del,p_max_del, p_min_ins,p_max_ins;

static __u32 *num_bonds, *num_cluster, *size_giant;
static __u64 *sec_cs_moment,*four_cs_moment;

/******************************************************************************
 *****************************************************************************/
static void extract_observables(__u32 i) {
    static __u32 clust_cnt=0;
    static __u64 sum=0,sum_2=0;
    static __u32 maxc=0;
    static s_tree *c;
    num_bonds[i] = nte + te;
    c = levels[0]->head;
    while(c) {
        clust_cnt++;
        sum += pow(c->root->n,2);//*c->root->n;
        sum_2 += pow(c->root->n,4);
        if(c->root->n > maxc)
            maxc = c->root->n;
        c = c->next;
    }
    num_cluster[i] = clust_cnt;
    size_giant[i] = maxc;
    sec_cs_moment[i] = sum;
    four_cs_moment[i] = sum_2;
    sum=sum_2=0;
    maxc = 0;
    clust_cnt=0;
}

/******************************************************************************
 *****************************************************************************/

static inline void mcStep(void) {
    // get a randomly selected edge
    cEdge = &adjMatrix[gsl_rng_uniform_int(r,2*DX*DX)];

	rnd_num = gsl_rng_uniform(r);
 	if(ACTIVE_EDGE(cEdge)) { //edge is active hence delete it
	  dB = -1;
	  if(rnd_num < p_min_del) {
	    delete_dc(cEdge);
	  }
	  else {
	      if(rnd_num < p_max_del) {
		dN = delete_dc(cEdge);
		if(rnd_num >= rcWeights[dB == -1? dN : 2-dN])
		  insert_dc(cEdge);
	      }

	  }

	}
      else {
        // Insert cedge. If adjacent vertices are already connected, the new
        // edge will be a non-tree edge, i.e. dN = 0
        dB = 1;
        if(rnd_num < p_min_ins) {
	  insert_dc(cEdge);
	}
	else {
	  if(rnd_num < p_max_ins) {
	    dN = connected_dc_e(cEdge,1) ? 0: -1;
	    dB = 1;
	    if(rnd_num < rcWeights[dB == -1? dN : 2-dN])
		insert_dc(cEdge);
	  }

	}
      }

}
/******************************************************************************
 *****************************************************************************/
static inline void sweep() {
    for(j=0;j<2*DX*DX;j++)mcStep();
}
/******************************************************************************
 *****************************************************************************/
char init_sweeny_dc(double _q,unsigned int _l,double _beta,double _coupl,
        unsigned int _cutoff,unsigned _tslength,unsigned int rng_seed,
        void *ts_0,void *ts_1,void * ts_2,void *ts_3, void *ts_4) {

    q = _q;
    DX = _l;
    beta = _beta;
    coupling = _coupl;
    cutoff = _cutoff;
    steps = _tslength;
    seed = rng_seed;
    num_bonds = (__u32 *)ts_0;
    num_cluster = (__u32 *)ts_1;
    size_giant = (__u32 *)ts_2;
    sec_cs_moment = (__u64 *)ts_3;
    four_cs_moment = (__u64 *)ts_4;

    r = gsl_rng_alloc(gsl_rng_mt19937);
    gsl_rng_set(r,seed);
    v = exp(coupling*beta) - 1.; // >0
    rcWeights[0] = pow(v,-1); //db==-1,dN==0
    rcWeights[1] = rcWeights[0]*q; //db==-1,dN=1
    rcWeights[2] = v; //db==1,dN==0
    rcWeights[3] = v*pow(q,-1); //db==1,dN==-1
    p_min_del = MIN(rcWeights[0],rcWeights[1]);
    p_max_del = MAX(rcWeights[0],rcWeights[1]);
    p_min_ins = MIN(rcWeights[2],rcWeights[3]);
    p_max_ins = MAX(rcWeights[2],rcWeights[3]);
    return setup=(r && init_dc(DX));

}

/******************************************************************************
 *****************************************************************************/
void destroy_sweeny_dc(void) {
    if(setup) {
        gsl_rng_free(r);
        destroy_dc();
    }
    setup=0;
}
/******************************************************************************
 *****************************************************************************/
char simulate_sweeny_dc(void) {
    if(!setup) return 0;
    __u32 i;
    for(i=0;i<cutoff;i++)sweep();
    if(verbose)
        printf("Equilibration done!\n");
    equilibration=0;
    for(i=0;i<steps;i++){
        sweep();
        extract_observables(i);
    }
    return 1;
}
/******************************************************************************
 *****************************************************************************/

max_stars_repo_name: ernmeel/sweeny

max_issues_repo_name: ernmeel/sweeny


text:
/**
 * @file testmatmultc.c
 * @author Jeff Daily, PNNL jeff.daily@pnl.gov
 *
 * This is the C version of testmatmul.F.  The port was as direct as possible,
 * hopefully.
 */

#if HAVE_CONFIG_H
#   include "config.h"
#endif

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

#if defined(USE_ELEMENTAL)
#include <cblas.h>
#include "El.h"
int me, nproc;
ElGlobalArrays_d eldga;
#else
#include "mp3.h"
#include "ga.h"
#include "macdecls.h"
#include "xgemm.h"
#endif

void load_ga(int handle, double *f, int dim1, int dim2);
void verify_ga_dgemm(char xt1, char xt2, int num_m, int num_n, int num_k,
        double alpha, int g_a, int g_b, double beta, int g_c,
        double *tmpa, double *tmpb, double *tmpc);

#define dgemm_verify 1
#define nummax 1024
#define howmany 2
#define ntrans 4

#if defined(USE_ELEMENTAL)
#define MP_TIMER MPI_Wtime
void GA_Error (const char *str, int item)
{
    printf ("%s: %d\n", str, item);
    MPI_Abort(MPI_COMM_WORLD, -101);
}
#endif

/*
 * test ga_dgemm
 * Note: - change nummax for large arrays
 *       - turn off "dgemm_verify" for large arrays due to memory
 *         limitations, as dgemm_verify=1 for large arrays produces
 *         segfault, dumps core,or any crap.
 */
int main(int argc, char **argv)
{
    int num_m;
    int num_n;
    int num_k;
    int i;
    int ii;
    double *h0;
    int g_c;
    int g_b;
    int g_a;
    double a;
    double t1;
    double mf;
    double avg_t[ntrans];
    double avg_mf[ntrans];
    int itime;
    int ntimes;
    int nums_m[/*howmany*/] = {512,1024};
    int nums_n[/*howmany*/] = {512,1024};
    int nums_k[/*howmany*/] = {512,1024};
    char transa[/*ntrans*/] = "ntnt";
    char transb[/*ntrans*/] = "nntt";
    char ta;
    char tb;
    double *tmpa;
    double *tmpb;
    double *tmpc;
    int ndim;
    int dims[2];
#ifdef BLOCK_CYCLIC
    int block_size[2];
#endif

#if defined(USE_ELEMENTAL)
    // initialize Elemental (which will initialize MPI)
    ElInitialize( &argc, &argv );
    ElMPICommRank( MPI_COMM_WORLD, &me );
    ElMPICommSize( MPI_COMM_WORLD, &nproc );
    // instantiate el::global array
    ElGlobalArraysConstruct_d( &eldga );
    // initialize global arrays
    ElGlobalArraysInitialize_d( eldga );
#else
    MP_INIT(argc,argv);
    if (!MA_init(MT_DBL,1,20000000)) {
        GA_Error("failed: ma_init(MT_DBL,1,20000000)",10);
    }
    GA_INIT(argc,argv);
    me = GA_Nodeid();
#endif

    h0 = (double*)malloc(sizeof(double) * nummax*nummax);
    tmpa = (double*)malloc(sizeof(double) * nummax*nummax);
    tmpb = (double*)malloc(sizeof(double) * nummax*nummax);
    tmpc = (double*)malloc(sizeof(double) * nummax*nummax);

    ii = 0;
    for (i=0; i<nummax*nummax; i++) {
        ii = ii + 1;
        if (ii > nummax) {
            ii = 0;
        }
        h0[i] = ii;
    }

    /* Compute times assuming 500 mflops and 5 second target time */
    /* ntimes = max(3.0d0,5.0d0/(4.0d-9*num**3)); */
    ntimes = 5;

    for (ii=0; ii<howmany; ii++) {
        num_m = nums_m[ii];
        num_n = nums_n[ii];
        num_k = nums_k[ii];
        a = 0.5/(num_m*num_n);
        if (num_m > nummax || num_n > nummax || num_k > nummax) {
            GA_Error("Insufficient memory: check nummax", 1);
        }

#ifndef BLOCK_CYCLIC
        ndim = 2;

	/*
        dims[0] = num_m;
        dims[1] = num_n;
	*/
        dims[1] = num_m;
        dims[0] = num_n;

#if defined(USE_ELEMENTAL)
        ElGlobalArraysCreate_d( eldga, ndim, dims, "g_c", NULL, &g_c );
#else
        if (!((g_c = NGA_Create(MT_DBL,ndim,dims,"g_c",NULL)))) {
            GA_Error("failed: create g_c",20);
        }
#endif
	/*
        dims[0] = num_k;
        dims[1] = num_n;
	*/
        dims[1] = num_k;
        dims[0] = num_n;
#if defined(USE_ELEMENTAL)
        ElGlobalArraysCreate_d( eldga, ndim, dims, "g_b", NULL, &g_b );
#else
        if (!((g_b = NGA_Create(MT_DBL,ndim,dims,"g_b",NULL)))) {
            GA_Error("failed: create g_b",30);
        }
#endif
	/*
        dims[0] = num_m;
        dims[1] = num_k;
	*/
        dims[1] = num_m;
        dims[0] = num_k;
#if defined(USE_ELEMENTAL)
        ElGlobalArraysCreate_d( eldga, ndim, dims, "g_a", NULL, &g_a );
#else
        if (!((g_a = NGA_Create(MT_DBL,ndim,dims,"g_a",NULL)))) {
            GA_Error("failed: create g_a",40);
        }
#endif
#else
        ndim = 2;
        block_size[0] = 128;
        block_size[1] = 128;

        dims[0] = num_m;
        dims[1] = num_n;
        g_c = GA_Create_handle();
        GA_Set_data(g_c,ndim,dims,MT_DBL);
        GA_Set_array_name(g_c,"g_c");
        GA_Set_block_cyclic(g_c,block_size);
        if (!GA_Allocate(g_c)) {
            GA_Error("failed: create g_c",40);
        }

        dims[0] = num_k;
        dims[1] = num_n;
        g_b = GA_Create_handle();
        GA_Set_data(g_b,ndim,dims,MT_DBL);
        GA_Set_array_name(g_b,"g_b");
        GA_Set_block_cyclic(g_b,block_size);
        if (!ga_allocate(g_b)) {
            GA_Error("failed: create g_b",40);
        }

        dims[0] = num_m;
        dims[1] = num_k;
        g_a = GA_Create_handle();
        GA_Set_data(g_a,ndim,dims,MT_DBL);
        GA_Set_array_name(g_a,"g_a");
        GA_Set_block_cyclic(g_a,block_size);
        if (!ga_allocate(g_a)) {
            GA_Error('failed: create g_a',40);
        }
#endif

        /* Initialize matrices A and B */
        if (me == 0) {
            load_ga(g_a, h0, num_m, num_k);
            load_ga(g_b, h0, num_k, num_n);
        }
#if defined(USE_ELEMENTAL)
        double zero = 0.0;
        ElGlobalArraysFill_d( eldga, g_c, &zero );
	ElGlobalArraysSync_d( eldga );
#else
        GA_Zero(g_c);
        GA_Sync();
#endif
#if defined(USE_ELEMENTAL)
        if (me == 0) {
#else
        if (GA_Nodeid() == 0) {
#endif
            printf("\nMatrix Multiplication on C = A[%ld,%ld]xB[%ld,%ld]\n",
                    (long)num_m, (long)num_k, (long)num_k, (long)num_n);
            fflush(stdout);
        }

        for (i=0; i<ntrans; i++) {
            avg_t[i]  = 0.0;
            avg_mf[i] = 0.0;
        }

        for (itime=0; itime<ntimes; itime++) {
            for (i=0; i<ntrans; i++) {
#if defined(USE_ELEMENTAL)
	        ElGlobalArraysSync_d( eldga );
#else
                GA_Sync();
#endif
                ta = transa[i];
                tb = transb[i];
                t1 = MP_TIMER();
#if defined(USE_ELEMENTAL)
		ElGlobalArraysDgemm_d( eldga, ta, tb, num_m, num_n, num_k, 1.0, g_a, g_b, 0.0, g_c );
#else
                GA_Dgemm(ta,tb,num_m,num_n,num_k,1.0, g_a, g_b, 0.0, g_c);
#endif
                t1 = MP_TIMER() - t1;
#if defined(USE_ELEMENTAL)
                if (me == 0) {
#else
                if (GA_Nodeid() == 0) {
#endif
#if defined(USE_ELEMENTAL)
                    mf = 2e0*num_m*num_n*num_k/t1*1e-6/nproc;
#else
                    mf = 2e0*num_m*num_n*num_k/t1*1e-6/GA_Nnodes();
#endif
                    avg_t[i]  = avg_t[i]+t1;
                    avg_mf[i] = avg_mf[i] + mf;
                    printf("%15s%2d: %12.4f seconds %12.1f mflops/proc  %c %c\n",
                            "Run#", itime, t1, mf, ta, tb);
                    fflush(stdout);
                    if (dgemm_verify && itime == 0) {
                        /* recall the C API swaps the matrix order */
                        /* we swap it here for the Fortran-based verify */
                        verify_ga_dgemm(tb, ta, num_n, num_m, num_k, 1.0,
                                g_b, g_a, 0.0, g_c, tmpb, tmpa, tmpc);
                    }
                }
            }
        }
#if defined(USE_ELEMENTAL)
        if (me == 0) {
#else
        if (GA_Nodeid() == 0) {
#endif
            printf("\n");
            for (i=0; i<ntrans; i++) {
                printf("%17s: %12.4f seconds %12.1f mflops/proc  %c %c\n",
                        "Average", avg_t[i]/ntimes, avg_mf[i]/ntimes,
                        transa[i], transb[i]);
            }
            if(dgemm_verify) {
                printf("All GA_Dgemms are verified...O.K.\n");
            }
            fflush(stdout);
        }

        /*
           GA_Print(g_a);
           GA_Print(g_b);
           GA_Print(g_c);
           */
#if defined(USE_ELEMENTAL)
        ElGlobalArraysDestroy_d( eldga, g_a );
        ElGlobalArraysDestroy_d( eldga, g_b );
        ElGlobalArraysDestroy_d( eldga, g_c );
#else
        GA_Destroy(g_c);
        GA_Destroy(g_b);
        GA_Destroy(g_a);
#endif
    }

    /* ???
       format(a15, i2, ': ', e12.4, ' seconds ',f12.1,
       .     ' mflops/proc ', 3a2)
       */
#if defined(USE_ELEMENTAL)
    if (me == 0) {
#else
    if (GA_Nodeid() == 0) {
#endif
        printf("All tests successful\n");
    }

    free(h0);
    free(tmpa);
    free(tmpb);
    free(tmpc);
#if defined(USE_ELEMENTAL)
    // call el::global arrays destructor
    ElGlobalArraysTerminate_d( eldga );
    ElGlobalArraysDestruct_d( eldga );
    ElFinalize();
#else
    GA_Terminate();
    MP_FINALIZE();
#endif
    return 0;
}

/*
 * Verify for correctness. Process 0 computes BLAS dgemm
 * locally. For larger arrays, disbale this test as memory
 * might not be sufficient
 */
void verify_ga_dgemm(char xt1, char xt2, int num_m, int num_n, int num_k,
        double alpha, int g_a, int g_b, double beta, int g_c,
        double *tmpa, double *tmpb, double *tmpc)
{
    int i,j,type,ndim,dims[2],lo[2],hi[2];
    double abs_value;

    for (i=0; i<num_n; i++) {
        for (j=0; j<num_m; j++) {
            tmpc[j+i*num_m] = -1.0;
            tmpa[j+i*num_m] = -2.0;
        }
    }

#if defined(USE_ELEMENTAL)
    ElGlobalArraysInquire_d( eldga, g_a, &ndim, dims );
#else
    NGA_Inquire(g_a, &type, &ndim, dims);
#endif
    lo[0] = 0;
    lo[1] = 0;
    hi[0] = dims[0]-1;
    hi[1] = dims[1]-1;
#if defined(USE_ELEMENTAL)
    ElGlobalArraysGet_d( eldga, g_a, lo, hi, tmpa, &dims[1] );
#else
    NGA_Get(g_a, lo, hi, tmpa, &dims[1]);
#endif

#if defined(USE_ELEMENTAL)
    ElGlobalArraysInquire_d( eldga, g_a, &ndim, dims );
#else
    NGA_Inquire(g_a, &type, &ndim, dims);
#endif
    lo[0] = 0;
    lo[1] = 0;
    hi[0] = dims[0]-1;
    hi[1] = dims[1]-1;
#if defined(USE_ELEMENTAL)
    ElGlobalArraysGet_d( eldga, g_b, lo, hi, tmpb, &dims[1] );
#else
    NGA_Get(g_b, lo, hi, tmpb, &dims[1]);
#endif

    /* compute dgemm sequentially */
#if defined(USE_ELEMENTAL)
    cblas_dgemm ( CblasRowMajor, ( xt1 == 'n'? CblasNoTrans: CblasTrans ),
	    ( xt2 == 'n'? CblasNoTrans: CblasTrans ),
	    num_m /* M */, num_n /* N */, num_k /* K */,
	    alpha, tmpa, num_m, /* lda */
	    tmpb, num_k, /* ldb */ beta,
	    tmpc, num_m /* ldc */);
#else
    xb_dgemm(&xt1, &xt2, &num_m, &num_n, &num_k,
            &alpha, tmpa, &num_m,
            tmpb, &num_k, &beta,
            tmpc, &num_m);
#endif

    /* after computing c locally, verify it with the values in g_c */

#if defined(USE_ELEMENTAL)
    ElGlobalArraysInquire_d( eldga, g_a, &ndim, dims );
#else
    NGA_Inquire(g_a, &type, &ndim, dims);
#endif
    lo[0] = 0;
    lo[1] = 0;
    hi[0] = dims[0]-1;
    hi[1] = dims[1]-1;
#if defined(USE_ELEMENTAL)
    ElGlobalArraysGet_d( eldga, g_c, lo, hi, tmpa, &dims[1] );
#else
    NGA_Get(g_c, lo, hi, tmpa, &dims[1]);
#endif

    for (i=0; i<num_n; i++) {
        for (j=0; j<num_m; j++) {
            abs_value = fabs(tmpc[j+i*num_m]-tmpa[j+i*num_m]);
            if(abs_value > 1.0 || abs_value < -1.0) {
                printf("Values are = %f %f\n",
                        tmpc[j+i*num_m], tmpa[j+i*num_m]);
                printf("Values are = %f %f\n",
                        fabs(tmpc[j+i*num_m]-tmpa[j*i*num_m]), abs_value);
                fflush(stdout);
                GA_Error("verify ga_dgemm failed", 1);
            }
        }
    }
}

/**
 * called by process '0' (or your master process )
 */
void load_ga(int handle, double *f, int dim1, int dim2)
{
      int lo[2], hi[2];

      if (dim1 < 0 || dim2 < 0) {
          return;
      }

      lo[0] = 0;
      lo[1] = 0;
      hi[0] = dim1-1;
      hi[1] = dim2-1;
#if defined(USE_ELEMENTAL)
      ElGlobalArraysPut_d( eldga, handle, lo, hi, f, &dim1 );
#else
      NGA_Put(handle, lo, hi, f, &dim1);
#endif
}

/*
c
c-----------------------------------------------------------------------
c     must be called by all processors, if you need to fillup the
c     entire array
c
      subroutine load_ga_from_square(handle,num,f,ndim)
      implicit none
      integer handle, memhandle
      integer num,ndim
      real*8 f(ndim,ndim)
      integer ilo, ihi, jlo, jhi, nx, ny, ibuff
      integer ga_nodeid, i1, i2, i, j, ix, jx
#include "mafdecls.fh"

      call ga_distribution(handle, ga_nodeid(), ilo, ihi, jlo, jhi)

      if(ihi.le.0)return
      if(jhi.le.0)return

c     nx = ihi - ilo + 1
c     ny = jhi - jlo + 1

      do i = ilo,ihi,ndim
         do j = jlo,jhi,ndim
            call ga_put(handle,i,min(ihi,i+ndim),j,min(jhi,j+ndim),
     &                    f,ndim)
         enddo
      enddo

      return
      end
*/

/*
c
c-----------------------------------------------------------------------
c     must be called by all processors, if you need to fillup the
c     entire array
c
      subroutine load_ga_from_triangle(handle,f,ndim)
      implicit none
      integer handle, memhandle
      real*8 f(*)
      integer ndim
      integer ilo, ihi, jlo, jhi, nx, ny, ibuff
      integer ga_nodeid, i1, i2, i, j, ix, jx
#include "mafdecls.fh"

      call ga_distribution(handle, ga_nodeid(), ilo, ihi, jlo, jhi)

      if(ihi.le.0)return
      if(jhi.le.0)return

      nx = ihi - ilo + 1
      ny = jhi - jlo + 1

      if (.not.ma_alloc_get(MT_DBL,nx*ny,'flap',memhandle,ibuff)) then
         call ga_error('failed: allocate triangle',100)
      endif

      do i = 1,nx
         do j = 1,ny
            ix = i + ilo - 1
            jx = j + jlo - 1
            i1 = min(ix,jx)
            i2 = max(ix,jx)
            dbl_mb(ibuff + nx*(j-1) + (i-1) ) = f(i2*(i2-1)/2 + i1)
         enddo
      enddo

      call ga_put(handle,ilo,ihi,jlo,jhi,
     &              dbl_mb(ibuff),nx)

      if (.not.ma_free_heap(memhandle)) then
         call ga_error('failed: free triangle',100)
      endif

      return
      end
*/