senne.vb

一个lbm的程序

Imports MyIO
Imports ss.Sne2
Imports ss.Collide
Imports ss.Stream
Public Class Senne
    Private num, cnt As Int32
    Private st(1) As Byte
    Private emptyArray As New List(Of Int32())
    Private fillArray, tmpArray As New List(Of Int32())
    Private en, exmass, vtotal, sitaw As Double
    Private avgrho, avgu, avgv As Double
    Private vi(q) As Double
    Event returnStep(ByVal num As Int32)
    Event returnComplete()
    Event returnWrite(ByVal str As String)

    Public Sub start()
        RaiseEvent returnWrite("started..." & strNextLine)
        RaiseEvent returnWrite("Itf number:" & itfArray.Count & strNextLine)
        Dim t As Int32
        FileWritter.writeTcp(strWrite & t, u, v, mass, "u", "v", "mass")
        For num = 1 To steps
            preStream()
            streaming()
            poststream()
            collision()
            reconstruction()
            RaiseEvent returnStep(num)
            If num Mod interval = 0 Then
                t += 1
                FileWritter.writeTcp(strWrite & t, u, v, mass, "u", "v", "mass")
                RaiseEvent returnWrite(num & " Done" & strNextLine & "Itf number:" & itfArray.Count & strNextLine)
            End If
        Next
        RaiseEvent returnComplete()
    End Sub

    Private Sub preStream()
        Dim i, j As Int32
        For i = 0 To m
            For j = 0 To n
                Select Case flag(i, j)
                    Case CONST_WALL
                        dealWallPrestreamAt(i, j)
                        Exit Select
                    Case CONST_IOLET
                        dealIOletPrestreamAt(i, j)
                        Exit Select
                    Case CONST_INTERFACE
                        Exit Select
                End Select
            Next
        Next
        dealItfPrestream()
    End Sub

    Private Sub dealWallPrestreamAt(ByVal i As Int32, ByVal j As Int32)
        Dim k As Byte
        Dim wcell As WallCell = field(i, j)
        For k = 1 To q
            If isIn(i + ei(k, 0), j + ei(k, 1)) Then
                If flag(i + ei(k, 0), j + ei(k, 1)) >= CONST_INTERFACE Then
                    wcell.f(k) = field(i + ei(k, 0), j + ei(k, 1)).f(oppk(k))
                End If
            End If
        Next
    End Sub

    Private Sub dealIOletPrestreamAt(ByVal i As Int32, ByVal j As Int32)
        Dim k As Byte
        Dim f() As Double
        Dim iocell As IOletCell = field(i, j)
        If iocell.flag2 = CONST_EMPTY Then
            For k = 1 To q
                If isIn(i + ei(k, 0), j + ei(k, 1)) Then
                    If flag(i + ei(k, 0), j + ei(k, 1)) = CONST_FLUID Then
                        f = field(i + ei(k, 0), j + ei(k, 1)).f.Clone
                        field(i + ei(k, 0), j + ei(k, 1)) = New InterfaceCell
                        field(i + ei(k, 0), j + ei(k, 1)).f = f.Clone
                        flag(i, j) = CONST_INTERFACE
                        itfArray.Add(New Int32() {i, j})
                    End If
                End If
            Next
        End If
    End Sub

    Private Sub dealItfPrestream()
        Dim k As Byte
        Dim a, b As Int32
        Dim location(1) As Int32
        Dim icell As InterfaceCell
        Dim iocell As IOletCell
        For Each location In itfArray
            icell = field(location(0), location(1))
            icell.mass = mass(location(0), location(1))
            If hasCellType(location(0), location(1), CONST_WALL) Then
                icell.normalDirection(0) = u(location(0), location(1))
                icell.normalDirection(1) = v(location(0), location(1))
            Else
                icell.normalDirection = getNormalDirection(location(0), location(1))
            End If

            For k = 1 To q
                a = location(0) + ei(k, 0)
                b = location(1) + ei(k, 1)
                Select Case flag(a, b)
                    Case CONST_EMPTY
                        icell.tf(oppk(k)) = getfe2(location(0), location(1), k) + getfe2(location(0), location(1), oppk(k)) - icell.f(k)
                        Exit Select
                    Case CONST_FLUID
                        icell.tf(oppk(k)) = field(a, b).f(oppk(k))
                        icell.mass += field(a, b).f(oppk(k)) - icell.f(k)
                        Exit Select
                    Case CONST_INTERFACE
                        icell.tf(oppk(k)) = field(a, b).f(oppk(k))
                        dtable(location(0), location(1), a, b, k, icell)
                        Exit Select
                    Case CONST_WALL
                        icell.tf(oppk(k)) = field(a, b).f(oppk(k))
                        Exit Select
                    Case CONST_IOLET
                        iocell = field(a, b)
                        If iocell.flag2 = CONST_EMPTY Then
                            icell.tf(oppk(k)) = getfe2(location(0), location(1), k) + getfe2(location(0), location(1), oppk(k)) - icell.f(k)
                        ElseIf iocell.flag2 = CONST_FLUID Then
                            icell.tf(oppk(k)) = field(a, b).f(oppk(k))
                            icell.mass += field(a, b).f(oppk(k)) - icell.f(k)
                        Else
                            icell.tf(oppk(k)) = field(a, b).f(oppk(k))
                            icell.mass += (field(a, b).f(oppk(k) - icell.f(k))) * (eps(a, b) + eps(location(0), location(1))) / 2
                        End If
                        Exit Select

                End Select
            Next
        Next
    End Sub

    Private Sub dtable(ByVal i As Int32, ByVal j As Int32, ByVal a As Int32, ByVal b As Int32, ByVal k As Byte, ByRef icell As InterfaceCell)
        st(0) = getState(i, j)
        st(1) = getState(a, b)
        If st(0) = 1 Then
            Select Case st(1)
                Case 1
                    icell.mass += (field(a, b).f(oppk(k)) - icell.f(k)) * (eps(a, b) + eps(i, j)) / 2
                Case 0
                    icell.mass += field(a, b).f(oppk(k)) * (eps(a, b) + eps(i, j)) / 2
                Case 2
                    icell.mass += -1 * icell.f(k) * (eps(a, b) + eps(i, j)) / 2
            End Select
        ElseIf st(0) = 0 Then
            Select Case st(1)
                Case 1
                    icell.mass += -1 * icell.f(k) * (eps(a, b) + eps(i, j)) / 2
                Case 0
                    icell.mass += (field(a, b).f(oppk(k)) - icell.f(k)) * (eps(a, b) + eps(i, j)) / 2
                Case 2
                    icell.mass += -1 * icell.f(k) * (eps(a, b) + eps(i, j)) / 2
            End Select
        Else
            Select Case st(1)
                Case 1
                    icell.mass += field(a, b).f(oppk(k)) * (eps(a, b) + eps(i, j)) / 2
                Case 0
                    icell.mass += field(a, b).f(oppk(k)) * (eps(a, b) + eps(i, j)) / 2
                Case 2
                    icell.mass += (field(a, b).f(oppk(k)) - icell.f(k)) * (eps(a, b) + eps(i, j)) / 2
            End Select
        End If
    End Sub

    Private Sub poststream()
        Dim location(1) As Int32
        Dim itfcell As InterfaceCell
        For Each location In itfArray
            itfcell = field(location(0), location(1))
            mass(location(0), location(1)) = itfcell.mass
            itfcell.reNewf()
        Next
    End Sub

    Private Sub reconstruction()
        Dim location(1) As Int32
        fillArray.Clear()
        emptyArray.Clear()
        cnt = 0
        Do While cnt < itfArray.Count
            location = itfArray.Item(cnt)
            If mass(location(0), location(1)) > (1 + kapa) * rho(location(0), location(1)) Then
                fillArray.Add(location.Clone)
                itfArray.RemoveAt(cnt)
            ElseIf mass(location(0), location(1)) < (0 - kapa) * rho(location(0), location(1)) Then
                emptyArray.Add(location.Clone)
                flag(location(0), location(1)) = CONST_MEDIUM
                itfArray.RemoveAt(cnt)
            Else
                cnt += 1
            End If
        Loop
        dfill()
        dempty()
        dcolor()
        ddel()
        deps()
    End Sub

    Private Sub dfill()
        Dim a, b As Int32
        Dim location(1) As Int32
        Dim k As Byte
        For Each location In fillArray
            For k = 1 To q
                a = location(0) + ei(k, 0)
                b = location(1) + ei(k, 1)
                If flag(a, b) = CONST_EMPTY Then
                    emptyToItf(a, b)
                End If