package scientific.linear

import std.unittest.*
import std.unittest.testmacro.*

import scientific.numbers.*

/*
 * Compute the value of 1-norm, Frobenius norm, infinity-norm
 * or the largest absolute value of a general rectangular matrix.
 */
foreign func LAPACKE_dlange(
    matrix_layout: IntNative,  // row or column major
    norm: UInt8,               // type of norm, one of 'F', 'I', 'l', 'M'
    m: IntNative,              // number of rows in A
    n: IntNative,              // number of columns in A
    a: CPointer<Unit>,         // matrix A (Float64)
    lda: IntNative             // leading dimension of A
): Float64

/*
 * Compute the value of 1-norm, Frobenius norm, infinity-norm
 * or the largest absolute value of a general rectangular matrix.
 */
foreign func LAPACKE_zlange(
    matrix_layout: IntNative,  // row or column major
    norm: UInt8,               // type of norm, one of 'F', 'I', 'l', 'M'
    m: IntNative,              // number of rows in A
    n: IntNative,              // number of columns in A
    a: CPointer<Unit>,         // matrix A (Complex64)
    lda: IntNative             // leading dimension of A
): Float64

public func norm(a: Matrix<Float64>, norm: Rune): Float64 {
    let LAPACK_COL_MAJOR: Int64 = 102
    let cnorm: UInt32 = UInt32(norm)
    let m = a.getRows()
    let n = a.getCols()
    let lda = m
    var work = zeros<Float64>(m)
    let DLANGE = unsafe{ LAPACKE_dlange(IntNative(LAPACK_COL_MAJOR), UInt8(cnorm),
                         IntNative(m), IntNative(n), a.ptr, IntNative(lda)) }
    return DLANGE
}

public func norm(a: Matrix<Complex64>, norm: Rune): Float64 {
    let LAPACK_COL_MAJOR: Int64 = 102
    let cnorm: UInt32 = UInt32(norm)
    let m = a.getRows()
    let n = a.getCols()
    let lda = m
    var work = zeros<Float64>(m)
    let DLANGE = unsafe{ LAPACKE_zlange(IntNative(LAPACK_COL_MAJOR), UInt8(cnorm),
                         IntNative(m), IntNative(n), a.ptr, IntNative(lda)) }
    return DLANGE
}

func testDlange() {
    let LAPACK_COL_MAJOR: Int64 = 102
    var a = matrix<Float64>(
        [[-4.0, -3.0, -2.0],
         [-1.0, 0.0, 1.0],
         [2.0, 3.0, 4.0]]
    )
    let norm: Rune = r'F'
    let cnorm: UInt32 = UInt32(norm)
    let m = a.getRows()
    let n = a.getCols()
    let lda = m
    var work = zeros<Float64>(m)
    let frobenius_norm = unsafe{ LAPACKE_dlange(IntNative(LAPACK_COL_MAJOR), UInt8(cnorm),
                                    IntNative(m), IntNative(n), a.ptr, IntNative(lda)) }
    let expected_frobenius_norm = 7.745967
    @Assert(approxEqual(frobenius_norm, expected_frobenius_norm, atol:1e-6))
}

@Test
public class TestNorm {
    @TestCase
    func testMatrixNorm(): Unit {
        var a = matrix<Float64>(
            [[-4.0, -3.0, -2.0],
             [-1.0, 0.0, 1.0],
             [2.0, 3.0, 4.0]]
        )
        let frobenius_norm = norm(a, r'F')
        let infinity_norm = norm(a, r'I')
        let one_norm = norm(a, r'1')
        let max_norm = norm(a, r'M')

        let expected_frobenius_norm = 7.745967
        let expected_infinity_norm = 9.0
        let expected_one_norm = 7.0
        let expected_max_norm = 4.0
        @Assert(approxEqual(frobenius_norm, expected_frobenius_norm, atol:1e-6))
        @Assert(approxEqual(infinity_norm, expected_infinity_norm, atol:1e-6))
        @Assert(approxEqual(one_norm, expected_one_norm, atol:1e-6))
        @Assert(approxEqual(max_norm, expected_max_norm, atol:1e-6))
    }

    @TestCase
    func testComplex64MatrixNorm(): Unit {
        var a = matrix<Complex64>([
            [Complex64(0.0, -4.0), Complex64(0.0, -3.0), Complex64(0.0, -2.0)],
            [Complex64(0.0, -1.0), Complex64(0.0,  0.0), Complex64(0.0,  1.0)],
            [Complex64(0.0,  2.0), Complex64(0.0,  3.0), Complex64(0.0,  4.0)]
        ])
        let frobenius_norm = norm(a, r'F')
        let infinity_norm = norm(a, r'I')
        let one_norm = norm(a, r'1')
        let max_norm = norm(a, r'M')

        let expected_frobenius_norm = 7.745967
        let expected_infinity_norm = 9.0
        let expected_one_norm = 7.0
        let expected_max_norm = 4.0
        @Assert(approxEqual(frobenius_norm, expected_frobenius_norm, atol:1e-6))
        @Assert(approxEqual(infinity_norm, expected_infinity_norm, atol:1e-6))
        @Assert(approxEqual(one_norm, expected_one_norm, atol:1e-6))
        @Assert(approxEqual(max_norm, expected_max_norm, atol:1e-6))
    }

    @TestCase
    func testLapackDlange(): Unit {
        testDlange()
    }
}