Attribute VB_Name = "MATRIX_ARITHM_MULT_LIBR"

Option Explicit     'Requires that all variables to be declared explicitly.
Option Base 1       'The "Option Base" statement allows to specify 0 or 1 as the
                    'default first index of arrays.


'************************************************************************************
'************************************************************************************
'© Copyright NicoSystem 2009. All rights reserved by Rafael Nicolas Fermin Cota,
'San Francisco, CA. USA  www.rnfc.org
'nfermincota.hba2005@ivey.ca
'************************************************************************************
'************************************************************************************
'FUNCTION      : VECTOR_ELEMENTS_3D_PRODUCT_FUNC
'DESCRIPTION   : Return vector product (only 3 dimension)
'LIBRARY       : MATRIX
'GROUP         : ARITHMETIC_MULT
'ID            : 00X
'AUTHOR        : RAFAEL NICOLAS FERMIN COTA
'LAST UPDATE   : 01/22/2009
'************************************************************************************
'************************************************************************************

Function VECTOR_ELEMENTS_3D_PRODUCT_FUNC(ByRef ADATA_RNG As Variant, _
ByRef BDATA_RNG As Variant)

Dim ADATA_VECTOR As Variant
Dim BDATA_VECTOR As Variant

Dim TEMP_MATRIX As Variant

On Error GoTo ERROR_LABEL

ADATA_VECTOR = ADATA_RNG
    If UBound(ADATA_VECTOR, 1) = 1 Then: _
        ADATA_VECTOR = MATRIX_TRANSPOSE_FUNC(ADATA_VECTOR)

BDATA_VECTOR = BDATA_RNG
    If UBound(BDATA_VECTOR, 1) = 1 Then: _
        BDATA_VECTOR = MATRIX_TRANSPOSE_FUNC(BDATA_VECTOR)

If UBound(ADATA_VECTOR, 1) <> UBound(BDATA_VECTOR, 1) Or _
   UBound(ADATA_VECTOR, 1) <> 3 Then: GoTo ERROR_LABEL
    
ReDim TEMP_MATRIX(1 To 3, 1 To 1)
    
TEMP_MATRIX(1, 1) = ADATA_VECTOR(2, 1) * BDATA_VECTOR(3, 1) - _
                    BDATA_VECTOR(2, 1) * ADATA_VECTOR(3, 1)
    
TEMP_MATRIX(2, 1) = ADATA_VECTOR(3, 1) * BDATA_VECTOR(1, 1) - _
                    ADATA_VECTOR(1, 1) * BDATA_VECTOR(3, 1)
    
TEMP_MATRIX(3, 1) = ADATA_VECTOR(1, 1) * BDATA_VECTOR(2, 1) - _
                    BDATA_VECTOR(1, 1) * ADATA_VECTOR(2, 1)

VECTOR_ELEMENTS_3D_PRODUCT_FUNC = TEMP_MATRIX

Exit Function
ERROR_LABEL:
VECTOR_ELEMENTS_3D_PRODUCT_FUNC = Err.number
End Function

'************************************************************************************
'************************************************************************************
'© Copyright NicoSystem 2009. All rights reserved by Rafael Nicolas Fermin Cota,
'San Francisco, CA. USA  www.rnfc.org
'nfermincota.hba2005@ivey.ca
'************************************************************************************
'************************************************************************************
'FUNCTION      : VECTOR_ELEMENTS_DOT_PRODUCT_FUNC
'DESCRIPTION   : Dot Product Mult
'LIBRARY       : MATRIX
'GROUP         : ARITHMETIC_MULT
'ID            : 00X
'AUTHOR        : RAFAEL NICOLAS FERMIN COTA
'LAST UPDATE   : 01/22/2009
'************************************************************************************
'************************************************************************************

Function VECTOR_ELEMENTS_DOT_PRODUCT_FUNC(ByRef ADATA_RNG As Variant, _
ByRef BDATA_RNG As Variant)
    
Dim i As Long
Dim NROWS As Long

Dim TEMP_SUM As Double

Dim ADATA_VECTOR As Variant
Dim BDATA_VECTOR As Variant
        
On Error GoTo ERROR_LABEL
           
ADATA_VECTOR = ADATA_RNG
BDATA_VECTOR = BDATA_RNG

If IS_MATRIX_FUNC(ADATA_VECTOR) And IS_MATRIX_FUNC(BDATA_VECTOR) Then
    NROWS = UBound(ADATA_VECTOR, 1)
    TEMP_SUM = 0
    For i = 1 To NROWS
        TEMP_SUM = TEMP_SUM + ADATA_VECTOR(i, 1) * BDATA_VECTOR(i, 1)
    Next i
ElseIf IS_ARRAY_FUNC(ADATA_VECTOR) And IS_ARRAY_FUNC(BDATA_VECTOR) Then
    NROWS = UBound(ADATA_VECTOR)
    TEMP_SUM = 0
    For i = 1 To NROWS
        TEMP_SUM = TEMP_SUM + ADATA_VECTOR(i) * BDATA_VECTOR(i)
    Next i
Else
    GoTo ERROR_LABEL
End If

VECTOR_ELEMENTS_DOT_PRODUCT_FUNC = TEMP_SUM

Exit Function
ERROR_LABEL:
VECTOR_ELEMENTS_DOT_PRODUCT_FUNC = Err.number
End Function



'************************************************************************************
'************************************************************************************
'© Copyright NicoSystem 2009. All rights reserved by Rafael Nicolas Fermin Cota,
'San Francisco, CA. USA  www.rnfc.org
'nfermincota.hba2005@ivey.ca
'************************************************************************************
'************************************************************************************
'FUNCTION      : VECTOR_ELEMENTS_MULT_SCALAR_FUNC
'DESCRIPTION   : Multiplies all the numbers in a vector by a given SCALAR_VAL
'LIBRARY       : MATRIX
'GROUP         : ARITHMETIC_MULT
'ID            : 00X
'AUTHOR        : RAFAEL NICOLAS FERMIN COTA
'LAST UPDATE   : 01/22/2009
'************************************************************************************
'************************************************************************************

Function VECTOR_ELEMENTS_MULT_SCALAR_FUNC(ByRef DATA_RNG As Variant, _
ByVal SCALAR_VAL As Double, _
Optional ByVal VERSION As Integer = 0)
    
    Dim i As Long
    Dim NCOLUMNS As Long
    Dim DATA_VECTOR As Variant
    Dim TEMP_VECTOR As Variant
    
    On Error GoTo ERROR_LABEL
    
    DATA_VECTOR = DATA_RNG
    
    Select Case VERSION
        Case 0
            NCOLUMNS = UBound(DATA_VECTOR, 2)
            ReDim TEMP_VECTOR(1 To 1, 1 To NCOLUMNS)
            For i = 1 To NCOLUMNS
                TEMP_VECTOR(1, i) = DATA_VECTOR(1, i) * SCALAR_VAL
            Next i
        Case 1
            NCOLUMNS = UBound(DATA_VECTOR, 2)
            ReDim TEMP_VECTOR(1 To NCOLUMNS, 1 To 1)
            For i = 1 To NCOLUMNS
                TEMP_VECTOR(i, 1) = DATA_VECTOR(1, i) * SCALAR_VAL
            Next i
        Case 2
            NCOLUMNS = UBound(DATA_VECTOR, 1)
            ReDim TEMP_VECTOR(1 To 1, 1 To NCOLUMNS)
            For i = 1 To NCOLUMNS
                TEMP_VECTOR(1, i) = DATA_VECTOR(i, 1) * SCALAR_VAL
            Next i
        Case Else
            NCOLUMNS = UBound(DATA_VECTOR, 1)
            ReDim TEMP_VECTOR(1 To NCOLUMNS, 1 To 1)
            For i = 1 To NCOLUMNS
                TEMP_VECTOR(i, 1) = DATA_VECTOR(i, 1) * SCALAR_VAL
            Next i
    End Select

    VECTOR_ELEMENTS_MULT_SCALAR_FUNC = TEMP_VECTOR

Exit Function
ERROR_LABEL:
VECTOR_ELEMENTS_MULT_SCALAR_FUNC = Err.number
End Function


'************************************************************************************
'************************************************************************************
'© Copyright NicoSystem 2009. All rights reserved by Rafael Nicolas Fermin Cota,
'San Francisco, CA. USA  www.rnfc.org
'nfermincota.hba2005@ivey.ca
'************************************************************************************
'************************************************************************************
'FUNCTION      : VECTOR_ELEMENTS_MULT_FUNC
'DESCRIPTION   : This routine multiplies two vectors
'LIBRARY       : MATRIX
'GROUP         : ARITHMETIC_MULT
'ID            : 00X
'AUTHOR        : RAFAEL NICOLAS FERMIN COTA
'LAST UPDATE   : 01/22/2009
'************************************************************************************
'************************************************************************************

Function VECTOR_ELEMENTS_MULT_FUNC(ByRef ADATA_RNG As Variant, _
ByRef BDATA_RNG As Variant)
    
Dim i As Long
Dim NROWS As Long

Dim TEMP_VECTOR As Variant
Dim ADATA_VECTOR As Variant
Dim BDATA_VECTOR As Variant

On Error GoTo ERROR_LABEL
           
ADATA_VECTOR = ADATA_RNG
If UBound(ADATA_VECTOR, 1) = 1 Then
    ADATA_VECTOR = MATRIX_TRANSPOSE_FUNC(ADATA_VECTOR)
End If

BDATA_VECTOR = BDATA_RNG
If UBound(BDATA_VECTOR, 1) = 1 Then
    BDATA_VECTOR = MATRIX_TRANSPOSE_FUNC(BDATA_VECTOR)
End If
If UBound(ADATA_VECTOR, 1) <> UBound(BDATA_VECTOR, 1) Then: GoTo ERROR_LABEL
NROWS = UBound(ADATA_VECTOR, 1)
ReDim TEMP_VECTOR(1 To NROWS, 1 To 1)

For i = 1 To NROWS
    TEMP_VECTOR(i, 1) = ADATA_VECTOR(i, 1) * BDATA_VECTOR(i, 1)
Next i
    
VECTOR_ELEMENTS_MULT_FUNC = TEMP_VECTOR

Exit Function
ERROR_LABEL:
VECTOR_ELEMENTS_MULT_FUNC = Err.number
End Function



'************************************************************************************
'************************************************************************************
'© Copyright NicoSystem 2009. All rights reserved by Rafael Nicolas Fermin Cota,
'San Francisco, CA. USA  www.rnfc.org
'nfermincota.hba2005@ivey.ca
'************************************************************************************
'************************************************************************************
'FUNCTION      : VECTOR_ELEMENTS_MATRIX_VECTOR_MULT_FUNC
'DESCRIPTION   : This routine multiplies a vector into a square matrix
'and then into the transpose of the vector to yield a SCALAR_VAL
'outcome = xMx'
'LIBRARY       : MATRIX
'GROUP         : ARITHMETIC_MULT
'ID            : 00X
'AUTHOR        : RAFAEL NICOLAS FERMIN COTA
'LAST UPDATE   : 01/22/2009
'************************************************************************************
'************************************************************************************

Function VECTOR_ELEMENTS_MATRIX_VECTOR_MULT_FUNC(ByRef VECTOR_RNG As Variant, _
ByRef MATRIX_RNG As Variant)

Dim i As Long
Dim NROWS As Long

Dim TEMP_SUM As Double

Dim DATA_VECTOR As Variant
Dim DATA_MATRIX As Variant

On Error GoTo ERROR_LABEL

DATA_VECTOR = VECTOR_RNG
If UBound(DATA_VECTOR, 1) = 1 Then
    DATA_VECTOR = MATRIX_TRANSPOSE_FUNC(DATA_VECTOR)
End If

DATA_MATRIX = MATRIX_RNG
NROWS = UBound(DATA_VECTOR, 1)
If NROWS <> UBound(DATA_MATRIX, 2) Then: GoTo ERROR_LABEL

DATA_MATRIX = MMULT_FUNC(DATA_MATRIX, DATA_VECTOR)

TEMP_SUM = 0
For i = 1 To NROWS
    TEMP_SUM = TEMP_SUM + DATA_MATRIX(i, 1) * DATA_VECTOR(i, 1)
Next i

VECTOR_ELEMENTS_MATRIX_VECTOR_MULT_FUNC = TEMP_SUM

Exit Function
ERROR_LABEL:
VECTOR_ELEMENTS_MATRIX_VECTOR_MULT_FUNC = Err.number
End Function


'************************************************************************************
'************************************************************************************
'© Copyright NicoSystem 2009. All rights reserved by Rafael Nicolas Fermin Cota,
'San Francisco, CA. USA  www.rnfc.org
'nfermincota.hba2005@ivey.ca
'************************************************************************************
'************************************************************************************
'FUNCTION      : MATRIX_ELEMENTS_VECTOR_MULT_FUNC
'DESCRIPTION   : This routine multiplies a matrix with a vector to yield
'a vector.
'LIBRARY       : MATRIX
'GROUP         : ARITHMETIC_MULT
'ID            : 00X
'AUTHOR        : RAFAEL NICOLAS FERMIN COTA
'LAST UPDATE   : 01/22/2009
'************************************************************************************
'************************************************************************************

Function MATRIX_ELEMENTS_VECTOR_MULT_FUNC(ByRef VECTOR_RNG As Variant, _
ByRef MATRIX_RNG As Variant)

Dim NROWS As Long

Dim DATA_VECTOR As Variant
Dim DATA_MATRIX As Variant

On Error GoTo ERROR_LABEL

DATA_VECTOR = VECTOR_RNG
If UBound(DATA_VECTOR, 1) = 1 Then
    DATA_VECTOR = MATRIX_TRANSPOSE_FUNC(DATA_VECTOR)
End If
DATA_MATRIX = MATRIX_RNG

NROWS = UBound(DATA_VECTOR, 1)
If NROWS <> UBound(DATA_MATRIX, 2) Then: GoTo ERROR_LABEL
DATA_MATRIX = MMULT_FUNC(DATA_MATRIX, DATA_VECTOR)

MATRIX_ELEMENTS_VECTOR_MULT_FUNC = DATA_MATRIX

Exit Function
ERROR_LABEL:
MATRIX_ELEMENTS_VECTOR_MULT_FUNC = Err.number
End Function


'************************************************************************************
'************************************************************************************
'© Copyright NicoSystem 2009. All rights reserved by Rafael Nicolas Fermin Cota,
'San Francisco, CA. USA  www.rnfc.org
'nfermincota.hba2005@ivey.ca
'************************************************************************************
'************************************************************************************
'FUNCTION      : VECTOR_ELEMENTS_MATRIX_MULT_FUNC
'DESCRIPTION   : n x 1 Vector into (n x n) mult matrix
'LIBRARY       : MATRIX
'GROUP         : ARITHMETIC_MULT
'ID            : 00X
'AUTHOR        : RAFAEL NICOLAS FERMIN COTA
'LAST UPDATE   : 01/22/2009
'************************************************************************************
'************************************************************************************

Function VECTOR_ELEMENTS_MATRIX_MULT_FUNC(ByRef DATA_RNG As Variant)

Dim i As Long
Dim j As Long

Dim NROWS As Long
Dim DATA_VECTOR As Variant
Dim TEMP_MATRIX As Variant

On Error GoTo ERROR_LABEL

DATA_VECTOR = DATA_RNG

If UBound(DATA_VECTOR, 1) = 1 Then
    DATA_VECTOR = MATRIX_TRANSPOSE_FUNC(DATA_VECTOR)
End If
NROWS = UBound(DATA_VECTOR, 1)

ReDim TEMP_MATRIX(1 To NROWS, 1 To NROWS)

For i = 1 To NROWS
    TEMP_MATRIX(i, i) = DATA_VECTOR(i, 1) * DATA_VECTOR(i, 1)
    For j = 1 To i - 1
        TEMP_MATRIX(i, j) = DATA_VECTOR(i, 1) * DATA_VECTOR(j, 1)
        TEMP_MATRIX(j, i) = TEMP_MATRIX(i, j)
    Next j
Next i
VECTOR_ELEMENTS_MATRIX_MULT_FUNC = TEMP_MATRIX

Exit Function
ERROR_LABEL:
VECTOR_ELEMENTS_MATRIX_MULT_FUNC = Err.number
End Function

'************************************************************************************
'************************************************************************************
'© Copyright NicoSystem 2009. All rights reserved by Rafael Nicolas Fermin Cota,
'San Francisco, CA. USA  www.rnfc.org
'nfermincota.hba2005@ivey.ca
'************************************************************************************
'************************************************************************************
'FUNCTION      : MATRIX_ELEMENTS_PRODUCT_SUM_FUNC
'DESCRIPTION   : Matrix scalar product
'LIBRARY       : MATRIX
'GROUP         : ARITHMETIC_MULT
'ID            : 00X
'AUTHOR        : RAFAEL NICOLAS FERMIN COTA
'LAST UPDATE   : 01/22/2009
'************************************************************************************
'************************************************************************************

Function MATRIX_ELEMENTS_PRODUCT_SUM_FUNC(ByRef DATA_RNG As Variant, _
ByVal j As Long, _
ByVal k As Long)

Dim i As Long
Dim NSIZE As Long
Dim TEMP_SUM As Double
Dim DATA_MATRIX As Variant

On Error GoTo ERROR_LABEL
    
DATA_MATRIX = DATA_RNG
NSIZE = UBound(DATA_MATRIX)
TEMP_SUM = 0
For i = 1 To NSIZE
    TEMP_SUM = TEMP_SUM + DATA_MATRIX(i, j) * DATA_MATRIX(i, k)
Next i
MATRIX_ELEMENTS_PRODUCT_SUM_FUNC = TEMP_SUM

Exit Function
ERROR_LABEL:
MATRIX_ELEMENTS_PRODUCT_SUM_FUNC = Err.number
End Function


'************************************************************************************
'************************************************************************************
'© Copyright NicoSystem 2009. All rights reserved by Rafael Nicolas Fermin Cota,
'San Francisco, CA. USA  www.rnfc.org
'nfermincota.hba2005@ivey.ca
'************************************************************************************
'************************************************************************************
'FUNCTION      : MATRIX_ELEMENTS_MULT_SCALAR_FUNC
'DESCRIPTION   : Matrix SCALAR_VAL multiplication
'LIBRARY       : MATRIX
'GROUP         : ARITHMETIC_MULT
'ID            : 00X
'AUTHOR        : RAFAEL NICOLAS FERMIN COTA
'LAST UPDATE   : 01/22/2009
'************************************************************************************
'************************************************************************************

Function MATRIX_ELEMENTS_MULT_SCALAR_FUNC(ByRef DATA_RNG As Variant, _
ByVal SCALAR_VAL As Double)
  
  Dim i As Long
  Dim j As Long
  
  Dim NROWS As Long
  Dim NCOLUMNS As Long
  
  Dim DATA_MATRIX As Variant
  Dim TEMP_MATRIX As Variant
  
  On Error GoTo ERROR_LABEL
  
  DATA_MATRIX = DATA_RNG
  
  NROWS = UBound(DATA_MATRIX, 1)
  NCOLUMNS = UBound(DATA_MATRIX, 2)
  
  ReDim TEMP_MATRIX(1 To NROWS, 1 To NCOLUMNS)
  
  For i = 1 To NROWS
    For j = 1 To NCOLUMNS
          TEMP_MATRIX(i, j) = SCALAR_VAL * DATA_MATRIX(i, j)
    Next j
  Next i
  
MATRIX_ELEMENTS_MULT_SCALAR_FUNC = TEMP_MATRIX
Exit Function
ERROR_LABEL:
MATRIX_ELEMENTS_MULT_SCALAR_FUNC = Err.number
End Function


'************************************************************************************
'************************************************************************************
'© Copyright NicoSystem 2009. All rights reserved by Rafael Nicolas Fermin Cota,
'San Francisco, CA. USA  www.rnfc.org
'nfermincota.hba2005@ivey.ca
'************************************************************************************
'************************************************************************************
'FUNCTION      : MATRIX_ELEMENTS_MULT_FUNC
'DESCRIPTION   : Returns the M = aA x bB
'LIBRARY       : MATRIX
'GROUP         : ARITHMETIC_MULT
'ID            : 00X
'AUTHOR        : RAFAEL NICOLAS FERMIN COTA
'LAST UPDATE   : 01/22/2009
'************************************************************************************
'************************************************************************************

Function MATRIX_ELEMENTS_MULT_FUNC(ByRef ADATA_RNG As Variant, _
ByRef BDATA_RNG As Variant, _
Optional ByVal ASCALAR_VAL As Double = 1, _
Optional ByVal BSCALAR_VAL As Double = 1)
  
  Dim i As Long
  Dim j As Long
  
  Dim ANROWS As Long
  Dim ANCOLUMNS As Long
  
  Dim BNROWS As Long
  Dim BNCOLUMNS As Long
  
  Dim ADATA_MATRIX As Variant
  Dim BDATA_MATRIX As Variant
  
  Dim TEMP_MATRIX As Variant

  On Error GoTo ERROR_LABEL

  ADATA_MATRIX = ADATA_RNG
  BDATA_MATRIX = BDATA_RNG

  ANROWS = UBound(ADATA_MATRIX, 1)
  BNROWS = UBound(BDATA_MATRIX, 1)
  
  ANCOLUMNS = UBound(ADATA_MATRIX, 2)
  BNCOLUMNS = UBound(BDATA_MATRIX, 2)
   
'  If (ANROWS <> BNROWS) Or (ANCOLUMNS <> BNCOLUMNS) Then: GoTo ERROR_LABEL

    ReDim TEMP_MATRIX(1 To ANROWS, 1 To ANCOLUMNS)
    For i = 1 To ANROWS
      For j = 1 To ANCOLUMNS
        TEMP_MATRIX(i, j) = (ASCALAR_VAL * ADATA_MATRIX(i, j)) * _
            (BSCALAR_VAL * BDATA_MATRIX(i, j))
      Next j
    Next i
    MATRIX_ELEMENTS_MULT_FUNC = TEMP_MATRIX
  
Exit Function
ERROR_LABEL:
MATRIX_ELEMENTS_MULT_FUNC = Err.number
End Function


'************************************************************************************
'************************************************************************************
'© Copyright NicoSystem 2009. All rights reserved by Rafael Nicolas Fermin Cota,
'San Francisco, CA. USA  www.rnfc.org
'nfermincota.hba2005@ivey.ca
'************************************************************************************
'************************************************************************************
'FUNCTION      : MMULT_FUNC
'DESCRIPTION   : Fast matrix multiplication without size limitation
'LIBRARY       : MATRIX
'GROUP         : ARITHMETIC_MULT
'ID            : 00X
'AUTHOR        : RAFAEL NICOLAS FERMIN COTA
'LAST UPDATE   : 01/22/2009
'************************************************************************************
'************************************************************************************


Function MMULT_FUNC(ByRef ADATA_RNG As Variant, _
ByRef BDATA_RNG As Variant, _
Optional ByVal THRESHOLD As Single = 70)

Dim i As Single
Dim j As Single
Dim k As Single

Dim ii As Single
Dim jj As Single
Dim kk As Single

Dim iii As Single
Dim jjj As Single
Dim kkk As Single

Dim IMAX_VAL As Single
Dim IMIN_VAL As Single

Dim JMAX_VAL As Single
Dim JMIN_VAL As Single

Dim KMAX_VAL As Single
Dim KMIN_VAL As Single

Dim NSIZE As Single
Dim MSIZE As Single
Dim PSIZE As Single

Dim ANROWS As Single
Dim ANCOLUMNS As Single
Dim BNCOLUMNS As Single

Dim ATEMP_VECTOR As Variant
Dim BTEMP_VECTOR As Variant
Dim CTEMP_VECTOR As Variant

Dim TEMP_MATRIX As Variant

Dim ADATA_MATRIX As Variant
Dim BDATA_MATRIX As Variant

On Error GoTo ERROR_LABEL

ADATA_MATRIX = ADATA_RNG
BDATA_MATRIX = BDATA_RNG
  
If UBound(ADATA_MATRIX, 2) <> UBound(BDATA_MATRIX, 1) Then: GoTo ERROR_LABEL
  
ANROWS = UBound(ADATA_MATRIX, 1)  'rows of ADATA_MATRIX
ANCOLUMNS = UBound(ADATA_MATRIX, 2)
  'columns of ADATA_MATRIX = rows of BDATA_MATRIX
BNCOLUMNS = UBound(BDATA_MATRIX, 2)  'columns of BDATA_MATRIX

If ANROWS <= THRESHOLD And BNCOLUMNS <= THRESHOLD Then
    'fast multiplication
    MMULT_FUNC = MMULT2_FUNC(ADATA_MATRIX, BDATA_MATRIX)
    Exit Function
End If

'sub-matrix multiplication begins
NSIZE = Int(ANROWS / THRESHOLD)   'row-blocks of ADATA_MATRIX
PSIZE = Int(ANCOLUMNS / THRESHOLD)
'column-blocks of ADATA_MATRIX = row-blocks of BDATA_MATRIX
MSIZE = Int(BNCOLUMNS / THRESHOLD)   'column-blocks of BDATA_MATRIX

If NSIZE * THRESHOLD < ANROWS Then NSIZE = NSIZE + 1
If PSIZE * THRESHOLD < ANCOLUMNS Then PSIZE = PSIZE + 1
If MSIZE * THRESHOLD < BNCOLUMNS Then MSIZE = MSIZE + 1

ReDim TEMP_MATRIX(1 To ANROWS, 1 To BNCOLUMNS)

For ii = 1 To NSIZE
    For jj = 1 To MSIZE
        For kk = 1 To PSIZE

            'extract the sub-matrix ADATA_MATRIX(ii, kk) -> ATEMP_VECTOR
            IMIN_VAL = THRESHOLD * (ii - 1) + 1
            IMAX_VAL = THRESHOLD * ii
    
            If IMAX_VAL > ANROWS Then IMAX_VAL = ANROWS
            KMIN_VAL = THRESHOLD * (kk - 1) + 1
            KMAX_VAL = THRESHOLD * kk
    
            If KMAX_VAL > ANCOLUMNS Then KMAX_VAL = ANCOLUMNS
            iii = IMAX_VAL - IMIN_VAL + 1
            jjj = KMAX_VAL - KMIN_VAL + 1
    
            ReDim ATEMP_VECTOR(1 To iii, 1 To jjj)
            For i = 1 To UBound(ATEMP_VECTOR, 1)
                For k = 1 To UBound(ATEMP_VECTOR, 2)
                    ATEMP_VECTOR(i, k) = ADATA_MATRIX(i + IMIN_VAL - 1, _
                        k + KMIN_VAL - 1)
                Next k
            Next i

            'extract the sub-matrix BDATA_MATRIX(kk, jj) -> BTEMP_VECTOR
            KMIN_VAL = THRESHOLD * (kk - 1) + 1
            KMAX_VAL = THRESHOLD * kk
            If KMAX_VAL > ANCOLUMNS Then KMAX_VAL = ANCOLUMNS
            JMIN_VAL = THRESHOLD * (jj - 1) + 1
            JMAX_VAL = THRESHOLD * jj
            If JMAX_VAL > BNCOLUMNS Then JMAX_VAL = BNCOLUMNS
            jjj = KMAX_VAL - KMIN_VAL + 1
            kkk = JMAX_VAL - JMIN_VAL + 1
    
            ReDim BTEMP_VECTOR(1 To jjj, 1 To kkk)
            For k = 1 To UBound(BTEMP_VECTOR, 1)
                For j = 1 To UBound(BTEMP_VECTOR, 2)
                    BTEMP_VECTOR(k, j) = BDATA_MATRIX(k + KMIN_VAL - 1, _
                        j + JMIN_VAL - 1)
                Next j
            Next k
    
            'performs the multiplication of the sub-matrices
            CTEMP_VECTOR = MMULT2_FUNC(ATEMP_VECTOR, BTEMP_VECTOR)
    
            IMIN_VAL = THRESHOLD * (ii - 1) + 1
            JMIN_VAL = THRESHOLD * (jj - 1) + 1
    
            'accumulate the sub-matrix result
            For i = 1 To UBound(CTEMP_VECTOR, 1)
                For j = 1 To UBound(CTEMP_VECTOR, 2)
                    iii = i + IMIN_VAL - 1
                    jjj = j + JMIN_VAL - 1
                    TEMP_MATRIX(iii, jjj) = TEMP_MATRIX(iii, jjj) + CTEMP_VECTOR(i, j)
                Next j
            Next i
    
        Next kk
    Next jj
Next ii

MMULT_FUNC = TEMP_MATRIX
Exit Function
ERROR_LABEL:
MMULT_FUNC = Err.number
End Function

'************************************************************************************
'************************************************************************************
'© Copyright NicoSystem 2009. All rights reserved by Rafael Nicolas Fermin Cota,
'San Francisco, CA. USA  www.rnfc.org
'nfermincota.hba2005@ivey.ca
'************************************************************************************
'************************************************************************************
'FUNCTION      : MMULT2_FUNC
'DESCRIPTION   : Returns the matrix product of two arrays. The result is an array
'with the same number of rows as array1 and the same number of columns as array.
'LIBRARY       : MATRIX
'GROUP         : ARITHMETIC_MULT
'ID            : 015
'AUTHOR        : RAFAEL NICOLAS FERMIN COTA
'LAST UPDATE   : 01/22/2009
'************************************************************************************
'************************************************************************************

Function MMULT2_FUNC(ByRef ADATA_RNG As Variant, _
ByRef BDATA_RNG As Variant)

Dim i As Long
Dim j As Long
Dim k As Long

Dim ANROWS As Long
Dim BNROWS As Long

Dim ANCOLUMNS As Long
Dim BNCOLUMNS As Long

Dim ADATA_MATRIX As Variant
Dim BDATA_MATRIX As Variant

Dim TEMP_MATRIX As Variant

On Error GoTo ERROR_LABEL

ADATA_MATRIX = ADATA_RNG
BDATA_MATRIX = BDATA_RNG
  
ANROWS = UBound(ADATA_MATRIX, 1)
BNROWS = UBound(BDATA_MATRIX, 1)
  
ANCOLUMNS = UBound(ADATA_MATRIX, 2)
BNCOLUMNS = UBound(BDATA_MATRIX, 2)
  
If ANCOLUMNS <> BNROWS Then: GoTo ERROR_LABEL
    
ReDim TEMP_MATRIX(1 To ANROWS, 1 To BNCOLUMNS)

For i = 1 To ANROWS
    For j = 1 To BNCOLUMNS
        TEMP_MATRIX(i, j) = 0
        For k = 1 To ANCOLUMNS
            TEMP_MATRIX(i, j) = TEMP_MATRIX(i, j) + ADATA_MATRIX(i, k) * _
                                BDATA_MATRIX(k, j)
        Next k
    Next j
Next i
    
MMULT2_FUNC = TEMP_MATRIX

Exit Function
ERROR_LABEL:
MMULT2_FUNC = Err.number
End Function

'************************************************************************************
'************************************************************************************
'© Copyright NicoSystem 2009. All rights reserved by Rafael Nicolas Fermin Cota,
'San Francisco, CA. USA  www.rnfc.org
'nfermincota.hba2005@ivey.ca
'************************************************************************************
'************************************************************************************
'FUNCTION      : VECTOR_VANDERMONDE_MATRIX_FUNC
'DESCRIPTION   : Returns the Vandermonde's matrix for a given vector:
'x = (x1, x2, ...xn)
'LIBRARY       : MATRIX
'GROUP         : ARITHMETIC_MULT
'ID            : 00X
'AUTHOR        : RAFAEL NICOLAS FERMIN COTA
'LAST UPDATE   : 01/22/2009
'************************************************************************************
'************************************************************************************

Function VECTOR_VANDERMONDE_MATRIX_FUNC(ByRef DATA_RNG As Variant)

Dim i As Long
Dim j As Long

Dim NROWS As Long
Dim NCOLUMNS As Long

Dim TEMP_MATRIX As Variant
Dim DATA_VECTOR As Variant

On Error GoTo ERROR_LABEL

DATA_VECTOR = DATA_RNG
If UBound(DATA_VECTOR, 1) = 1 Then
    DATA_VECTOR = MATRIX_TRANSPOSE_FUNC(DATA_VECTOR)
End If
NROWS = UBound(DATA_VECTOR, 1)
NCOLUMNS = UBound(DATA_VECTOR, 2)
'If NROWS > 1 And NCOLUMNS > 1 Then: GoTo ERROR_LABEL

ReDim TEMP_MATRIX(1 To NROWS, 1 To NROWS)

For i = 1 To NROWS
    TEMP_MATRIX(i, 1) = 1
    For j = 2 To NROWS
        TEMP_MATRIX(i, j) = TEMP_MATRIX(i, j - 1) * DATA_VECTOR(i, 1)
    Next j
Next i

VECTOR_VANDERMONDE_MATRIX_FUNC = TEMP_MATRIX

Exit Function
ERROR_LABEL:
VECTOR_VANDERMONDE_MATRIX_FUNC = Err.number
End Function