form1.frm

基于VB语言的微粒群算法

'########求lnr3(求α、β相中活度系数lnr3)##########
For i = 1 To DataNum
    Y31(i) = (x1b(i) * a(2) * G(2) + x2b(i) * a(4) * G(4)) * (x1b(i) * G(2) + x2b(i) * G(4)) / (x1b(i) * G(2) + x2b(i) * G(4) + x3b(i)) ^ 2
    Y32(i) = x1b(i) * G(5) * (a(5) * x1b(i) + a(5) * x2b(i) * G(3) - x2b(i) * G(3) * a(3)) / (x1b(i) + x2b(i) * G(3) + x3b(i) * G(5)) ^ 2
    Y33(i) = x2b(i) * G(6) * (x2b(i) * a(6) + x1b(i) * a(6) * G(6) - x1b(i) * G(1) * a(1)) / (x1b(i) * G(1) + x2b(i) + x3b(i) * G(6)) ^ 2
    Y34(i) = (x1a(i) * a(2) * G(2) + x2a(i) * G(4) * a(4)) * (x1a(i) * G(2) + x2a(i) * G(4)) / (x1a(i) * G(2) + x2a(i) * G(4) + x3a(i)) ^ 2
    Y35(i) = x1a(i) * G(5) * (a(5) * x1a(i) + a(5) * x2a(i) * G(3) - x2a(i) * a(3) * G(3)) / (x1a(i) + x2a(i) * G(3) + x3a(i) * G(5)) ^ 2
    Y36(i) = x2a(i) * G(6) * (x2a(i) * a(6) + x1a(i) * a(6) * G(6) - x1a(i) * a(1) * G(1)) / (x1a(i) * G(1) + x2a(i) + x3a(i) * G(6)) ^ 2
    '表达式3  lnγ3I=Y31+Y32+Y33   lnγ3II=Y34+Y35+Y36  Y3=ln(γ3I/γ3II)
Next i
For i = 1 To DataNum
    rb(i, 3) = Exp(Y31(i) + Y32(i) + Y33(i))
    ra(i, 3) = Exp(Y34(i) + Y35(i) + Y36(i))
Next i
'########求lnr3(求α、β相中活度系数lnr2)##########

For i = 1 To DataNum
    For j = 1 To 3
        KReal(i, j) = x1a(j) / x1b(j)
        KCal(i, j) = rb(i, j) / ra(i, j)
    Next j
Next i
FunctionObject = 0  '另初值为零，得注意这里
For i = 1 To DataNum
    For j = 1 To 3
        FunctionObject = FunctionObject + ((KReal(i, j) - KCal(i, j)) / KReal(i, j)) ^ 2  '目标函数  值越小越好
    Next j
Next i
FunctionObject = FunctionObject / 3 / DataNum
'以上目标函数变时需要修改
End Function

Public Function Pso(amin() As Double, amax() As Double, VaryNum As Integer, Popsize As Integer, MaxGeneration As Integer, ObjectValue As Integer)
'amin  变量的取值下限 amax 变量的取值上限
'VaryNum 变量数  MaxGeneration 最大迭代数目
'Popsize 微粒群规模
Dim RealValue(DataNum) As Double
Dim i%, j%, Generation%, PNum%
Dim c1#, c2#, Rnd1#, Rnd2#, w#
'##########重定义数组变量##############
Dim InBird() As Individual, IndividualBest() As Individual, GlobalBest As Individual
ReDim GlobalBest.v(VaryNum) As Double
ReDim GlobalBest.x(VaryNum) As Double
ReDim InBird(MaxGeneration, Popsize) As Individual, IndividualBest(Popsize) As Individual
For i = 0 To MaxGeneration
    For j = 1 To Popsize
        ReDim InBird(i, j).v(VaryNum) As Double
        ReDim InBird(i, j).x(VaryNum) As Double
    Next j
Next i
For i = 1 To Popsize
    ReDim IndividualBest(i).v(VaryNum) As Double
    ReDim IndividualBest(i).x(VaryNum) As Double
Next i
'##########重定义数组变量##############

Randomize
'##########初始化微粒群##############
For i = 1 To Popsize
    For j = 1 To VaryNum
        InBird(0, i).x(j) = amin(j) + Rnd * (amax(j) - amin(j)) '初始化变量值
    Next j
Next i
For i = 1 To Popsize
    InBird(0, i).Value = FunctionObject(InBird(0, i).x(), 0)      '通过调用函数计算目标函数值
    IndividualBest(i) = InBird(0, i)                            '初始化时将初始值作为最优值
Next i
GlobalBest = IndividualBest(1)
For i = 1 To Popsize
    If ObjectValue = 1 Then        '目标函数取最大值为目的
       If GlobalBest.Value < IndividualBest(i).Value Then
          GlobalBest = IndividualBest(i)                           '寻求初始群体中的最优位置
       End If
    ElseIf ObjectValue = 0 Then     '目标函数取最小值为目的
       If GlobalBest.Value > IndividualBest(i).Value Then
          GlobalBest = IndividualBest(i)                           '寻求初始群体中的最优位置
       End If
    End If
Next i
'##########初始化微粒群##############

For Generation = 1 To MaxGeneration
    w = 0.4 - Generation * (0.8 - 0.4) / 1000      '线性调整权重因子w的值
    c1 = 2: c2 = 2                                '比例参数
    For i = 1 To Popsize
        For j = 1 To VaryNum
Recalculation:                                    '重新计算
            Rnd1 = Rnd: Rnd2 = Rnd                '产生随机数
            '产生下一代微粒移动的距离
            InBird(Generation, i).v(j) = w * InBird(Generation - 1, i).v(j) + c1 * Rnd1 * (IndividualBest(i).x(j) - InBird(Generation - 1, i).x(j)) + c2 * Rnd2 * (GlobalBest.x(j) - InBird(Generation - 1, i).x(j))
            InBird(Generation, i).x(j) = InBird(Generation - 1, i).x(j) + InBird(Generation, i).v(j)
            If InBird(Generation, i).x(j) > amax(j) Or InBird(Generation, i).x(j) < amin(j) Then
               PNum = PNum + 1
               If PNum <= 10 Then
                  GoTo Recalculation
               Else
                  InBird(Generation, i).v(j) = Rnd
                  InBird(Generation, i).x(j) = amin(j) + (amax(j) - amin(j)) * InBird(Generation, i).v(j)
               End If
            End If
            PNum = 0
        Next j
    Next i
    For i = 1 To Popsize
        InBird(Generation, i).Value = FunctionObject(InBird(Generation, i).x(), Generation)
    Next i
    For i = 1 To Popsize   '寻求每一个微粒群个体运行至某一代时发现的最优位置
        If ObjectValue = 1 Then        '目标函数取最大值为目的
           If IndividualBest(i).Value < InBird(Generation, i).Value Then
              IndividualBest(i) = InBird(Generation, i)
           End If
        ElseIf ObjectValue = 0 Then        '目标函数取最小值为目的
           If IndividualBest(i).Value > InBird(Generation, i).Value Then
              IndividualBest(i) = InBird(Generation, i)
           End If
        End If
    Next i
    For i = 1 To Popsize
        If ObjectValue = 1 Then        '目标函数取最大值为目的
           If GlobalBest.Value < IndividualBest(i).Value Then
              GlobalBest = IndividualBest(i) '寻求整个微粒群体运行至某一代时发现的最优位置
           End If
        ElseIf ObjectValue = 0 Then       '目标函数取最小值为目的
           If GlobalBest.Value > IndividualBest(i).Value Then
              GlobalBest = IndividualBest(i) '寻求整个微粒群体运行至某一代时发现的最优位置
           End If
        End If
    Next i
    Debug.Print "当前运算到第"; Generation; " 代    ";
    Debug.Print "GlobalBest.Value="; GlobalBest.Value
Next Generation
    For i = 1 To VaryNum
        Debug.Print "GlobalBest.x("; i; ")="; GlobalBest.x(i)
    Next i
    Call SubValue(GlobalBest.x(), RealValue())
End Function


Public Sub SubValue(Vary() As Double, fValue() As Double)
  Dim Dimension As Integer, x() As Double, i%
  Dimension = UBound(Vary)
  ReDim x(Dimension) As Double
  For i = 1 To Dimension
      x(i) = Vary(i)
  Next i
    
   '以下目标函数变时需要修改
  Dim j%
  For j = 1 To DataNum
      fValue(j) = x(1) + x(2) * X0(j, 1) + x(3) * X0(j, 2) + x(4) * X0(j, 3) + x(5) * X0(j, 1) * X0(j, 2) + x(6) * X0(j, 1) * X0(j, 3) + x(7) * X0(j, 2) * X0(j, 3) + x(8) * X0(j, 1) ^ 2 + x(9) * X0(j, 2) ^ 2 + x(10) * X0(j, 3) ^ 2
  Next j
  '以上目标函数变时需要修改
    For i = 1 To DataNum
         Debug.Print "fValue("; i; ")="; fValue(i)
    Next i
End Sub