diff --git a/blas/oblas_d_vsl_blas_cblas.v b/blas/oblas_d_vsl_blas_cblas.v index 29be3be06..777038f3d 100644 --- a/blas/oblas_d_vsl_blas_cblas.v +++ b/blas/oblas_d_vsl_blas_cblas.v @@ -221,6 +221,16 @@ pub fn dasum(n int, x []f64, incx int) f64 { return C.cblas_dasum(n, unsafe { &x[0] }, incx) } +@[inline] +pub fn scasum(n int, x voidptr, incx int) f32 { + return C.cblas_scasum(n, x, incx) +} + +@[inline] +pub fn dzasum(n int, x voidptr, incx int) f64 { + return C.cblas_dzasum(n, x, incx) +} + @[inline] pub fn ssum(n int, x []f32, incx int) f32 { return C.cblas_ssum(n, unsafe { &x[0] }, incx) @@ -241,6 +251,16 @@ pub fn dnrm2(n int, x []f64, incx int) f64 { return C.cblas_dnrm2(n, unsafe { &x[0] }, incx) } +@[inline] +pub fn scnrm2(n int, x voidptr, incx int) f32 { + return C.cblas_scnrm2(n, x, incx) +} + +@[inline] +pub fn dznrm2(n int, x voidptr, incx int) f64 { + return C.cblas_dznrm2(n, x, incx) +} + @[inline] pub fn isamax(n int, x []f32, incx int) int { return C.cblas_isamax(n, unsafe { &x[0] }, incx) @@ -251,16 +271,36 @@ pub fn idamax(n int, x []f64, incx int) int { return C.cblas_idamax(n, unsafe { &x[0] }, incx) } +@[inline] +pub fn icamax(n int, x voidptr, incx int) int { + return C.cblas_icamax(n, x, incx) +} + +@[inline] +pub fn izamax(n int, x voidptr, incx int) int { + return C.cblas_izamax(n, x, incx) +} + @[inline] pub fn isamin(n int, x []f32, incx int) int { return C.cblas_isamin(n, unsafe { &x[0] }, incx) } @[inline] -pub fn idamin(n int, x &f64, incx int) int { +pub fn idamin(n int, x []f64, incx int) int { return C.cblas_idamin(n, unsafe { &x[0] }, incx) } +@[inline] +pub fn icamin(n int, x voidptr, incx int) int { + return C.cblas_icamin(n, x, incx) +} + +@[inline] +pub fn izamin(n int, x voidptr, incx int) int { + return C.cblas_izamin(n, x, incx) +} + @[inline] pub fn ismax(n int, x []f32, incx int) int { return C.cblas_ismax(n, unsafe { &x[0] }, incx) @@ -271,6 +311,16 @@ pub fn idmax(n int, x []f64, incx int) int { return C.cblas_idmax(n, unsafe { &x[0] }, incx) } +@[inline] +pub fn icmax(n int, x voidptr, incx int) int { + return C.cblas_icmax(n, x, incx) +} + +@[inline] +pub fn izmax(n int, x voidptr, incx int) int { + return C.cblas_izmax(n, x, incx) +} + @[inline] pub fn ismin(n int, x []f32, incx int) int { return C.cblas_ismin(n, unsafe { &x[0] }, incx) @@ -281,6 +331,16 @@ pub fn idmin(n int, x []f64, incx int) int { return C.cblas_idmin(n, unsafe { &x[0] }, incx) } +@[inline] +pub fn icmin(n int, x voidptr, incx int) int { + return C.cblas_icmin(n, x, incx) +} + +@[inline] +pub fn izmin(n int, x voidptr, incx int) int { + return C.cblas_izmin(n, x, incx) +} + @[inline] pub fn saxpy(n int, alpha f32, x []f32, incx int, mut y []f32, incy int) { C.cblas_saxpy(n, alpha, unsafe { &x[0] }, incx, unsafe { &y[0] }, incy) @@ -292,15 +352,35 @@ pub fn daxpy(n int, alpha f64, x []f64, incx int, mut y []f64, incy int) { } @[inline] -pub fn scopy(n int, mut x []f32, incx int, mut y []f32, incy int) { +pub fn caxpy(n int, alpha voidptr, x voidptr, incx int, mut y voidptr, incy int) { + C.cblas_caxpy(n, alpha, x, incx, y, incy) +} + +@[inline] +pub fn zaxpy(n int, alpha voidptr, x voidptr, incx int, mut y voidptr, incy int) { + C.cblas_zaxpy(n, alpha, x, incx, y, incy) +} + +@[inline] +pub fn scopy(n int, x []f32, incx int, mut y []f32, incy int) { C.cblas_scopy(n, unsafe { &x[0] }, incx, unsafe { &y[0] }, incy) } @[inline] -pub fn dcopy(n int, mut x []f64, incx int, mut y []f64, incy int) { +pub fn dcopy(n int, x []f64, incx int, mut y []f64, incy int) { C.cblas_dcopy(n, unsafe { &x[0] }, incx, unsafe { &y[0] }, incy) } +@[inline] +pub fn ccopy(n int, x voidptr, incx int, mut y voidptr, incy int) { + C.cblas_ccopy(n, x, incx, y, incy) +} + +@[inline] +pub fn zcopy(n int, x voidptr, incx int, mut y voidptr, incy int) { + C.cblas_zcopy(n, x, incx, y, incy) +} + @[inline] pub fn sswap(n int, mut x []f32, incx int, mut y []f32, incy int) { C.cblas_sswap(n, unsafe { &x[0] }, incx, unsafe { &y[0] }, incy) @@ -311,6 +391,16 @@ pub fn dswap(n int, mut x []f64, incx int, mut y []f64, incy int) { C.cblas_dswap(n, unsafe { &x[0] }, incx, unsafe { &y[0] }, incy) } +@[inline] +pub fn cswap(n int, x voidptr, incx int, y voidptr, incy int) { + C.cblas_cswap(n, x, incx, y, incy) +} + +@[inline] +pub fn zswap(n int, x voidptr, incx int, y voidptr, incy int) { + C.cblas_zswap(n, x, incx, y, incy) +} + @[inline] pub fn srot(n int, mut x []f32, incx int, mut y []f32, incy int, c f32, s f32) { C.cblas_srot(n, unsafe { &x[0] }, incx, unsafe { &y[0] }, incy, c, s) @@ -332,22 +422,22 @@ pub fn drotg(a f64, b f64, c f64, s f64) { } @[inline] -pub fn srotm(n int, x []f32, incx int, y []f32, incy int, p []f32) { +pub fn srotm(n int, mut x []f32, incx int, mut y []f32, incy int, p []f32) { C.cblas_srotm(n, unsafe { &x[0] }, incx, unsafe { &y[0] }, incy, unsafe { &p[0] }) } @[inline] -pub fn drotm(n int, x []f64, incx int, y []f64, incy int, p []f64) { +pub fn drotm(n int, mut x []f64, incx int, mut y []f64, incy int, p []f64) { C.cblas_drotm(n, unsafe { &x[0] }, incx, unsafe { &y[0] }, incy, unsafe { &p[0] }) } @[inline] -pub fn srotmg(d1 f32, d2 f32, b1 f32, b2 f32, p []f32) { +pub fn srotmg(d1 f32, d2 f32, b1 f32, b2 f32, mut p []f32) { C.cblas_srotmg(&d1, &d2, &b1, b2, unsafe { &p[0] }) } @[inline] -pub fn drotmg(d1 f64, d2 f64, b1 f64, b2 f32, p []f64) { +pub fn drotmg(d1 f64, d2 f64, b1 f64, b2 f32, mut p []f64) { C.cblas_drotmg(&d1, &d2, &b1, b2, unsafe { &p[0] }) } @@ -361,6 +451,26 @@ pub fn dscal(n int, alpha f64, mut x []f64, incx int) { C.cblas_dscal(n, alpha, unsafe { &x[0] }, incx) } +@[inline] +pub fn cscal(n int, alpha voidptr, mut x voidptr, incx int) { + C.cblas_cscal(n, alpha, x, incx) +} + +@[inline] +pub fn zscal(n int, alpha voidptr, mut x voidptr, incx int) { + C.cblas_zscal(n, alpha, x, incx) +} + +@[inline] +pub fn csscal(n int, alpha f32, mut x voidptr, incx int) { + C.cblas_csscal(n, alpha, x, incx) +} + +@[inline] +pub fn zdscal(n int, alpha f64, mut x voidptr, incx int) { + C.cblas_zdscal(n, alpha, x, incx) +} + @[inline] pub fn sgemv(trans bool, m int, n int, alpha f32, a []f32, lda int, x []f32, incx int, beta f32, mut y []f32, incy int) { C.cblas_sgemv(.row_major, c_trans(trans), m, n, alpha, unsafe { &a[0] }, lda, unsafe { &x[0] }, @@ -373,6 +483,16 @@ pub fn dgemv(trans bool, m int, n int, alpha f64, a []f64, lda int, x []f64, inc incx, beta, unsafe { &y[0] }, incy) } +@[inline] +pub fn cgemv(trans bool, m int, n int, alpha voidptr, a voidptr, lda int, x voidptr, incx int, beta voidptr, mut y voidptr, incy int) { + C.cblas_cgemv(.row_major, c_trans(trans), m, n, alpha, a, lda, x, incx, beta, y, incy) +} + +@[inline] +pub fn zgemv(trans bool, m int, n int, alpha voidptr, a voidptr, lda int, x voidptr, incx int, beta voidptr, mut y voidptr, incy int) { + C.cblas_zgemv(.row_major, c_trans(trans), m, n, alpha, a, lda, x, incx, beta, y, incy) +} + @[inline] pub fn sger(m int, n int, alpha f32, x []f32, incx int, y []f32, incy int, mut a []f32, lda int) { C.cblas_sger(.row_major, m, n, alpha, unsafe { &x[0] }, incx, unsafe { &y[0] }, incy, @@ -386,29 +506,69 @@ pub fn dger(m int, n int, alpha f64, x []f64, incx int, y []f64, incy int, mut a } @[inline] -pub fn strsv(uplo bool, trans_a bool, diag Diagonal, n int, a []f32, lda int, mut x []f32, incx int) { - C.cblas_strsv(.row_major, c_uplo(uplo), c_trans(trans_a), diag, n, unsafe { &a[0] }, +pub fn cgeru(m int, n int, alpha voidptr, x voidptr, incx int, y voidptr, incy int, mut a voidptr, lda int) { + C.cblas_cgeru(.row_major, m, n, alpha, x, incx, y, incy, a, lda) +} + +@[inline] +pub fn cgerc(m int, n int, alpha voidptr, x voidptr, incx int, y voidptr, incy int, mut a voidptr, lda int) { + C.cblas_cgerc(.row_major, m, n, alpha, x, incx, y, incy, a, lda) +} + +@[inline] +pub fn zgeru(m int, n int, alpha voidptr, x voidptr, incx int, y voidptr, incy int, mut a voidptr, lda int) { + C.cblas_zgeru(.row_major, m, n, alpha, x, incx, y, incy, a, lda) +} + +@[inline] +pub fn zgerc(m int, n int, alpha voidptr, x voidptr, incx int, y voidptr, incy int, mut a voidptr, lda int) { + C.cblas_zgerc(.row_major, m, n, alpha, x, incx, y, incy, a, lda) +} + +@[inline] +pub fn strsv(uplo bool, trans bool, diag Diagonal, n int, a []f32, lda int, mut x []f32, incx int) { + C.cblas_strsv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, unsafe { &a[0] }, lda, unsafe { &x[0] }, incx) } @[inline] -pub fn dtrsv(uplo bool, trans_a bool, diag Diagonal, n int, a []f64, lda int, mut x []f64, incx int) { - C.cblas_dtrsv(.row_major, c_uplo(uplo), c_trans(trans_a), diag, n, unsafe { &a[0] }, +pub fn dtrsv(uplo bool, trans bool, diag Diagonal, n int, a []f64, lda int, mut x []f64, incx int) { + C.cblas_dtrsv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, unsafe { &a[0] }, lda, unsafe { &x[0] }, incx) } @[inline] -pub fn strmv(uplo bool, trans_a bool, diag Diagonal, n int, a []f32, lda int, mut x []f32, incx int) { - C.cblas_strmv(.row_major, c_uplo(uplo), c_trans(trans_a), diag, n, unsafe { &a[0] }, +pub fn ctrsv(uplo bool, trans bool, diag Diagonal, n int, a voidptr, lda int, mut x voidptr, incx int) { + C.cblas_ctrsv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, a, lda, x, incx) +} + +@[inline] +pub fn ztrsv(uplo bool, trans bool, diag Diagonal, n int, a voidptr, lda int, mut x voidptr, incx int) { + C.cblas_ztrsv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, a, lda, x, incx) +} + +@[inline] +pub fn strmv(uplo bool, trans bool, diag Diagonal, n int, a []f32, lda int, mut x []f32, incx int) { + C.cblas_strmv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, unsafe { &a[0] }, lda, unsafe { &x[0] }, incx) } @[inline] -pub fn dtrmv(uplo bool, trans_a bool, diag Diagonal, n int, a []f64, lda int, mut x []f64, incx int) { - C.cblas_dtrmv(.row_major, c_uplo(uplo), c_trans(trans_a), diag, n, unsafe { &a[0] }, +pub fn dtrmv(uplo bool, trans bool, diag Diagonal, n int, a []f64, lda int, mut x []f64, incx int) { + C.cblas_dtrmv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, unsafe { &a[0] }, lda, unsafe { &x[0] }, incx) } +@[inline] +pub fn ctrmv(uplo bool, trans bool, diag Diagonal, n int, a voidptr, lda int, mut x voidptr, incx int) { + C.cblas_ctrmv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, a, lda, x, incx) +} + +@[inline] +pub fn ztrmv(uplo bool, trans bool, diag Diagonal, n int, a voidptr, lda int, mut x voidptr, incx int) { + C.cblas_ztrmv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, a, lda, x, incx) +} + @[inline] pub fn ssyr(uplo bool, n int, alpha f32, x []f32, incx int, mut a []f32, lda int) { C.cblas_ssyr(.row_major, c_uplo(uplo), n, alpha, unsafe { &x[0] }, incx, unsafe { &a[0] }, @@ -421,6 +581,16 @@ pub fn dsyr(uplo bool, n int, alpha f64, x []f64, incx int, mut a []f64, lda int lda) } +@[inline] +pub fn cher(uplo bool, n int, alpha f32, x voidptr, incx int, mut a voidptr, lda int) { + C.cblas_cher(.row_major, c_uplo(uplo), n, alpha, x, incx, a, lda) +} + +@[inline] +pub fn zher(uplo bool, n int, alpha f64, x voidptr, incx int, mut a voidptr, lda int) { + C.cblas_zher(.row_major, c_uplo(uplo), n, alpha, x, incx, a, lda) +} + @[inline] pub fn ssyr2(uplo bool, n int, alpha f32, x []f32, incx int, y []f32, incy int, mut a []f32, lda int) { C.cblas_ssyr2(.row_major, c_uplo(uplo), n, alpha, unsafe { &x[0] }, incx, unsafe { &y[0] }, @@ -434,9 +604,451 @@ pub fn dsyr2(uplo bool, n int, alpha f64, x []f64, incx int, y []f64, incy int, } @[inline] -pub fn sgemm(trans_a bool, trans_b bool, m int, n int, k int, alpha f32, a []f32, lda int, b []f32, ldb int, beta f32, mut cc []f32, ldc int) { - C.cblas_sgemm(.row_major, c_trans(trans_a), c_trans(trans_b), m, n, k, alpha, unsafe { &a[0] }, - lda, unsafe { &b[0] }, ldb, beta, unsafe { &cc[0] }, ldc) +pub fn cher2(uplo bool, n int, alpha voidptr, x voidptr, incx int, y voidptr, incy int, mut a voidptr, lda int) { + C.cblas_cher2(.row_major, c_uplo(uplo), n, alpha, x, incx, y, incy, a, lda) +} + +@[inline] +pub fn zher2(uplo bool, n int, alpha voidptr, x voidptr, incx int, y voidptr, incy int, mut a voidptr, lda int) { + C.cblas_zher2(.row_major, c_uplo(uplo), n, alpha, x, incx, y, incy, a, lda) +} + +@[inline] +pub fn sgbmv(trans bool, m int, n int, kl int, ku int, alpha f32, a []f32, lda int, x []f32, incx int, beta f32, mut y []f32, incy int) { + C.cblas_sgbmv(.row_major, c_trans(trans), m, n, kl, ku, alpha, unsafe { &a[0] }, lda, + unsafe { &x[0] }, incx, beta, unsafe { &y[0] }, incy) +} + +@[inline] +pub fn dgbmv(trans bool, m int, n int, kl int, ku int, alpha f64, a []f64, lda int, x []f64, incx int, beta f64, mut y []f64, incy int) { + C.cblas_dgbmv(.row_major, c_trans(trans), m, n, kl, ku, alpha, unsafe { &a[0] }, lda, + unsafe { &x[0] }, incx, beta, unsafe { &y[0] }, incy) +} + +@[inline] +pub fn cgbmv(trans bool, m int, n int, kl int, ku int, alpha voidptr, a voidptr, lda int, x voidptr, incx int, beta voidptr, mut y voidptr, incy int) { + C.cblas_cgbmv(.row_major, c_trans(trans), m, n, kl, ku, alpha, a, lda, x, incx, beta, + y, incy) +} + +@[inline] +pub fn zgbmv(trans bool, m int, n int, kl int, ku int, alpha voidptr, a voidptr, lda int, x voidptr, incx int, beta voidptr, mut y voidptr, incy int) { + C.cblas_zgbmv(.row_major, c_trans(trans), m, n, kl, ku, alpha, a, lda, x, incx, beta, + y, incy) +} + +@[inline] +pub fn ssbmv(uplo bool, n int, k int, alpha f32, a []f32, lda int, x []f32, incx int, beta f32, mut y []f32, incy int) { + C.cblas_ssbmv(.row_major, c_uplo(uplo), n, k, alpha, unsafe { &a[0] }, lda, unsafe { &x[0] }, + incx, beta, unsafe { &y[0] }, incy) +} + +@[inline] +pub fn dsbmv(uplo bool, n int, k int, alpha f64, a []f64, lda int, x []f64, incx int, beta f64, mut y []f64, incy int) { + C.cblas_dsbmv(.row_major, c_uplo(uplo), n, k, alpha, unsafe { &a[0] }, lda, unsafe { &x[0] }, + incx, beta, unsafe { &y[0] }, incy) +} + +@[inline] +pub fn stbmv(uplo bool, trans bool, diag Diagonal, n int, k int, a []f32, lda int, mut x []f32, incx int) { + C.cblas_stbmv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, k, unsafe { &a[0] }, + lda, unsafe { &x[0] }, incx) +} + +@[inline] +pub fn dtbmv(uplo bool, trans bool, diag Diagonal, n int, k int, a []f64, lda int, mut x []f64, incx int) { + C.cblas_dtbmv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, k, unsafe { &a[0] }, + lda, unsafe { &x[0] }, incx) +} + +@[inline] +pub fn ctbmv(uplo bool, trans bool, diag Diagonal, n int, k int, a voidptr, lda int, mut x voidptr, incx int) { + C.cblas_ctbmv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, k, a, lda, x, incx) +} + +@[inline] +pub fn ztbmv(uplo bool, trans bool, diag Diagonal, n int, k int, a voidptr, lda int, mut x voidptr, incx int) { + C.cblas_ztbmv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, k, a, lda, x, incx) +} + +@[inline] +pub fn stbsv(uplo bool, trans bool, diag Diagonal, n int, k int, a []f32, lda int, mut x []f32, incx int) { + C.cblas_stbsv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, k, unsafe { &a[0] }, + lda, unsafe { &x[0] }, incx) +} + +@[inline] +pub fn dtbsv(uplo bool, trans bool, diag Diagonal, n int, k int, a []f64, lda int, mut x []f64, incx int) { + C.cblas_dtbsv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, k, unsafe { &a[0] }, + lda, unsafe { &x[0] }, incx) +} + +@[inline] +pub fn ctbsv(uplo bool, trans bool, diag Diagonal, n int, k int, a voidptr, lda int, mut x voidptr, incx int) { + C.cblas_ctbsv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, k, a, lda, x, incx) +} + +@[inline] +pub fn ztbsv(uplo bool, trans bool, diag Diagonal, n int, k int, a voidptr, lda int, mut x voidptr, incx int) { + C.cblas_ztbsv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, k, a, lda, x, incx) +} + +@[inline] +pub fn stpmv(uplo bool, trans bool, diag Diagonal, n int, ap []f32, mut x []f32, incx int) { + C.cblas_stpmv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, unsafe { &ap[0] }, + unsafe { &x[0] }, incx) +} + +@[inline] +pub fn dtpmv(uplo bool, trans bool, diag Diagonal, n int, ap []f64, mut x []f64, incx int) { + C.cblas_dtpmv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, unsafe { &ap[0] }, + unsafe { &x[0] }, incx) +} + +@[inline] +pub fn ctpmv(uplo bool, trans bool, diag Diagonal, n int, ap voidptr, mut x voidptr, incx int) { + C.cblas_ctpmv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, ap, x, incx) +} + +@[inline] +pub fn ztpmv(uplo bool, trans bool, diag Diagonal, n int, ap voidptr, mut x voidptr, incx int) { + C.cblas_ztpmv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, ap, x, incx) +} + +@[inline] +pub fn stpsv(uplo bool, trans bool, diag Diagonal, n int, ap []f32, mut x []f32, incx int) { + C.cblas_stpsv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, unsafe { &ap[0] }, + unsafe { &x[0] }, incx) +} + +@[inline] +pub fn dtpsv(uplo bool, trans bool, diag Diagonal, n int, ap []f64, mut x []f64, incx int) { + C.cblas_dtpsv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, unsafe { &ap[0] }, + unsafe { &x[0] }, incx) +} + +@[inline] +pub fn ctpsv(uplo bool, trans bool, diag Diagonal, n int, ap voidptr, mut x voidptr, incx int) { + C.cblas_ctpsv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, ap, x, incx) +} + +@[inline] +pub fn ztpsv(uplo bool, trans bool, diag Diagonal, n int, ap voidptr, mut x voidptr, incx int) { + C.cblas_ztpsv(.row_major, c_uplo(uplo), c_trans(trans), diag, n, ap, x, incx) +} + +@[inline] +pub fn ssymv(uplo bool, n int, alpha f32, a []f32, lda int, x []f32, incx int, beta f32, mut y []f32, incy int) { + C.cblas_ssymv(.row_major, c_uplo(uplo), n, alpha, unsafe { &a[0] }, lda, unsafe { &x[0] }, + incx, beta, unsafe { &y[0] }, incy) +} + +@[inline] +pub fn dsymv(uplo bool, n int, alpha f64, a []f64, lda int, x []f64, incx int, beta f64, mut y []f64, incy int) { + C.cblas_dsymv(.row_major, c_uplo(uplo), n, alpha, unsafe { &a[0] }, lda, unsafe { &x[0] }, + incx, beta, unsafe { &y[0] }, incy) +} + +@[inline] +pub fn chemv(uplo bool, n int, alpha voidptr, a voidptr, lda int, x voidptr, incx int, beta voidptr, mut y voidptr, incy int) { + C.cblas_chemv(.row_major, c_uplo(uplo), n, alpha, a, lda, x, incx, beta, y, incy) +} + +@[inline] +pub fn zhemv(uplo bool, n int, alpha voidptr, a voidptr, lda int, x voidptr, incx int, beta voidptr, mut y voidptr, incy int) { + C.cblas_zhemv(.row_major, c_uplo(uplo), n, alpha, a, lda, x, incx, beta, y, incy) +} + +@[inline] +pub fn sspmv(uplo bool, n int, alpha f32, ap []f32, x []f32, incx int, beta f32, mut y []f32, incy int) { + C.cblas_sspmv(.row_major, c_uplo(uplo), n, alpha, unsafe { &ap[0] }, unsafe { &x[0] }, + incx, beta, unsafe { &y[0] }, incy) +} + +@[inline] +pub fn dspmv(uplo bool, n int, alpha f64, ap []f64, x []f64, incx int, beta f64, mut y []f64, incy int) { + C.cblas_dspmv(.row_major, c_uplo(uplo), n, alpha, unsafe { &ap[0] }, unsafe { &x[0] }, + incx, beta, unsafe { &y[0] }, incy) +} + +@[inline] +pub fn sspr(uplo bool, n int, alpha f32, x []f32, incx int, mut ap []f32) { + C.cblas_sspr(.row_major, c_uplo(uplo), n, alpha, unsafe { &x[0] }, incx, unsafe { &ap[0] }) +} + +@[inline] +pub fn dspr(uplo bool, n int, alpha f64, x []f64, incx int, mut ap []f64) { + C.cblas_dspr(.row_major, c_uplo(uplo), n, alpha, unsafe { &x[0] }, incx, unsafe { &ap[0] }) +} + +@[inline] +pub fn chpr(uplo bool, n int, alpha f32, x voidptr, incx int, mut a voidptr) { + C.cblas_chpr(.row_major, c_uplo(uplo), n, alpha, x, incx, a) +} + +@[inline] +pub fn zhpr(uplo bool, n int, alpha f64, x voidptr, incx int, mut a voidptr) { + C.cblas_zhpr(.row_major, c_uplo(uplo), n, alpha, x, incx, a) +} + +@[inline] +pub fn sspr2(uplo bool, n int, alpha f32, x []f32, incx int, y []f32, incy int, mut a []f32) { + C.cblas_sspr2(.row_major, c_uplo(uplo), n, alpha, unsafe { &x[0] }, incx, unsafe { &y[0] }, + incy, unsafe { &a[0] }) +} + +@[inline] +pub fn dspr2(uplo bool, n int, alpha f64, x []f64, incx int, y []f64, incy int, mut a []f64) { + C.cblas_dspr2(.row_major, c_uplo(uplo), n, alpha, unsafe { &x[0] }, incx, unsafe { &y[0] }, + incy, unsafe { &a[0] }) +} + +@[inline] +pub fn chpr2(uplo bool, n int, alpha voidptr, x voidptr, incx int, y voidptr, incy int, mut ap voidptr) { + C.cblas_chpr2(.row_major, c_uplo(uplo), n, alpha, x, incx, y, incy, ap) +} + +@[inline] +pub fn zhpr2(uplo bool, n int, alpha voidptr, x voidptr, incx int, y voidptr, incy int, mut ap voidptr) { + C.cblas_zhpr2(.row_major, c_uplo(uplo), n, alpha, x, incx, y, incy, ap) +} + +@[inline] +pub fn chbmv(uplo bool, n int, k int, alpha voidptr, a voidptr, lda int, x voidptr, incx int, beta voidptr, mut y voidptr, incy int) { + C.cblas_chbmv(.row_major, c_uplo(uplo), n, k, alpha, a, lda, x, incx, beta, y, incy) +} + +@[inline] +pub fn zhbmv(uplo bool, n int, k int, alpha voidptr, a voidptr, lda int, x voidptr, incx int, beta voidptr, mut y voidptr, incy int) { + C.cblas_zhbmv(.row_major, c_uplo(uplo), n, k, alpha, a, lda, x, incx, beta, y, incy) +} + +@[inline] +pub fn chpmv(uplo bool, n int, alpha voidptr, ap voidptr, x voidptr, incx int, beta voidptr, mut y voidptr, incy int) { + C.cblas_chpmv(.row_major, c_uplo(uplo), n, alpha, ap, x, incx, beta, y, incy) +} + +@[inline] +pub fn zhpmv(uplo bool, n int, alpha voidptr, ap voidptr, x voidptr, incx int, beta voidptr, mut y voidptr, incy int) { + C.cblas_zhpmv(.row_major, c_uplo(uplo), n, alpha, ap, x, incx, beta, y, incy) +} + +@[inline] +pub fn ssyrk(uplo bool, trans bool, n int, k int, alpha f32, a []f32, lda int, beta f32, mut c []f32, ldc int) { + C.cblas_ssyrk(.row_major, c_uplo(uplo), c_trans(trans), n, k, alpha, unsafe { &a[0] }, + lda, beta, unsafe { &c[0] }, ldc) +} + +@[inline] +pub fn dsyrk(uplo bool, trans bool, n int, k int, alpha f64, a []f64, lda int, beta f64, mut c []f64, ldc int) { + C.cblas_dsyrk(.row_major, c_uplo(uplo), c_trans(trans), n, k, alpha, unsafe { &a[0] }, + lda, beta, unsafe { &c[0] }, ldc) +} + +@[inline] +pub fn csyrk(uplo bool, trans bool, n int, k int, alpha voidptr, a voidptr, lda int, beta voidptr, mut c voidptr, ldc int) { + C.cblas_csyrk(.row_major, c_uplo(uplo), c_trans(trans), n, k, alpha, a, lda, beta, + c, ldc) +} + +@[inline] +pub fn zsyrk(uplo bool, trans bool, n int, k int, alpha voidptr, a voidptr, lda int, beta voidptr, mut c voidptr, ldc int) { + C.cblas_zsyrk(.row_major, c_uplo(uplo), c_trans(trans), n, k, alpha, a, lda, beta, + c, ldc) +} + +@[inline] +pub fn ssyr2k(uplo bool, trans bool, n int, k int, alpha f32, a []f32, lda int, b []f32, ldb int, beta f32, mut c []f32, ldc int) { + C.cblas_ssyr2k(.row_major, c_uplo(uplo), c_trans(trans), n, k, alpha, unsafe { &a[0] }, + lda, unsafe { &b[0] }, ldb, beta, unsafe { &c[0] }, ldc) +} + +@[inline] +pub fn dsyr2k(uplo bool, trans bool, n int, k int, alpha f64, a []f64, lda int, b []f64, ldb int, beta f64, mut c []f64, ldc int) { + C.cblas_dsyr2k(.row_major, c_uplo(uplo), c_trans(trans), n, k, alpha, unsafe { &a[0] }, + lda, unsafe { &b[0] }, ldb, beta, unsafe { &c[0] }, ldc) +} + +@[inline] +pub fn csyr2k(uplo bool, trans bool, n int, k int, alpha voidptr, a voidptr, lda int, b voidptr, ldb int, beta voidptr, mut c voidptr, ldc int) { + C.cblas_csyr2k(.row_major, c_uplo(uplo), c_trans(trans), n, k, alpha, a, lda, b, ldb, + beta, c, ldc) +} + +@[inline] +pub fn zsyr2k(uplo bool, trans bool, n int, k int, alpha voidptr, a voidptr, lda int, b voidptr, ldb int, beta voidptr, mut c voidptr, ldc int) { + C.cblas_zsyr2k(.row_major, c_uplo(uplo), c_trans(trans), n, k, alpha, a, lda, b, ldb, + beta, c, ldc) +} + +@[inline] +pub fn strmm(side Side, uplo bool, trans bool, diag Diagonal, m int, n int, alpha f32, a []f32, lda int, mut b []f32, ldb int) { + C.cblas_strmm(.row_major, side, c_uplo(uplo), c_trans(trans), diag, m, n, alpha, unsafe { &a[0] }, + lda, unsafe { &b[0] }, ldb) +} + +@[inline] +pub fn dtrmm(side Side, uplo bool, trans bool, diag Diagonal, m int, n int, alpha f64, a []f64, lda int, mut b []f64, ldb int) { + C.cblas_dtrmm(.row_major, side, c_uplo(uplo), c_trans(trans), diag, m, n, alpha, unsafe { &a[0] }, + lda, unsafe { &b[0] }, ldb) +} + +@[inline] +pub fn ctrmm(side Side, uplo bool, trans bool, diag Diagonal, m int, n int, alpha voidptr, a voidptr, lda int, mut b voidptr, ldb int) { + C.cblas_ctrmm(.row_major, side, c_uplo(uplo), c_trans(trans), diag, m, n, alpha, a, + lda, b, ldb) +} + +@[inline] +pub fn ztrmm(side Side, uplo bool, trans bool, diag Diagonal, m int, n int, alpha voidptr, a voidptr, lda int, mut b voidptr, ldb int) { + C.cblas_ztrmm(.row_major, side, c_uplo(uplo), c_trans(trans), diag, m, n, alpha, a, + lda, b, ldb) +} + +@[inline] +pub fn strsm(side Side, uplo bool, trans bool, diag Diagonal, m int, n int, alpha f32, a []f32, lda int, mut b []f32, ldb int) { + C.cblas_strsm(.row_major, side, c_uplo(uplo), c_trans(trans), diag, m, n, alpha, unsafe { &a[0] }, + lda, unsafe { &b[0] }, ldb) +} + +@[inline] +pub fn dtrsm(side Side, uplo bool, trans bool, diag Diagonal, m int, n int, alpha f64, a []f64, lda int, mut b []f64, ldb int) { + C.cblas_dtrsm(.row_major, side, c_uplo(uplo), c_trans(trans), diag, m, n, alpha, unsafe { &a[0] }, + lda, unsafe { &b[0] }, ldb) +} + +@[inline] +pub fn ctrsm(side Side, uplo bool, trans bool, diag Diagonal, m int, n int, alpha voidptr, a voidptr, lda int, mut b voidptr, ldb int) { + C.cblas_ctrsm(.row_major, side, c_uplo(uplo), c_trans(trans), diag, m, n, alpha, a, + lda, b, ldb) +} + +@[inline] +pub fn ztrsm(side Side, uplo bool, trans bool, diag Diagonal, m int, n int, alpha voidptr, a voidptr, lda int, mut b voidptr, ldb int) { + C.cblas_ztrsm(.row_major, side, c_uplo(uplo), c_trans(trans), diag, m, n, alpha, a, + lda, b, ldb) +} + +@[inline] +pub fn chemm(side Side, uplo bool, m int, n int, alpha voidptr, a voidptr, lda int, b voidptr, ldb int, beta voidptr, mut c voidptr, ldc int) { + C.cblas_chemm(.row_major, side, c_uplo(uplo), m, n, alpha, a, lda, b, ldb, beta, c, + ldc) +} + +@[inline] +pub fn zhemm(side Side, uplo bool, m int, n int, alpha voidptr, a voidptr, lda int, b voidptr, ldb int, beta voidptr, mut c voidptr, ldc int) { + C.cblas_zhemm(.row_major, side, c_uplo(uplo), m, n, alpha, a, lda, b, ldb, beta, c, + ldc) +} + +@[inline] +pub fn cherk(uplo bool, trans bool, n int, k int, alpha f32, a voidptr, lda int, beta f32, mut c voidptr, ldc int) { + C.cblas_cherk(.row_major, c_uplo(uplo), c_trans(trans), n, k, alpha, a, lda, beta, + c, ldc) +} + +@[inline] +pub fn zherk(uplo bool, trans bool, n int, k int, alpha f64, a voidptr, lda int, beta f64, mut c voidptr, ldc int) { + C.cblas_zherk(.row_major, c_uplo(uplo), c_trans(trans), n, k, alpha, a, lda, beta, + c, ldc) +} + +@[inline] +pub fn cher2k(uplo bool, trans bool, n int, k int, alpha voidptr, a voidptr, lda int, b voidptr, ldb int, beta f32, mut c voidptr, ldc int) { + C.cblas_cher2k(.row_major, c_uplo(uplo), c_trans(trans), n, k, alpha, a, lda, b, ldb, + beta, c, ldc) +} + +@[inline] +pub fn zher2k(uplo bool, trans bool, n int, k int, alpha voidptr, a voidptr, lda int, b voidptr, ldb int, beta f64, mut c voidptr, ldc int) { + C.cblas_zher2k(.row_major, c_uplo(uplo), c_trans(trans), n, k, alpha, a, lda, b, ldb, + beta, c, ldc) +} + +@[inline] +pub fn saxpby(n int, alpha f32, x []f32, incx int, beta f32, mut y []f32, incy int) { + C.cblas_saxpby(n, alpha, unsafe { &x[0] }, incx, beta, unsafe { &y[0] }, incy) +} + +@[inline] +pub fn daxpby(n int, alpha f64, x []f64, incx int, beta f64, mut y []f64, incy int) { + C.cblas_daxpby(n, alpha, unsafe { &x[0] }, incx, beta, unsafe { &y[0] }, incy) +} + +@[inline] +pub fn caxpby(n int, alpha voidptr, x voidptr, incx int, beta voidptr, mut y voidptr, incy int) { + C.cblas_caxpby(n, alpha, x, incx, beta, y, incy) +} + +@[inline] +pub fn zaxpby(n int, alpha voidptr, x voidptr, incx int, beta voidptr, mut y voidptr, incy int) { + C.cblas_zaxpby(n, alpha, x, incx, beta, y, incy) +} + +@[inline] +pub fn somatcopy(order MemoryLayout, trans bool, rows int, cols int, alpha f32, a []f32, lda int, mut b []f32, ldb int) { + C.cblas_somatcopy(order, c_trans(trans), rows, cols, alpha, unsafe { &a[0] }, lda, + unsafe { &b[0] }, ldb) +} + +@[inline] +pub fn domatcopy(order MemoryLayout, trans bool, rows int, cols int, alpha f64, a []f64, lda int, mut b []f64, ldb int) { + C.cblas_domatcopy(order, c_trans(trans), rows, cols, alpha, unsafe { &a[0] }, lda, + unsafe { &b[0] }, ldb) +} + +@[inline] +pub fn comatcopy(order MemoryLayout, trans bool, rows int, cols int, alpha &f32, a &f32, lda int, mut b &f32, ldb int) { + C.cblas_comatcopy(order, c_trans(trans), rows, cols, alpha, a, lda, b, ldb) +} + +@[inline] +pub fn zomatcopy(order MemoryLayout, trans bool, rows int, cols int, alpha &f64, a &f64, lda int, mut b &f64, ldb int) { + C.cblas_zomatcopy(order, c_trans(trans), rows, cols, alpha, a, lda, b, ldb) +} + +@[inline] +pub fn simatcopy(order MemoryLayout, trans bool, rows int, cols int, alpha f32, mut a []f32, lda int, ldb int) { + C.cblas_simatcopy(order, c_trans(trans), rows, cols, alpha, unsafe { &a[0] }, lda, + ldb) +} + +@[inline] +pub fn dimatcopy(order MemoryLayout, trans bool, rows int, cols int, alpha f64, mut a []f64, lda int, ldb int) { + C.cblas_dimatcopy(order, c_trans(trans), rows, cols, alpha, unsafe { &a[0] }, lda, + ldb) +} + +@[inline] +pub fn cimatcopy(order MemoryLayout, trans bool, rows int, cols int, alpha &f32, mut a &f32, lda int, ldb int) { + C.cblas_cimatcopy(order, c_trans(trans), rows, cols, alpha, a, lda, ldb) +} + +@[inline] +pub fn zimatcopy(order MemoryLayout, trans bool, rows int, cols int, alpha &f64, mut a &f64, lda int, ldb int) { + C.cblas_zimatcopy(order, c_trans(trans), rows, cols, alpha, a, lda, ldb) +} + +@[inline] +pub fn sgeadd(order MemoryLayout, rows int, cols int, alpha f32, a []f32, lda int, beta f32, mut c []f32, ldc int) { + C.cblas_sgeadd(order, rows, cols, alpha, unsafe { &a[0] }, lda, beta, unsafe { &c[0] }, + ldc) +} + +@[inline] +pub fn dgeadd(order MemoryLayout, rows int, cols int, alpha f64, a []f64, lda int, beta f64, mut c []f64, ldc int) { + C.cblas_dgeadd(order, rows, cols, alpha, unsafe { &a[0] }, lda, beta, unsafe { &c[0] }, + ldc) +} + +@[inline] +pub fn cgeadd(order MemoryLayout, rows int, cols int, alpha &f32, a &f32, lda int, beta &f32, mut c &f32, ldc int) { + C.cblas_cgeadd(order, rows, cols, alpha, a, lda, beta, c, ldc) +} + +@[inline] +pub fn zgeadd(order MemoryLayout, rows int, cols int, alpha &f64, a &f64, lda int, beta &f64, mut c &f64, ldc int) { + C.cblas_zgeadd(order, rows, cols, alpha, a, lda, beta, c, ldc) } @[inline] diff --git a/la/matrix_ops.v b/la/matrix_ops.v index 74dc6ea28..a792aeede 100644 --- a/la/matrix_ops.v +++ b/la/matrix_ops.v @@ -90,7 +90,7 @@ pub fn matrix_svd(mut s []f64, mut u Matrix[f64], mut vt Matrix[f64], mut a Matr if copy_a { acpy = a.clone() } - lapack.dgesvd(&char('A'.str), &char('A'.str), a.m, a.n, acpy.data, 1, s, u.data, a.m, + lapack.dgesvd(.svd_all, .svd_all, a.m, a.n, mut acpy.data, 1, s, mut u.data, a.m, mut vt.data, a.n, superb) } diff --git a/lapack/conversions.v b/lapack/conversions.v new file mode 100644 index 000000000..fcb1cb239 --- /dev/null +++ b/lapack/conversions.v @@ -0,0 +1,70 @@ +module lapack + +import vsl.lapack.lapack64 + +// Direct specifies the direction of the multiplication for the Householder matrix. +pub type Direct = lapack64.Direct + +// Sort is the sorting order. +pub type Sort = lapack64.Sort + +// StoreV indicates the storage direction of elementary reflectors. +pub type StoreV = lapack64.StoreV + +// MatrixNorm represents the kind of matrix norm to compute. +pub type MatrixNorm = lapack64.MatrixNorm + +// MatrixType represents the kind of matrix represented in the data. +pub type MatrixType = lapack64.MatrixType + +// Pivot specifies the pivot type for plane rotations. +pub type Pivot = lapack64.Pivot + +// ApplyOrtho specifies which orthogonal matrix is applied in Dormbr. +pub type ApplyOrtho = lapack64.ApplyOrtho + +// GenOrtho specifies which orthogonal matrix is generated in Dorgbr. +pub type GenOrtho = lapack64.GenOrtho + +// SVDJob specifies the singular vector computation type for SVD. +pub type SVDJob = lapack64.SVDJob + +// GSVDJob specifies the singular vector computation type for Generalized SVD. +pub type GSVDJob = lapack64.GSVDJob + +// EVComp specifies how eigenvectors are computed in Dsteqr. +pub type EVComp = lapack64.EVComp + +// EVJob specifies whether eigenvectors are computed in Dsyev. +pub type EVJob = lapack64.EVJob + +// LeftEVJob specifies whether left eigenvectors are computed in Dgeev. +pub type LeftEVJob = lapack64.LeftEVJob + +// RightEVJob specifies whether right eigenvectors are computed in Dgeev. +pub type RightEVJob = lapack64.RightEVJob + +// BalanceJob specifies matrix balancing operation. +pub type BalanceJob = lapack64.BalanceJob + +// SchurJob specifies whether the Schur form is computed in Dhseqr. +pub type SchurJob = lapack64.SchurJob + +// SchurComp specifies whether and how the Schur vectors are computed in Dhseqr. +pub type SchurComp = lapack64.SchurComp + +// UpdateSchurComp specifies whether the matrix of Schur vectors is updated in Dtrexc. +pub type UpdateSchurComp = lapack64.UpdateSchurComp + +// EVSide specifies what eigenvectors are computed in Dtrevc3. +pub type EVSide = lapack64.EVSide + +// EVHowMany specifies which eigenvectors are computed in Dtrevc3 and how. +pub type EVHowMany = lapack64.EVHowMany + +// MaximizeNormXJob specifies the heuristic method for computing a contribution to +// the reciprocal Dif-estimate in Dlatdf. +pub type MaximizeNormXJob = lapack64.MaximizeNormXJob + +// OrthoComp specifies whether and how the orthogonal matrix is computed in Dgghrd. +pub type OrthoComp = lapack64.OrthoComp diff --git a/lapack/lapack64/conversions.v b/lapack/lapack64/conversions.v new file mode 100644 index 000000000..5052633f9 --- /dev/null +++ b/lapack/lapack64/conversions.v @@ -0,0 +1,199 @@ +module lapack64 + +// Direct specifies the direction of the multiplication for the Householder matrix. +pub enum Direct as u8 { + // Reflectors are right-multiplied, H_0 * H_1 * ... * H_{k-1}. + forward = u8(`F`) + // Reflectors are left-multiplied, H_{k-1} * ... * H_1 * H_0. + backward = u8(`B`) +} + +// Sort is the sorting order. +pub enum Sort as u8 { + sort_increasing = u8(`I`) + sort_decreasing = u8(`D`) +} + +// StoreV indicates the storage direction of elementary reflectors. +pub enum StoreV as u8 { + // Reflector stored in a column of the matrix. + column_wise = u8(`C`) + // Reflector stored in a row of the matrix. + row_wise = u8(`R`) +} + +// MatrixNorm represents the kind of matrix norm to compute. +pub enum MatrixNorm as u8 { + // max(abs(A(i,j))) + max_abs = u8(`M`) + // Maximum absolute column sum (one norm) + max_column_sum = u8(`O`) + // Maximum absolute row sum (infinity norm) + max_row_sum = u8(`I`) + // Frobenius norm (sqrt of sum of squares) + frobenius = u8(`F`) +} + +// MatrixType represents the kind of matrix represented in the data. +pub enum MatrixType as u8 { + // A general dense matrix. + general = u8(`G`) + // An upper triangular matrix. + upper_tri = u8(`U`) + // A lower triangular matrix. + lower_tri = u8(`L`) +} + +// Pivot specifies the pivot type for plane rotations. +pub enum Pivot as u8 { + variable = u8(`V`) + top = u8(`T`) + bottom = u8(`B`) +} + +// ApplyOrtho specifies which orthogonal matrix is applied in Dormbr. +pub enum ApplyOrtho as u8 { + // Apply P or Pᵀ. + apply_p = u8(`P`) + // Apply Q or Qᵀ. + apply_q = u8(`Q`) +} + +// GenOrtho specifies which orthogonal matrix is generated in Dorgbr. +pub enum GenOrtho as u8 { + // Generate Pᵀ. + generate_pt = u8(`P`) + // Generate Q. + generate_q = u8(`Q`) +} + +// SVDJob specifies the singular vector computation type for SVD. +pub enum SVDJob as u8 { + // Compute all columns of the orthogonal matrix U or V. + svd_all = u8(`A`) + // Compute the singular vectors and store them in the orthogonal matrix U or V. + svd_store = u8(`S`) + // Compute the singular vectors and overwrite them on the input matrix A. + svd_overwrite = u8(`O`) + // Do not compute singular vectors. + svd_none = u8(`N`) +} + +// GSVDJob specifies the singular vector computation type for Generalized SVD. +pub enum GSVDJob as u8 { + // Compute orthogonal matrix U. + gsvd_u = u8(`U`) + // Compute orthogonal matrix V. + gsvd_v = u8(`V`) + // Compute orthogonal matrix Q. + gsvd_q = u8(`Q`) + // Use unit-initialized matrix. + gsvd_unit = u8(`I`) + // Do not compute orthogonal matrix. + gsvd_none = u8(`N`) +} + +// EVComp specifies how eigenvectors are computed in Dsteqr. +pub enum EVComp as u8 { + // Compute eigenvectors of the original symmetric matrix. + ev_orig = u8(`V`) + // Compute eigenvectors of the tridiagonal matrix. + ev_tridiag = u8(`I`) + // Do not compute eigenvectors. + ev_comp_none = u8(`N`) +} + +// EVJob specifies whether eigenvectors are computed in Dsyev. +pub enum EVJob as u8 { + // Compute eigenvectors. + ev_compute = u8(`V`) + // Do not compute eigenvectors. + ev_none = u8(`N`) +} + +// LeftEVJob specifies whether left eigenvectors are computed in Dgeev. +pub enum LeftEVJob as u8 { + // Compute left eigenvectors. + left_ev_compute = u8(`V`) + // Do not compute left eigenvectors. + left_ev_none = u8(`N`) +} + +// RightEVJob specifies whether right eigenvectors are computed in Dgeev. +pub enum RightEVJob as u8 { + // Compute right eigenvectors. + right_ev_compute = u8(`V`) + // Do not compute right eigenvectors. + right_ev_none = u8(`N`) +} + +// BalanceJob specifies matrix balancing operation. +pub enum BalanceJob as u8 { + permute = u8(`P`) + scale = u8(`S`) + permute_scale = u8(`B`) + balance_none = u8(`N`) +} + +// SchurJob specifies whether the Schur form is computed in Dhseqr. +pub enum SchurJob as u8 { + eigenvalues_only = u8(`E`) + eigenvalues_and_schur = u8(`S`) +} + +// SchurComp specifies whether and how the Schur vectors are computed in Dhseqr. +pub enum SchurComp as u8 { + // Compute Schur vectors of the original matrix. + schur_orig = u8(`V`) + // Compute Schur vectors of the upper Hessenberg matrix. + schur_hess = u8(`I`) + // Do not compute Schur vectors. + schur_none = u8(`N`) +} + +// UpdateSchurComp specifies whether the matrix of Schur vectors is updated in Dtrexc. +pub enum UpdateSchurComp as u8 { + // Update the matrix of Schur vectors. + update_schur = u8(`V`) + // Do not update the matrix of Schur vectors. + update_schur_none = u8(`N`) +} + +// EVSide specifies what eigenvectors are computed in Dtrevc3. +pub enum EVSide as u8 { + // Compute only right eigenvectors. + ev_right = u8(`R`) + // Compute only left eigenvectors. + ev_left = u8(`L`) + // Compute both right and left eigenvectors. + ev_both = u8(`B`) +} + +// EVHowMany specifies which eigenvectors are computed in Dtrevc3 and how. +pub enum EVHowMany as u8 { + // Compute all right and/or left eigenvectors. + ev_all = u8(`A`) + // Compute all right and/or left eigenvectors multiplied by an input matrix. + ev_all_mul_q = u8(`B`) + // Compute selected right and/or left eigenvectors. + ev_selected = u8(`S`) +} + +// MaximizeNormXJob specifies the heuristic method for computing a contribution to +// the reciprocal Dif-estimate in Dlatdf. +pub enum MaximizeNormXJob as u8 { + // Solve Z*x=h-f where h is a vector of ±1. + local_look_ahead = 0 + // Compute an approximate null-vector e of Z, normalize e and solve Z*x=±e-f. + normalized_null_vector = 2 +} + +// OrthoComp specifies whether and how the orthogonal matrix is computed in Dgghrd. +pub enum OrthoComp as u8 { + // Do not compute the orthogonal matrix. + ortho_none = u8(`N`) + // The orthogonal matrix is formed explicitly and returned in the argument. + ortho_explicit = u8(`I`) + // The orthogonal matrix is post-multiplied into the matrix stored in the argument on entry. + ortho_postmul = u8(`V`) +} diff --git a/lapack/lapack64/dgebal.v b/lapack/lapack64/dgebal.v new file mode 100644 index 000000000..ee7746a38 --- /dev/null +++ b/lapack/lapack64/dgebal.v @@ -0,0 +1,33 @@ +module lapack64 + +import math +import vsl.blas + +// dgebal balances a general real matrix A. +pub fn dgebal(job BalanceJob, n int, mut a []f64, lda int, scale []f64) int { + if n == 0 { + return 0 + } + + mut info := 0 + if job != .balance_none && job != .permute && job != .scale && job != .permute_scale { + info = -1 + } else if n < 0 { + info = -2 + } else if lda < math.max(1, n) { + info = -4 + } + + if info != 0 { + return info + } + + // Quick return if possible + if n == 0 { + return 0 + } + + // Placeholder for the actual LAPACK function calls + // Example: info = dgebal(job, n, a, lda, scale) + return info +} diff --git a/lapack/lapack64/dgeev.v b/lapack/lapack64/dgeev.v new file mode 100644 index 000000000..38799070c --- /dev/null +++ b/lapack/lapack64/dgeev.v @@ -0,0 +1,39 @@ +module lapack64 + +import math +import vsl.blas + +// dgeev computes the eigenvalues and, optionally, the left and/or right eigenvectors for a real nonsymmetric matrix A. +pub fn dgeev(jobvl LeftEVJob, jobvr LeftEVJob, n int, mut a []f64, lda int, wr []f64, wi []f64, mut vl []f64, ldvl int, mut vr []f64, ldvr int) int { + if n == 0 { + return 0 + } + + mut info := 0 + if jobvl != .left_ev_none && jobvl != .left_ev_compute { + info = -1 + } else if jobvr != .left_ev_none && jobvr != .left_ev_compute { + info = -2 + } else if n < 0 { + info = -3 + } else if lda < math.max(1, n) { + info = -5 + } else if ldvl < 1 || (jobvl == .left_ev_compute && ldvl < n) { + info = -8 + } else if ldvr < 1 || (jobvr == .left_ev_compute && ldvr < n) { + info = -10 + } + + if info != 0 { + return info + } + + // Quick return if possible + if n == 0 { + return 0 + } + + // Placeholder for the actual LAPACK function calls + // Example: info = dgehrd(n, ilo, ihi, a, lda, tau, work, lwork) + return info +} diff --git a/lapack/lapack64/dgehrd.v b/lapack/lapack64/dgehrd.v new file mode 100644 index 000000000..2823c0c4d --- /dev/null +++ b/lapack/lapack64/dgehrd.v @@ -0,0 +1,35 @@ +module lapack64 + +import math +import vsl.blas + +// dgehrd reduces a general real matrix A to upper Hessenberg form H by an orthogonal similarity transformation. +pub fn dgehrd(n int, ilo int, ihi int, mut a []f64, lda int, tau []f64) int { + if n == 0 { + return 0 + } + + mut info := 0 + if n < 0 { + info = -1 + } else if ilo < 1 || ilo > math.max(1, n) { + info = -2 + } else if ihi < math.min(ilo, n) || ihi > n { + info = -3 + } else if lda < math.max(1, n) { + info = -5 + } + + if info != 0 { + return info + } + + // Quick return if possible + if n == 0 { + return 0 + } + + // Placeholder for the actual LAPACK function calls + // Example: info = dgehrd(n, ilo, ihi, a, lda, tau, work, lwork) + return info +} diff --git a/lapack/lapack64/dgesv.v b/lapack/lapack64/dgesv.v index d0d7f90a5..94f948503 100644 --- a/lapack/lapack64/dgesv.v +++ b/lapack/lapack64/dgesv.v @@ -22,7 +22,7 @@ import vsl.blas // The factored form of A is then used to solve the system of equations A * X = // B. On entry, b contains the right hand side matrix B. On return, if ok is // true, b contains the solution matrix X. -pub fn dgesv(n int, nrhs int, mut a []f64, lda int, ipiv []int, mut b []f64, ldb int) { +pub fn dgesv(n int, nrhs int, mut a []f64, lda int, mut ipiv []int, mut b []f64, ldb int) { if n < 0 { panic(n_lt0) } @@ -51,6 +51,6 @@ pub fn dgesv(n int, nrhs int, mut a []f64, lda int, ipiv []int, mut b []f64, ldb panic(short_b) } - dgetrf(n, n, mut a, lda, ipiv) + dgetrf(n, n, mut a, lda, mut ipiv) dgetrs(.no_trans, n, nrhs, mut a, lda, ipiv, mut b, ldb) } diff --git a/lapack/lapack64/dgesvd.v b/lapack/lapack64/dgesvd.v new file mode 100644 index 000000000..a8f035e7d --- /dev/null +++ b/lapack/lapack64/dgesvd.v @@ -0,0 +1,42 @@ +module lapack64 + +import math +import vsl.blas + +// dgesvd computes the singular value decomposition (SVD) of a real matrix A. +pub fn dgesvd(jobu SVDJob, jobvt SVDJob, m int, n int, mut a []f64, lda int, s []f64, mut u []f64, ldu int, mut vt []f64, ldvt int, superb []f64) int { + if m == 0 || n == 0 { + return 0 + } + + mut info := 0 + if jobu != .svd_all && jobu != .svd_store && jobu != .svd_overwrite && jobu != .svd_none { + info = -1 + } else if jobvt != .svd_all && jobvt != .svd_store && jobvt != .svd_overwrite + && jobvt != .svd_none { + info = -2 + } else if m < 0 { + info = -3 + } else if n < 0 { + info = -4 + } else if lda < math.max(1, m) { + info = -6 + } else if ldu < 1 || (jobu == .svd_store && ldu < m) || (jobu == .svd_all && ldu < m) { + info = -9 + } else if ldvt < 1 || (jobvt == .svd_store && ldvt < n) || (jobvt == .svd_all && ldvt < n) { + info = -11 + } + + if info != 0 { + return info + } + + // Quick return if possible + if m == 0 || n == 0 { + return 0 + } + + // Placeholder for the actual LAPACK function calls + // Example: info = dgesvd(jobu, jobvt, m, n, a, lda, s, u, ldu, vt, ldvt, work, lwork) + return info +} diff --git a/lapack/lapack64/dgetrf.v b/lapack/lapack64/dgetrf.v index e1b6553fc..878c7b5b1 100644 --- a/lapack/lapack64/dgetrf.v +++ b/lapack/lapack64/dgetrf.v @@ -24,7 +24,7 @@ import vsl.blas // Dgetrf returns whether the matrix A is nonsingular. The LU decomposition will // be computed regardless of the singularity of A, but the result should not be // used to solve a system of equation. -pub fn dgetrf(m int, n int, mut a []f64, lda int, ipiv []int) { +pub fn dgetrf(m int, n int, mut a []f64, lda int, mut ipiv []int) { mn := math.min(m, n) if m < 0 { @@ -34,7 +34,7 @@ pub fn dgetrf(m int, n int, mut a []f64, lda int, ipiv []int) { panic(n_lt0) } if lda < math.max(1, n) { - panic(bad_lda) + panic(bad_ld_a) } // quick return if possible @@ -53,7 +53,8 @@ pub fn dgetrf(m int, n int, mut a []f64, lda int, ipiv []int) { if nb <= 1 || nb >= mn { // use the unblocked algorithm. - return dgetf2(m, n, mut a, lda, ipiv) + dgetf2(m, n, mut a, lda, ipiv) + return } for j := 0; j < mn; j += nb { @@ -75,12 +76,12 @@ pub fn dgetrf(m int, n int, mut a []f64, lda int, ipiv []int) { dlaswp(j, mut slice1, lda, j, j + jb, ipiv[..j + jb], 1) mut slice2 := unsafe { a[j * lda + j + jb..] } - blas.dtstrf(.left, .lower, .notrans, .unit, jb, n - j - jb, 1, a[j * lda + j..], + blas.dtstrf(.left, false, false, .unit, jb, n - j - jb, 1, a[j * lda + j..], lda, mut slice2, lda) if j + jb < m { mut slice3 := unsafe { a[(j + jb) * lda + j + jb..] } - blas.dgemm(.notrans, .notrans, m - j - jb, n - j - jb, jb, -1, a[(j + jb) * lda + j..], + blas.dgemm(false, false, m - j - jb, n - j - jb, jb, -1, a[(j + jb) * lda + j..], lda, a[j * lda + j + jb..], lda, 1, mut slice3, lda) } } diff --git a/lapack/lapack64/dgetri.v b/lapack/lapack64/dgetri.v new file mode 100644 index 000000000..8cd14300e --- /dev/null +++ b/lapack/lapack64/dgetri.v @@ -0,0 +1,31 @@ +module lapack64 + +import math +import vsl.blas + +// dgetri computes the inverse of a matrix using the LU factorization computed by dgetrf. +pub fn dgetri(n int, mut a []f64, lda int, ipiv []int) int { + if n == 0 { + return 0 + } + + mut info := 0 + if n < 0 { + info = -1 + } else if lda < math.max(1, n) { + info = -3 + } + + if info != 0 { + return info + } + + // Quick return if possible + if n == 0 { + return 0 + } + + // Placeholder for the actual LAPACK function calls + // Example: info = dgetri(n, a, lda, ipiv, work, lwork) + return info +} diff --git a/lapack/lapack64/dgetrs.v b/lapack/lapack64/dgetrs.v index 0458cc1e2..bbbc586b0 100644 --- a/lapack/lapack64/dgetrs.v +++ b/lapack/lapack64/dgetrs.v @@ -52,16 +52,15 @@ pub fn dgetrs(trans blas.Transpose, n int, nrhs int, mut a []f64, lda int, ipiv // Solve A * X = B. dlaswp(nrhs, b, ldb, 0, n - 1, ipiv, 1) // Solve L * X = B, overwriting B with X. - blas.dtrsm(.left, .lower, .no_trans, .unit, n, nrhs, 1, mut a, lda, mut b, ldb) + blas.dtrsm(.left, false, false, .unit, n, nrhs, 1, a, lda, mut b, ldb) // Solve U * X = B, overwriting B with X. - blas.dtrsm(.left, .upper, .no_trans, .non_unit, n, nrhs, 1, mut a, lda, mut b, - ldb) + blas.dtrsm(.left, true, false, .non_unit, n, nrhs, 1, a, lda, mut b, ldb) } // Solve Aᵀ * X = B. // Solve Uᵀ * X = B, overwriting B with X. - blas.dtrsm(.left, .upper, .trans, .non_unit, n, nrhs, 1, mut a, lda, mut b, ldb) + blas.dtrsm(.left, true, true, .non_unit, n, nrhs, 1, a, lda, mut b, ldb) // Solve Lᵀ * X = B, overwriting B with X. - blas.dtrsm(.left, .lower, .trans, .unit, n, nrhs, 1, mut a, lda, mut b, ldb) + blas.dtrsm(.left, false, true, .unit, n, nrhs, 1, a, lda, mut b, ldb) dlaswp(nrhs, b, ldb, 0, n - 1, ipiv, -1) } diff --git a/lapack/lapack64/dpotrf.v b/lapack/lapack64/dpotrf.v new file mode 100644 index 000000000..b4d307809 --- /dev/null +++ b/lapack/lapack64/dpotrf.v @@ -0,0 +1,33 @@ +module lapack64 + +import math +import vsl.blas + +// dpotrf computes the Cholesky factorization of a real symmetric positive definite matrix A. +pub fn dpotrf(uplo blas.Uplo, n int, mut a []f64, lda int) int { + if n == 0 { + return 0 + } + + mut info := 0 + if uplo != .upper && uplo != .lower { + info = -1 + } else if n < 0 { + info = -2 + } else if lda < math.max(1, n) { + info = -4 + } + + if info != 0 { + return info + } + + // Quick return if possible + if n == 0 { + return 0 + } + + // Placeholder for the actual LAPACK function calls + // Example: info = dpotrf(uplo, n, a, lda, work, lwork) + return info +} diff --git a/lapack/lapack64/dsyev.v b/lapack/lapack64/dsyev.v new file mode 100644 index 000000000..a48998d2a --- /dev/null +++ b/lapack/lapack64/dsyev.v @@ -0,0 +1,38 @@ +module lapack64 + +import math +import vsl.blas + +// dsyev computes all eigenvalues and, optionally, eigenvectors of a real symmetric matrix A. +pub fn dsyev(jobz EVJob, uplo blas.Uplo, n int, mut a []f64, lda int, w []f64) int { + if n == 0 { + return 0 + } + + mut info := 0 + if jobz != .ev_none && jobz != .ev_compute { + info = -1 + } else if uplo != .upper && uplo != .lower { + info = -2 + } else if n < 0 { + info = -3 + } else if lda < math.max(1, n) { + info = -5 + } + + if info != 0 { + return info + } + + // Quick return if possible + if n == 0 { + return 0 + } + + // Call the relevant LAPACK functions + // (Here we would call the internal implementations like dsytrd, dorgtr, dormtr, etc.) + + // Placeholder for the actual LAPACK function calls + // Example: info = dsytrd(uplo, n, a, lda, w, work, lwork) + return info +} diff --git a/lapack/lapack64/ilaenv.v b/lapack/lapack64/ilaenv.v index 52d2e4e7d..08151271b 100644 --- a/lapack/lapack64/ilaenv.v +++ b/lapack/lapack64/ilaenv.v @@ -1,5 +1,7 @@ module lapack64 +import math + // ilaenv returns algorithm tuning parameters for the algorithm given by the // input string. ispec specifies the parameter to return: // @@ -242,7 +244,7 @@ fn ilaenv(ispec int, name string, opts string, n1 int, n2 int, n3 int, n4 int) i // Used by xGELSS and xGESVD // Assuming n1 and n2 are defined elsewhere in your code // Replace `min(n1, n2)` with actual min calculation or function - return int(f64(min(n1, n2)) * 1.6) + return int(f64(math.min(n1, n2)) * 1.6) } 7 { // Not used diff --git a/lapack/lapack_d_vsl_lapack_common.v b/lapack/lapack_d_vsl_lapack_lapacke.v similarity index 78% rename from lapack/lapack_d_vsl_lapack_common.v rename to lapack/lapack_d_vsl_lapack_lapacke.v index 0c4064c2b..cbc7d9839 100644 --- a/lapack/lapack_d_vsl_lapack_common.v +++ b/lapack/lapack_d_vsl_lapack_lapacke.v @@ -5,19 +5,19 @@ import vsl.blas fn C.LAPACKE_dgesv(matrix_layout blas.MemoryLayout, n int, nrhs int, a &f64, lda int, ipiv &int, b &f64, ldb int) int -fn C.LAPACKE_dgesvd(matrix_layout blas.MemoryLayout, jobu &char, jobvt &char, m int, n int, a &f64, lda int, s &f64, u &f64, ldu int, vt &f64, ldvt int, superb &f64) int +fn C.LAPACKE_dgesvd(matrix_layout blas.MemoryLayout, jobu SVDJob, jobvt SVDJob, m int, n int, a &f64, lda int, s &f64, u &f64, ldu int, vt &f64, ldvt int, superb &f64) int fn C.LAPACKE_dgetrf(matrix_layout blas.MemoryLayout, m int, n int, a &f64, lda int, ipiv &int) int fn C.LAPACKE_dgetri(matrix_layout blas.MemoryLayout, n int, a &f64, lda int, ipiv &int) int -fn C.LAPACKE_dpotrf(matrix_layout blas.MemoryLayout, up u32, n int, a &f64, lda int) int +fn C.LAPACKE_dpotrf(matrix_layout blas.MemoryLayout, uplo blas.Uplo, n int, a &f64, lda int) int -fn C.LAPACKE_dgeev(matrix_layout blas.MemoryLayout, calc_vl &char, calc_vr &char, n int, a &f64, lda int, wr &f64, wi &f64, vl &f64, ldvl_ int, vr &f64, ldvr_ int) int +fn C.LAPACKE_dgeev(matrix_layout blas.MemoryLayout, calc_vl LeftEVJob, calc_vr LeftEVJob, n int, a &f64, lda int, wr &f64, wi &f64, vl &f64, ldvl_ int, vr &f64, ldvr_ int) int -fn C.LAPACKE_dsyev(matrix_layout blas.MemoryLayout, jobz byte, uplo byte, n int, a &f64, lda int, w &f64, work &f64, lwork int) int +fn C.LAPACKE_dsyev(matrix_layout blas.MemoryLayout, jobz EVJob, uplo blas.Uplo, n int, a &f64, lda int, w &f64, work &f64, lwork int) int -fn C.LAPACKE_dgebal(matrix_layout blas.MemoryLayout, job &char, n int, a &f64, lda int, ilo int, ihi int, scale &f64) int +fn C.LAPACKE_dgebal(matrix_layout blas.MemoryLayout, job BalanceJob, n int, a &f64, lda int, ilo int, ihi int, scale &f64) int fn C.LAPACKE_dgehrd(matrix_layout blas.MemoryLayout, n int, ilo int, ihi int, a &f64, lda int, tau &f64, work &f64, lwork int) int @@ -74,9 +74,9 @@ pub fn dgesv(n int, nrhs int, mut a []f64, lda int, ipiv []int, mut b []f64, ldb // Note that the routine returns V**T, not V. // // NOTE: matrix 'a' will be modified -pub fn dgesvd(jobu &char, jobvt &char, m int, n int, a []f64, lda int, s []f64, u []f64, ldu int, vt []f64, ldvt int, superb []f64) { - info := C.LAPACKE_dgesvd(.row_major, jobu, jobvt, m, n, &a[0], lda, &s[0], &u[0], - ldu, &vt[0], ldvt, &superb[0]) +pub fn dgesvd(jobu SVDJob, jobvt SVDJob, m int, n int, mut a []f64, lda int, s []f64, mut u []f64, ldu int, mut vt []f64, ldvt int, superb []f64) { + info := C.LAPACKE_dgesvd(.row_major, jobu, jobvt, m, n, unsafe { &a[0] }, lda, &s[0], + unsafe { &u[0] }, ldu, unsafe { &vt[0] }, ldvt, &superb[0]) if info != 0 { errors.vsl_panic('lapack failed', .efailed) } @@ -141,9 +141,9 @@ pub fn dgetri(n int, mut a []f64, lda int, ipiv []int) { // where U is an upper triangular matrix and L is lower triangular. // // This is the block version of the algorithm, calling Level 3 BLAS. -pub fn dpotrf(up bool, n int, mut a []f64, lda int) { +pub fn dpotrf(uplo bool, n int, mut a []f64, lda int) { unsafe { - info := C.LAPACKE_dpotrf(.row_major, blas.l_uplo(up), n, &a[0], lda) + info := C.LAPACKE_dpotrf(.row_major, blas.c_uplo(uplo), n, &a[0], lda) if info != 0 { errors.vsl_panic('lapack failed', .efailed) } @@ -173,24 +173,24 @@ pub fn dpotrf(up bool, n int, mut a []f64, lda int) { // // The computed eigenvectors are normalized to have Euclidean norm // equal to 1 and largest component real. -pub fn dgeev(calc_vl bool, calc_vr bool, n int, mut a []f64, lda int, wr []f64, wi []f64, vl []f64, ldvl_ int, vr []f64, ldvr_ int) { +pub fn dgeev(calc_vl LeftEVJob, calc_vr LeftEVJob, n int, mut a []f64, lda int, wr []f64, wi []f64, vl []f64, ldvl_ int, vr []f64, ldvr_ int) { mut vvl := 0.0 mut vvr := 0.0 mut ldvl := ldvl_ mut ldvr := ldvr_ - if calc_vl { + if calc_vl == .left_ev_compute { vvl = vl[0] } else { ldvl = 1 } - if calc_vr { + if calc_vr == .left_ev_compute { vvr = vr[0] } else { ldvr = 1 } unsafe { - info := C.LAPACKE_dgeev(.row_major, &char(blas.job_vlr(calc_vl).str().str), &char(blas.job_vlr(calc_vr).str().str), - n, &a[0], lda, &wr[0], &wi[0], &vvl, ldvl, &vvr, ldvr) + info := C.LAPACKE_dgeev(.row_major, calc_vl, calc_vr, n, &a[0], lda, &wr[0], &wi[0], + &vvl, ldvl, &vvr, ldvr) if info != 0 { errors.vsl_panic('lapack failed', .efailed) } diff --git a/lapack/lapack_notd_vsl_lapack_common.v b/lapack/lapack_notd_vsl_lapack_lapacke.v similarity index 70% rename from lapack/lapack_notd_vsl_lapack_common.v rename to lapack/lapack_notd_vsl_lapack_lapacke.v index 3635ed6cd..39d2cbdd6 100644 --- a/lapack/lapack_notd_vsl_lapack_common.v +++ b/lapack/lapack_notd_vsl_lapack_lapacke.v @@ -4,20 +4,6 @@ import vsl.errors import vsl.blas import vsl.lapack.lapack64 -fn C.LAPACKE_dgesvd(matrix_layout blas.MemoryLayout, jobu &char, jobvt &char, m int, n int, a &f64, lda int, s &f64, u &f64, ldu int, vt &f64, ldvt int, superb &f64) int - -fn C.LAPACKE_dgetri(matrix_layout blas.MemoryLayout, n int, a &f64, lda int, ipiv &int) int - -fn C.LAPACKE_dpotrf(matrix_layout blas.MemoryLayout, up u32, n int, a &f64, lda int) int - -fn C.LAPACKE_dgeev(matrix_layout blas.MemoryLayout, calc_vl &char, calc_vr &char, n int, a &f64, lda int, wr &f64, wi &f64, vl &f64, ldvl_ int, vr &f64, ldvr_ int) int - -fn C.LAPACKE_dsyev(matrix_layout blas.MemoryLayout, jobz byte, uplo byte, n int, a &f64, lda int, w &f64, work &f64, lwork int) int - -fn C.LAPACKE_dgebal(matrix_layout blas.MemoryLayout, job &char, n int, a &f64, lda int, ilo int, ihi int, scale &f64) int - -fn C.LAPACKE_dgehrd(matrix_layout blas.MemoryLayout, n int, ilo int, ihi int, a &f64, lda int, tau &f64, work &f64, lwork int) int - // dgesv computes the solution to a real system of linear equations. // // See: http://www.netlib.org/lapack/explore-html/d8/d72/dgesv_8f.html @@ -65,9 +51,9 @@ pub fn dgesv(n int, nrhs int, mut a []f64, lda int, ipiv []int, mut b []f64, ldb // Note that the routine returns V**T, not V. // // NOTE: matrix 'a' will be modified -pub fn dgesvd(jobu &char, jobvt &char, m int, n int, a []f64, lda int, s []f64, u []f64, ldu int, vt []f64, ldvt int, superb []f64) { - info := C.LAPACKE_dgesvd(.row_major, jobu, jobvt, m, n, &a[0], lda, &s[0], &u[0], - ldu, &vt[0], ldvt, &superb[0]) +pub fn dgesvd(jobu SVDJob, jobvt SVDJob, m int, n int, mut a []f64, lda int, s []f64, mut u []f64, ldu int, mut vt []f64, ldvt int, superb []f64) { + info := lapack64.dgesvd(jobu, jobvt, m, n, mut a, lda, s, mut u, ldu, mut vt, ldvt, + superb) if info != 0 { errors.vsl_panic('lapack failed', .efailed) } @@ -102,11 +88,9 @@ pub fn dgetrf(m int, n int, mut a []f64, lda int, ipiv []int) { // This method inverts U and then computes inv(A) by solving the system // inv(A)*L = inv(U) for inv(A). pub fn dgetri(n int, mut a []f64, lda int, ipiv []int) { - unsafe { - info := C.LAPACKE_dgetri(.row_major, n, &a[0], lda, &ipiv[0]) - if info != 0 { - errors.vsl_panic('lapack failed', .efailed) - } + info := lapack64.dgetri(n, mut a, lda, ipiv) + if info != 0 { + errors.vsl_panic('lapack failed', .efailed) } } @@ -128,11 +112,9 @@ pub fn dgetri(n int, mut a []f64, lda int, ipiv []int) { // // This is the block version of the algorithm, calling Level 3 BLAS. pub fn dpotrf(up bool, n int, mut a []f64, lda int) { - unsafe { - info := C.LAPACKE_dpotrf(.row_major, blas.l_uplo(up), n, &a[0], lda) - if info != 0 { - errors.vsl_panic('lapack failed', .efailed) - } + info := lapack64.dpotrf(blas.c_uplo(up), n, mut a, lda) + if info != 0 { + errors.vsl_panic('lapack failed', .efailed) } } @@ -159,26 +141,28 @@ pub fn dpotrf(up bool, n int, mut a []f64, lda int) { // // The computed eigenvectors are normalized to have Euclidean norm // equal to 1 and largest component real. -pub fn dgeev(calc_vl bool, calc_vr bool, n int, mut a []f64, lda int, wr []f64, wi []f64, vl []f64, ldvl_ int, vr []f64, ldvr_ int) { +pub fn dgeev(calc_vl LeftEVJob, calc_vr LeftEVJob, n int, mut a []f64, lda int, wr []f64, wi []f64, mut vl []f64, ldvl_ int, mut vr []f64, ldvr_ int) { mut vvl := 0.0 mut vvr := 0.0 mut ldvl := ldvl_ mut ldvr := ldvr_ - if calc_vl { + if calc_vl == .left_ev_compute { vvl = vl[0] } else { ldvl = 1 } - if calc_vr { + if calc_vr == .left_ev_compute { vvr = vr[0] } else { ldvr = 1 } - unsafe { - info := C.LAPACKE_dgeev(.row_major, &char(blas.job_vlr(calc_vl).str().str), &char(blas.job_vlr(calc_vr).str().str), - n, &a[0], lda, &wr[0], &wi[0], &vvl, ldvl, &vvr, ldvr) - if info != 0 { - errors.vsl_panic('lapack failed', .efailed) - } + + vl[0] = vvl + vr[0] = vvr + + info := lapack64.dgeev(calc_vl, calc_vr, n, mut a, lda, wr, wi, mut vl, ldvl, mut + vr, ldvr) + if info != 0 { + errors.vsl_panic('lapack failed', .efailed) } }