Симуляция Монте-Карло, включая гауссовское исключение с использованием VBA - PullRequest
Новая
       6

Симуляция Монте-Карло, включая гауссовское исключение с использованием VBA

1 голос
/ 26 марта 2019

Я написал код, который решает систему линейных уравнений времени сервировки подряд во время симуляции Монте-Карло.

Для каждого прогона вход немного изменяется, и решение нужно вычислять каждый раз снова.Цель этого состоит в том, чтобы получить функцию распределения вероятностей результатов (решения линейной системы).

Итак, мой вопрос:
Существует ли способ решения только системы линейных уравненийодин раз и сохраните общее решение, чтобы при каждом запуске в Монте-Карло решения можно было рассчитать напрямую?

Это сэкономило бы много времени, так как для правильного моделирования мне нужно как минимум 20 000 прогонов и даже для небольших систем из трехнеизвестные это занимает много времени.Мой код решает эти линейные уравнения каждый раз, когда в новом исходном коде число переменных и, следовательно, число входных величин должно быть замкнутым, поэтому общие решения неизвестны.

Вот мой алгоритм исключения Гаусса. 

Function gaussian_elimination(w As Variant, mm As Variant, R As Variant, rb As Variant, n_iso As Integer) As Variant()

'initializing running indexes
    Dim i As Integer
    Dim j As Integer
    Dim h As Integer
    Dim n As Integer

    n = n_iso

'runing variables for Gauss elimination
    Dim ip As Integer
    Dim q As Integer
    Dim p As Integer
    Dim z As Double
    Dim temp1(1, 1) As Variant
    Dim temp2(1, 1) As Variant
    Dim sum As Variant


'initializing b vector
    Dim b() As Variant
    ReDim b(1 To n - 1, 0 To 1)

'initializing k vector
    Dim k() As Variant
    ReDim k(1 To n - 1, 0 To 1)

'initializing A matrix
    Dim a() As Variant
    ReDim a(1 To n - 1, 1 To n - 1)

'initializing X matrix
    Dim x() As Variant
    ReDim x(1 To n - 1, 1 To n)


' calculating b vector

    For i = 1 To (n - 1) Step 1
        b(i, 0) = mm(1, 0) / (w(i + n - 1, 0) * (rb(i, 0) - R(i * n, 0))) - mm(1, 0) / (w(i, 0) * (R(i, 0) - rb(i, 0)))
    Next i

'calculating A matrix
    For i = 1 To (n - 1) Step 1
        For j = 1 To (n - 1) Step 1
            a(i, j) = mm(j + 1, 0) * ((R(j + i * (n - 1), 0) / (w(i + n - 1, 0) * (rb(i, 0) - R(i * n, 0)))) - (R(j, 0) / (w(i, 0) * (R(i, 0) - rb(i, 0)))))
        Next j
    Next i



    'using on board solving routine

    Dim A_Inv As Variant
    Dim k_vec As Variant
    Dim b_dummy As Variant





'filling X matrix
    For i = 1 To n - 1 Step 1
        For j = 1 To n Step 1
                If j = (n) Then
                    x(i, j) = b(i, 0)
                Else: x(i, j) = a(i, j)
                End If
        Next j
    Next i

   'Gaussian elimination
    For i = 1 To (n - 2) Step 1
        For j = i + 1 To (n - 1) Step 1
            If (Abs(x(j, i)) > Abs(x(i, i))) Then
                For h = 1 To n
                    temp1(1, 1) = x(i, h)
                    temp2(1, 1) = x(j, h)
                    x(i, h) = temp2(1, 1)
                    x(j, h) = temp1(1, 1)
                Next h
            End If
        Next

        For p = i + 1 To n - 1
            z = x(p, i) / x(i, i)
            For q = i + 1 To n
                x(p, q) = x(p, q) - z * x(i, q)
            Next q
            x(p, i) = 0
        Next p
    Next i






'calculatiing k factors backwards
    If Abs(x(UBound(x, 1), UBound(x, 2) - 1)) <= 0 Then
        MsgBox "Equation system can not be solved! Solving for k factors faild", vbExclamation, "Warning!"
        Exit Function
    End If

    k((UBound(x, 1)), 0) = x((UBound(x, 1)), UBound(x, 2)) / x((UBound(x, 1)), (UBound(x, 2) - 1))

    For i = ((UBound(x, 1) - 1)) To (LBound(x, 1)) Step -1
        sum = x(i, UBound(x, 2))
        For j = i + 1 To (UBound(x, 2) - 1) Step 1
            sum = sum - x(i, j) * k(j, 0)
        Next j
        k(i, 0) = sum / x(i, i)
    Next i

    For i = 1 To n - 1
        k(i, 0) = (-1) * k(i, 0)
    Next i

    gaussian_elimination = k



End Function
...