Attribute VB_Name = "NCOLES"
Option Explicit

Private nOutFileNo As Integer
Private nMtrxSize As Integer
Private GAMtrx() As Double
Private GUMtrx() As Double
Private GDMtrx() As Double
Private GLMtrx() As Double
Private GFMtrx() As Double
Private GYMtrx() As Double
Private GXMtrx() As Double

Public Sub NCLSetData(nSize As Integer, AMtrx() As Double, FMtrx() As Double)
    Dim i, j As Integer
    
    nMtrxSize = nSize
    
    '計算に必要な領域を確保
    ReDim GAMtrx(nMtrxSize, nMtrxSize)
    ReDim GUMtrx(nMtrxSize, nMtrxSize)
    ReDim GDMtrx(nMtrxSize, nMtrxSize)
    ReDim GLMtrx(nMtrxSize, nMtrxSize)
    ReDim GFMtrx(nMtrxSize)
    ReDim GYMtrx(nMtrxSize)
    ReDim GXMtrx(nMtrxSize)
    
    'パラメータで渡されたマトリクスをコピー
    For i = 0 To (nMtrxSize - 1)
        For j = 0 To (nMtrxSize - 1)
            GAMtrx(i, j) = AMtrx(i, j)
        Next j
        GFMtrx(i) = FMtrx(i)
    Next i
    
End Sub

Public Sub NCLCalculate()
    Dim i, j, k As Integer
    Dim dblS1, dblS2, dblS3, dblS4, dblS5 As Double

    '三角分解
    For j = 0 To (nMtrxSize - 1)
        For i = 0 To (nMtrxSize - 1)
            '下三角行列
            If i > j Then
                dblS1 = 0
                If i < 2 Then
                Else
                    For k = 0 To (i - 2)
                        dblS1 = dblS1 + GLMtrx(i, k) * GUMtrx(k, j) * GDMtrx(k, k)
                    Next k
                End If
        
                GLMtrx(i, j) = (GAMtrx(i, j) - dblS1) / GDMtrx(j, j)
            '上三角行列
            ElseIf i < j Then
                dblS2 = 0
                If j < 2 Then
                Else
                    For k = 0 To (j - 2)
                        dblS2 = dblS2 + GLMtrx(i, k) * GUMtrx(k, j) * GDMtrx(k, k)
                    Next k
                End If
                
                GUMtrx(i, j) = (GAMtrx(i, j) - dblS2) / GDMtrx(i, i)
            '対角行列
            Else
                dblS3 = 0
                If i = 0 Then
                Else
                    For k = 0 To i
                        dblS3 = dblS3 + GLMtrx(i, k) * GUMtrx(k, i) * GDMtrx(k, k)
                    Next k
                End If
                
                GDMtrx(i, i) = GAMtrx(i, i) - dblS3
                GLMtrx(i, i) = 1
                GUMtrx(i, i) = 1
            End If
        Next i
    Next j
            
    GYMtrx(0) = GFMtrx(0)
    
    For i = 0 To (nMtrxSize - 1)
        dblS4 = 0
        For j = 0 To (i - 1)
            dblS4 = dblS4 + GLMtrx(i, j) * GYMtrx(j)
        Next j
        
        GYMtrx(i) = GFMtrx(i) - dblS4
    Next i
            
    GXMtrx(nMtrxSize - 1) = GYMtrx(nMtrxSize - 1) / GDMtrx(nMtrxSize - 1, nMtrxSize - 1)

    For i = nMtrxSize - 2 To 0 Step -1
        dblS5 = 0
        For j = i + 1 To (nMtrxSize - 1)
            dblS5 = dblS5 + GUMtrx(i, j) * GXMtrx(j)
        Next j
        GXMtrx(i) = GYMtrx(i) / GDMtrx(i, i) - dblS5
    Next i
    
    Call OutputData
    
End Sub
        
Private Sub OutputData()
    Dim strOutputFile As String
    
    nOutFileNo = 10
    strOutputFile = App.Path & "\NColes.txt"
'    strOutputFile = ActiveWorkbook.Path & "\NColes.txt"

    Open strOutputFile For Output As #nOutFileNo

    Print #nOutFileNo, "連立方程式の解法：コレスキー法（非対称可）" & Chr(13) & Chr(13)

    Call MtrxOutput1("系行列", GAMtrx)
    Call MtrxOutput2("定数", GFMtrx)
    Call MtrxOutput1("上三角行列", GUMtrx)
    Call MtrxOutput1("対角角行列", GDMtrx)
    Call MtrxOutput1("下三角行列", GLMtrx)
    Call MtrxOutput2("解", GXMtrx)

    Close #nOutFileNo

End Sub

Private Sub MtrxOutput1(strComment As String, Mtrx() As Double)
    Dim i, j As Integer
    Dim strLineData As String
    
    Print #nOutFileNo, strComment

    For i = 0 To (nMtrxSize - 1)
        strLineData = ""
        For j = 0 To (nMtrxSize - 1)
            strLineData = strLineData & Format(Mtrx(i, j), "0.00000")
        Next j
        Print #nOutFileNo, strLineData
    Next i
    Print #nOutFileNo, ""

End Sub

Private Sub MtrxOutput2(strComment As String, Mtrx() As Double)
    Dim i As Integer
    Dim strLineData As String
    
    Print #nOutFileNo, strComment

    strLineData = ""
    For i = 0 To (nMtrxSize - 1)
        strLineData = strLineData & Format(Mtrx(i), "0.00000")
    Next i
    Print #nOutFileNo, strLineData
    Print #nOutFileNo, ""

End Sub

Public Sub NCLGetData(nMtrxSize As Integer, Mtrx() As Double)
    Dim i As Integer

    For i = 0 To (nMtrxSize - 1)
        Mtrx(i) = GXMtrx(i)
    Next i

End Sub