Option Explicit On
Option Strict On
Public Class Form2
    Private KBD As Integer = 100   '最大境界数
    Private KTJ As Integer = 10000 '最大節点数
    Private KCM As Integer = 100   '節点に集まる最大要素数
    Private KTE As Integer
    Private LTE As Integer = 100  '多角形の最大辺数
    Private pi As Double = Math.PI

    Private Sub Form2_Load(ByVal sender As System.Object, ByVal e As System.EventArgs) Handles MyBase.Load
        Call MAIN()

        Dim fm3 As New Form3()
        fm3.Show()
        Me.Close()
    End Sub


    Private Sub MAIN()
        Dim ibex() As Integer     '個々の外部境界上の節点数
        Dim ibin() As Integer     '個々の内部境界上の節点数
        Dim ibno(,) As Integer    '外部境界上の節点番号
        Dim index() As Integer    '作業領域
        Dim px() As Double        '節点のx座標
        Dim py() As Double        '節点のy座標
        Dim mtj(,) As Integer     '要素-節点関係
        Dim jac(,) As Integer     '要素-隣接要素関係
        Dim idm() As Integer      '要素-領域関係
        Dim ifix() As Integer     '固定節点（=1）・可動節点（=0）の判別
        Dim delx() As Double      '各部分領域の格子間隔

        ReDim ibex(KBD)
        ReDim ibin(KBD)
        ReDim ibno(KBD, KTJ)
        ReDim index(KBD + 1)
        ReDim px(KTJ + 3)
        ReDim py(KTJ + 3)
        ReDim mtj(2 * KTJ + 1, 3)
        ReDim jac(2 * KTJ + 1, 3)
        ReDim idm(2 * KTJ + 1)
        ReDim ifix(KTJ)
        ReDim delx(KBD)

        Dim NEX As Integer '外部境界数
        Dim NIN As Integer '内部境界数
        Dim NOB As Integer '境界上の節点総数
        Dim NIB As Integer '境界内部の節点総数
        Dim NODE As Integer '節点総数
        Dim NELM As Integer '要素総数


        Call TRINPUT(NEX, NIN, ibex, ibin, ibno, NOB, NIB, px, py, delx)
        Call TRMODEL(NEX, NIN, ibex, ibin, ibno, NOB, NIB, index, NODE, px, py, NELM, mtj, jac, idm, ifix, delx)
        Call TRDATA(NEX, NODE, NELM, mtj, jac, idm, px, py, ifix)

    End Sub


    Private Sub TRINPUT(ByRef NEX As Integer, ByRef NIN As Integer, ByRef ibex() As Integer, ByRef ibin() As Integer, _
                        ByRef ibno(,) As Integer, ByRef NOB As Integer, ByRef NIB As Integer, _
                        ByRef px() As Double, ByRef py() As Double, ByRef delx() As Double)

        Dim sr As System.IO.StreamReader
        Dim dat As String
        Dim sbuf() As String
        Dim delim() As Char = {","c}
        Dim fname1 As String
        Dim i As Integer
        Dim j As Integer

        fname1 = My.Forms.Form1.Label1.Text
        sr = New System.IO.StreamReader(fname1, System.Text.Encoding.Default)
        dat = sr.ReadLine() : sbuf = dat.Split(delim)
        NEX = CInt(sbuf(0))
        NIN = CInt(sbuf(1))
        If 0 < NEX Then
            For i = 1 To NEX
                dat = sr.ReadLine() : sbuf = dat.Split(delim)
                ibex(i) = CInt(sbuf(0))
                delx(i) = CDbl(sbuf(1))
            Next i
        End If
        If 0 < NIN Then
            For i = 1 To NIN
                dat = sr.ReadLine() : sbuf = dat.Split(delim)
                ibin(i) = CInt(sbuf(0))
            Next i
        End If
        For i = 1 To NEX
            dat = sr.ReadLine() : sbuf = dat.Split(delim)
            For j = 1 To ibex(i)
                ibno(i, j) = CInt(sbuf(j - 1))
        Next j, i
        dat = sr.ReadLine() : sbuf = dat.Split(delim)
        NOB = CInt(sbuf(0))
        NIB = CInt(sbuf(1))
        For i = 1 To NOB + NIB
            dat = sr.ReadLine() : sbuf = dat.Split(delim)
            px(i) = CDbl(sbuf(0))
            py(i) = CDbl(sbuf(1))
        Next i
        sr.Close()

    End Sub


    Private Sub TRMODEL(ByRef NEX As Integer, ByRef NIN As Integer, ByRef ibex() As Integer, ByRef ibin() As Integer, _
                        ByRef ibno(,) As Integer, ByRef NOB As Integer, ByRef NIB As Integer, _
                        ByRef index() As Integer, ByRef NODE As Integer, ByRef px() As Double, ByRef py() As Double, _
                        ByRef NELM As Integer, ByRef mtj(,) As Integer, ByRef jac(,) As Integer, _
                        ByRef idm() As Integer, ByRef ifix() As Integer, ByRef delx() As Double)

        Dim i As Integer
        Dim j As Integer
        Dim ntp As Integer
        Dim xmin As Double
        Dim xmax As Double
        Dim ymin As Double
        Dim ymax As Double
        Dim rax As Double
        Dim ray As Double
        Dim dmax As Double
        Dim ia As Integer
        Dim ib As Integer
        Dim ic As Integer

        Dim jnb() As Integer
        Dim nei(,) As Integer
        Dim angl() As Double
        Dim nedg() As Integer
        Dim kstack() As Integer
        Dim list() As Integer
        Dim istack() As Integer
        Dim kv() As Integer
        Dim ibr() As Integer
        Dim iadres() As Integer
        Dim map() As Integer

        ReDim jnb(KTJ)
        ReDim nei(KTJ, KCM)
        ReDim ifix(KTJ)
        ReDim angl(2 * KTJ + 1)
        ReDim nedg(2 * KTJ + 1)
        ReDim kstack(2 * KTJ + 1)
        ReDim list(KTJ)
        ReDim istack(KTJ)
        ReDim kv(2 * KTJ + 1)
        ReDim ibr(2 * KTJ + 1)
        ReDim iadres(KTJ + 3)
        ReDim map(2 * KTJ + 1)

        '初期化
        For i = 1 To KBD + 1
            index(i) = 0
        Next i
        KTE = 2 * KTJ + 1
        NODE = 0
        NELM = 0
        For i = 1 To KTE
            idm(i) = 0
            For j = 1 To 3
                mtj(i, j) = 0
                jac(i, j) = 0
        Next j, i
        For i = 1 To KTJ
            jnb(i) = 0
            For j = 1 To KCM
                nei(i, j) = 0
        Next j, i
        For i = 1 To KTJ
            ifix(i) = 0
        Next i

        ntp = NOB + NIB
        xmin = px(1)
        xmax = xmin
        ymin = py(1)
        ymax = ymin
        For i = 2 To ntp
            If px(i) < xmin Then xmin = px(i)
            If xmax < px(i) Then xmax = px(i)
            If py(i) < ymin Then ymin = py(i)
            If ymax < py(i) Then ymax = py(i)
        Next i
        rax = xmax - xmin
        ray = ymax - ymin
        dmax = rax
        If dmax < ray Then dmax = ray
        For i = 1 To ntp
            px(i) = (px(i) - xmin) / dmax
            py(i) = (py(i) - ymin) / dmax
        Next i

        NELM = 1
        ia = KTJ + 1
        ib = KTJ + 2
        ic = KTJ + 3
        mtj(1, 1) = ia
        mtj(1, 2) = ib
        mtj(1, 3) = ic
        jac(1, 1) = 0
        jac(1, 2) = 0
        jac(1, 3) = 0

        px(ia) = -1.23
        py(ia) = -0.5
        px(ib) = 2.23
        py(ib) = -0.5
        px(ic) = 0.5
        py(ic) = 2.5

        Call ROUGH(NEX, NIN, ibex, ibin, ibno, ntp, NIB, NODE, px, py, jnb, nei, NELM, mtj, jac, idm, _
                       list, iadres, istack, kv, ibr, map)

        For i = 1 To NODE
            ifix(i) = 1
        Next i

        Call TRFINE(NEX, xmin, ymin, dmax, NODE, px, py, jnb, nei, NELM, mtj, jac, idm, iadres, istack, delx)
        If NELM <> 2 * NODE + 1 Then
            'ERRORメッセージ＆実行中止
            MessageBox.Show("ERROR IN SUB TRMODEL", "Error終了")
            Me.Close()
        End If

        Call REMOVE(NEX, index, NELM, mtj, jac, idm)
        Call CHECK(NELM, mtj, jac)

        For i = 1 To NODE
            px(i) = px(i) * dmax + xmin
            py(i) = py(i) * dmax + ymin
        Next i

    End Sub


    Private Sub ROUGH(ByRef NEX As Integer, ByRef NIN As Integer, ByRef ibex() As Integer, ByRef ibin() As Integer, _
                      ByRef ibno(,) As Integer, ByRef ntp As Integer, ByRef NIB As Integer, _
                      ByRef NODE As Integer, ByRef px() As Double, ByRef py() As Double, _
                      ByRef jnb() As Integer, ByRef nei(,) As Integer, _
                      ByRef NELM As Integer, ByRef mtj(,) As Integer, ByRef jac(,) As Integer, ByRef idm() As Integer, _
                      ByRef list() As Integer, ByRef iadres() As Integer, ByRef istack() As Integer, ByRef kv() As Integer, _
                      ByRef ibr() As Integer, ByRef map() As Integer)

        Dim i As Integer
        Dim j As Integer
        Dim k As Integer
        Dim nb As Integer
        Dim np As Integer
        Dim ma As Integer
        Dim mb As Integer
        Dim mc As Integer
        Dim ms As Integer
        Dim mp As Integer
        Dim js As Integer
        Dim jz As Integer
        Dim xs As Double
        Dim ys As Double
        Dim loc As Integer
        Dim xa As Double
        Dim ya As Double
        Dim it As Integer

        For i = 1 To NEX
            jz = 0
            nb = i
            np = ibex(nb)
            For j = 1 To np
                list(j) = ibno(nb, j)
            Next j
            Call BOUGEN(jz, np, list, ntp, px, py, jnb, nei, NELM, mtj, jac, idm, iadres, istack, kv, ibr, map, NODE)
            For k = 1 To NELM
                If idm(k) = jz Then
                    ma = iadres(mtj(k, 1))
                    mb = iadres(mtj(k, 2))
                    mc = iadres(mtj(k, 3))
                    ms = ma * mb * mc
                    mp = (mb - ma) * (mc - mb) * (ma - mc)
                    If (ms <> 0) And (0 < mp) Then
                        idm(k) = nb
                    End If
                End If
            Next k
        Next i

        For i = 1 To NIN
            nb = i
            np = ibin(nb)
            For j = 1 To np
                list(j) = NODE + j
            Next j
            js = list(1)
            xs = px(js)
            ys = py(js)
            Call LOCATE(xs, ys, px, py, mtj, jac, NELM, loc)
            jz = idm(loc)
            Call BOUGEN(jz, np, list, ntp, px, py, jnb, nei, NELM, mtj, jac, idm, iadres, istack, kv, ibr, map, NODE)

            For k = 1 To NELM
                If idm(k) = jz Then
                    ma = iadres(mtj(k, 1))
                    mb = iadres(mtj(k, 2))
                    mc = iadres(mtj(k, 3))
                    ms = ma * mb * mc
                    mp = (mb - ma) * (mc - mb) * (ma - mc)
                    If (ms <> 0) And (mp < 0) Then
                        idm(k) = 0
                    End If
                End If
            Next k
        Next i

        For i = 1 To KTJ + 3
            iadres(i) = 0
        Next i
        For i = 1 To NIB
            NODE = NODE + 1
            xa = px(NODE)
            ya = py(NODE)
            Call LOCATE(xa, ya, px, py, mtj, jac, NELM, it)
            jz = idm(it)
            Call DELAUN(jz, NODE, NODE, ntp, px, py, jnb, nei, NELM, mtj, jac, idm, iadres, istack)
        Next i

        If NODE <> ntp Then
            'ERRORメッセージ＆実行中止
            MessageBox.Show("ERROR IN SUB ROUGH", "Error終了")
            Me.Close()
        End If

    End Sub


    Private Sub BOUGEN(ByRef jz As Integer, ByRef np As Integer, ByRef list() As Integer, ByRef ntp As Integer, _
                       ByRef px() As Double, ByRef py() As Double, ByRef jnb() As Integer, ByRef nei(,) As Integer, _
                       ByRef NELM As Integer, ByRef mtj(,) As Integer, ByRef jac(,) As Integer, ByRef idm() As Integer, _
                       ByRef iadres() As Integer, ByRef istack() As Integer, ByRef kv() As Integer, ByRef ibr() As Integer, _
                       ByRef map() As Integer, ByRef NODE As Integer)

        Dim i As Integer
        Dim inp As Integer
        Dim js As Integer
        Dim ip As Integer
        Dim iq As Integer
        Dim iv As Integer

        inp = 0
        For i = 1 To KTJ + 3
            iadres(i) = 0
        Next i
        For i = 1 To np
            iadres(list(i)) = i
        Next i
        js = list(1)
        If NODE < js Then
            inp = inp + 1
            Call DELAUN(jz, js, js, ntp, px, py, jnb, nei, NELM, mtj, jac, idm, iadres, istack)
        End If

        For i = 1 To np
            ip = list((i Mod np) + 1)
            iq = list(i)
            If (NODE < ip) And (i <> np) Then
                inp = inp + 1
                Call DELAUN(jz, ip, ip, ntp, px, py, jnb, nei, NELM, mtj, jac, idm, iadres, istack)
            End If
            For j = 1 To jnb(iq)
                For k = 1 To 3
                    If mtj(nei(iq, j), k) = ip Then GoTo bougen30
                Next k
            Next j
            Call SEARCH(jz, ip, iq, jnb, nei, NELM, mtj, jac, idm, istack, iv, kv, iadres, ibr)
            Call POLY(iq, ip, iv, kv, px, py, NELM, mtj, jac, jnb, nei, map)
bougen30:
        Next i
        NODE = NODE + inp

    End Sub


    Private Sub DELAUN(ByRef jz As Integer, ByRef js As Integer, ByRef jg As Integer, ByRef ntp As Integer, ByRef px() As Double, ByRef py() As Double, _
                       ByRef jnb() As Integer, ByRef nei(,) As Integer, ByRef NELM As Integer, ByRef mtj(,) As Integer, _
                       ByRef jac(,) As Integer, ByRef idm() As Integer, ByRef iadres() As Integer, ByRef istack() As Integer)
        Dim itop As Integer
        Dim maxstk As Integer
        Dim ip As Integer
        Dim it As Integer
        Dim xp As Double
        Dim yp As Double
        Dim ia As Integer
        Dim ib As Integer
        Dim ic As Integer
        Dim iv1 As Integer
        Dim iv2 As Integer
        Dim iv3 As Integer
        Dim idf As Integer
        Dim ms As Integer
        Dim iedge As Integer
        Dim il As Integer
        Dim ir As Integer
        Dim iera As Integer
        Dim ierb As Integer
        Dim ierl As Integer
        Dim iswap As Integer

        itop = 0
        maxstk = ntp
        For i = js To jg
            ip = i
            xp = px(ip)
            yp = py(ip)
            Call LOCATE(xp, yp, px, py, mtj, jac, NELM, it)
            If idm(it) <> jz Then
                'ERRORメッセージ＆実行中止
                MessageBox.Show("ERROR IN SUB DELAUN", "Error終了")
                Me.Close()
            End If

            ia = jac(it, 1)
            ib = jac(it, 2)
            ic = jac(it, 3)
            iv1 = mtj(it, 1)
            iv2 = mtj(it, 2)
            iv3 = mtj(it, 3)
            mtj(it, 1) = ip
            mtj(it, 2) = iv1
            mtj(it, 3) = iv2
            jac(it, 1) = NELM + 2
            jac(it, 2) = ia
            jac(it, 3) = NELM + 1

            NELM = NELM + 1
            idm(NELM) = jz
            mtj(NELM, 1) = ip
            mtj(NELM, 2) = iv2
            mtj(NELM, 3) = iv3
            jac(NELM, 1) = it
            jac(NELM, 2) = ib
            jac(NELM, 3) = NELM + 1
            NELM = NELM + 1
            idm(NELM) = jz
            mtj(NELM, 1) = ip
            mtj(NELM, 2) = iv3
            mtj(NELM, 3) = iv1
            jac(NELM, 1) = NELM - 1
            jac(NELM, 2) = ic
            jac(NELM, 3) = it

            Call INCR(iv1, NELM, jnb, nei)
            Call INCR(iv2, NELM - 1, jnb, nei)
            If iv3 <= KTJ Then
                nei(iv3, NEIBOR(iv3, it, jnb, nei)) = NELM - 1
            End If
            Call INCR(iv3, NELM, jnb, nei)
            jnb(ip) = 3
            nei(ip, 1) = it
            nei(ip, 2) = NELM - 1
            nei(ip, 3) = NELM
            If ia <> 0 Then
                ms = iadres(iv1) * iadres(iv2)
                idf = Math.Abs(iadres(iv1) - iadres(iv2))
                If (idm(ia) = jz) And ((ms = 0) Or (idf <> 1)) Then
                    itop = itop + 1
                    istack(itop) = IPUSH(it, maxstk, itop, istack)
                End If
            End If
            If ib <> 0 Then
                Call EDGE(ib, it, jac, iedge)
                jac(ib, iedge) = NELM - 1
                ms = iadres(iv2) * iadres(iv3)
                idf = Math.Abs(iadres(iv2) - iadres(iv3))
                If (idm(ib) = jz) And ((ms = 0) Or (idf <> 1)) Then
                    itop = itop + 1
                    istack(itop) = IPUSH(NELM - 1, maxstk, itop, istack)
                End If
            End If
            If ic <> 0 Then
                Call EDGE(ic, it, jac, iedge)
                jac(ic, iedge) = NELM
                ms = iadres(iv3) * iadres(iv1)
                idf = Math.Abs(iadres(iv3) - iadres(iv1))
                If (idm(ic) = jz) And ((ms = 0) Or (idf <> 1)) Then
                    itop = itop + 1
                    istack(itop) = IPUSH(NELM, maxstk, itop, istack)
                End If
            End If

delaun50:
            If 0 < itop Then
                il = istack(itop)
                itop = itop - 1
                ir = jac(il, 2)
                Call EDGE(ir, il, jac, ierl)
                iera = (ierl Mod 3) + 1
                ierb = (iera Mod 3) + 1
                iv1 = mtj(ir, ierl)
                iv2 = mtj(ir, iera)
                iv3 = mtj(ir, ierb)
                Call SWAP(px(iv1), py(iv1), px(iv2), py(iv2), px(iv3), py(iv3), xp, yp, iswap)
                If iswap = 1 Then
                    ia = jac(ir, iera)
                    ib = jac(ir, ierb)
                    ic = jac(il, 3)
                    mtj(il, 3) = iv3
                    jac(il, 2) = ia
                    jac(il, 3) = ir
                    mtj(ir, 1) = ip
                    mtj(ir, 2) = iv3
                    mtj(ir, 3) = iv1
                    jac(ir, 1) = il
                    jac(ir, 2) = ib
                    jac(ir, 3) = ic
                    Call DECR(iv1, il, jnb, nei)
                    Call DECR(iv2, ir, jnb, nei)
                    Call INCR(ip, ir, jnb, nei)
                    Call INCR(iv3, il, jnb, nei)

                    If ia <> 0 Then
                        Call EDGE(ia, ir, jac, iedge)
                        jac(ia, iedge) = il
                        ms = iadres(iv2) * iadres(iv3)
                        idf = Math.Abs(iadres(iv2) - iadres(iv3))
                        If (idm(ia) = jz) And ((ms = 0) Or (idf <> 1)) Then
                            itop = itop + 1
                            istack(itop) = IPUSH(il, maxstk, itop, istack)
                        End If
                    End If
                    If ib <> 0 Then
                        ms = iadres(iv3) * iadres(iv1)
                        idf = Math.Abs(iadres(iv3) - iadres(iv1))
                        If (idm(ib) = jz) And ((ms = 0) Or (idf <> 1)) Then
                            itop = itop + 1
                            istack(itop) = IPUSH(ir, maxstk, itop, istack)
                        End If
                    End If
                    If ic <> 0 Then
                        Call EDGE(ic, il, jac, iedge)
                        jac(ic, iedge) = ir
                    End If
                End If
                GoTo delaun50
            End If
        Next i

    End Sub


    Private Sub SEARCH(ByRef iz As Integer, ByRef ip As Integer, ByRef iq As Integer, ByRef jnb() As Integer, _
                       ByRef nei(,) As Integer, ByRef NELM As Integer, ByRef mtj(,) As Integer, ByRef jac(,) As Integer, _
                       ByRef idm() As Integer, ByRef istack() As Integer, ByRef iv As Integer, _
                       ByRef kv() As Integer, ByRef iadres() As Integer, ByRef ibr() As Integer)
        Dim i As Integer
        Dim j As Integer
        Dim ia As Integer
        Dim ib As Integer
        Dim ja As Integer
        Dim jb As Integer
        Dim mstk As Integer
        Dim nbr As Integer
        Dim idf As Integer
        Dim ms As Integer
        Dim ielm As Integer
        Dim jelm As Integer
        Dim kelm As Integer
        Dim min As Integer
        Dim jr As Integer

        iv = 0
        mstk = 0
        nbr = 0
        For i = 1 To NELM
            kv(i) = 0
            ibr(i) = 0
        Next i
        For i = 1 To KTJ
            istack(i) = 0
        Next i
        For i = 1 To jnb(iq)
            nbr = nbr + 1
            ibr(nei(iq, i)) = nbr
        Next i
        For i = 1 To jnb(iq)
            ielm = nei(iq, i)
            j = IVERT(ielm, iq, mtj)
            ja = (j Mod 3) + 1
            jb = (ja Mod 3) + 1
            ia = mtj(ielm, ja)
            ib = mtj(ielm, jb)
            idf = Math.Abs(iadres(ia) - iadres(ib))
            ms = iadres(ia) * iadres(ib)
            If (idf = 1) And (ms <> 0) Then GoTo search40
            jelm = jac(ielm, ja)
            If jelm = 0 Then GoTo search40
            If idm(jelm) <> iz Then GoTo search40
            nbr = nbr + 1
            ibr(jelm) = nbr
            If mtj(jelm, 1) = ip Then GoTo search80
            If mtj(jelm, 2) = ip Then GoTo search80
            If mtj(jelm, 3) = ip Then GoTo search80
            mstk = mstk + 1
            istack(mstk) = jelm
search40:
        Next i

search50:
        If mstk = 0 Then
            'ERRORメッセージ＆実行中止
            MessageBox.Show("ERROR IN SUB SEARCH", "Error終了")
            Me.Close()
        End If

        ielm = istack(1)
        mstk = mstk - 1
        For i = 1 To mstk
            istack(i) = istack(i + 1)
        Next i
        istack(mstk + 1) = 0
        For j = 1 To 3
            ja = (j Mod 3) + 1
            ia = mtj(ielm, j)
            ib = mtj(ielm, ja)
            idf = Math.Abs(iadres(ia) - iadres(ib))
            ms = iadres(ia) * iadres(ib)
            If (idf = 1) And (ms <> 0) Then GoTo search70
            jelm = jac(ielm, j)
            If jelm = 0 Then GoTo search70
            If ibr(jelm) <> 0 Then GoTo search70
            If idm(jelm) <> iz Then GoTo search70
            nbr = nbr + 1
            ibr(jelm) = nbr
            If mtj(jelm, 1) = ip Then GoTo search80
            If mtj(jelm, 2) = ip Then GoTo search80
            If mtj(jelm, 3) = ip Then GoTo search80
            mstk = mstk + 1
            istack(mstk) = jelm
search70:
        Next j
        GoTo search50

search80:
        iv = iv + 1
        kv(iv) = jelm
        If ibr(jelm) <= jnb(iq) Then GoTo search100
        min = KTE + 1
        For j = 1 To 3
            jr = jac(jelm, j)
            If jr = 0 Then GoTo search90
            If ibr(jr) = 0 Then GoTo search90
            ja = (j Mod 3) + 1
            ia = mtj(jelm, j)
            ib = mtj(jelm, ja)
            idf = Math.Abs(iadres(ia) - iadres(ib))
            ms = iadres(ia) * iadres(ib)
            If (idf = 1) And (ms <> 0) Then GoTo search90
            If ibr(jr) < min Then
                kelm = jr
                min = ibr(jr)
            End If
search90:
        Next j
        If min = KTE + 1 Then
            'ERRORメッセージ＆実行中止
            MessageBox.Show("ERROR IN SUB SEARCH", "Error終了")
            Me.Close()
        End If

        jelm = kelm
        GoTo search80
search100:

    End Sub


    Private Sub POLY(ByRef iq As Integer, ByRef ip As Integer, ByRef iv As Integer, ByRef kv() As Integer, _
                     ByRef px() As Double, ByRef py() As Double, ByRef NELM As Integer, ByRef mtj(,) As Integer, _
                     ByRef jac(,) As Integer, ByRef jnb() As Integer, ByRef nei(,) As Integer, ByRef map() As Integer)
        Dim i As Integer
        Dim j As Integer
        Dim k As Integer
        Dim l As Integer
        Dim ia As Integer
        Dim ja As Integer
        Dim ir As Integer
        Dim jr As Integer
        Dim iv1 As Integer
        Dim iv2 As Integer
        Dim iedge As Integer
        Dim npa As Integer
        Dim npb As Integer
        Dim nta As Integer
        Dim ntb As Integer
        Dim ielm As Integer
        Dim jelm As Integer
        Dim kelm As Integer
        Dim ivx As Integer
        Dim ips As Integer
        Dim ipg As Integer
        Dim iva As Integer
        Dim ivb As Integer

        Dim iena(,) As Integer
        Dim ienb(,) As Integer
        Dim jeea(,) As Integer
        Dim jeeb(,) As Integer
        Dim nsra() As Integer
        Dim nsrb() As Integer
        Dim ihen(,) As Integer
        Dim jhen() As Integer
        Dim iad() As Integer
        Dim jstack() As Integer

        ReDim iena(LTE, 3)
        ReDim ienb(LTE, 3)
        ReDim jeea(LTE, 3)
        ReDim jeeb(LTE, 3)
        ReDim nsra(LTE + 2)
        ReDim nsrb(LTE + 2)
        ReDim ihen(2 * LTE + 1, 2)
        ReDim jhen(2 * LTE + 1)
        ReDim iad(2 * LTE + 1)
        ReDim jstack(LTE + 2)

        If iv = 2 Then
            Call EDGE(kv(1), kv(2), jac, ia)
            Call EDGE(kv(2), kv(1), jac, ja)

            ir = jac(kv(1), (ia Mod 3) + 1)
            jr = jac(kv(2), (ja Mod 3) + 1)

            mtj(kv(1), (ia Mod 3) + 1) = iq
            jac(kv(1), ia) = jr
            jac(kv(1), (ja Mod 3) + 1) = kv(2)

            mtj(kv(2), (ja Mod 3) + 1) = ip
            jac(kv(2), ja) = ir
            jac(kv(2), (ja Mod 3) + 1) = kv(1)

            If ir <> 0 Then
                Call EDGE(ir, kv(1), jac, iedge)
                jac(ir, iedge) = kv(2)
            End If
            If jr <> 0 Then
                Call EDGE(jr, kv(2), jac, iedge)
                jac(jr, iedge) = kv(1)
            End If

            iv1 = mtj(kv(1), ia)
            iv2 = mtj(kv(2), ja)
            Call DECR(iv1, kv(2), jnb, nei)
            Call DECR(iv2, kv(1), jnb, nei)
            Call INCR(iq, kv(1), jnb, nei)
            Call INCR(ip, kv(2), jnb, nei)
        Else
            npa = 0
            npb = 0
            For i = 1 To LTE + 2
                nsra(i) = 0
                nsrb(i) = 0
            Next i
            For i = 1 To NELM
                map(i) = 0
            Next i
            For i = 1 To iv
                map(kv(1)) = 1
            Next i
            For i = 1 To iv
                ielm = kv(i)
                For j = 1 To 3
                    ivx = mtj(ielm, j)
                    Call DECR(ivx, ielm, jnb, nei)
                Next j
            Next i

            Call PICK(iq, ip, iv, kv, mtj, jac, map, npa, npb, nsra, nsrb)
            Call SUBDIV(npa, nsra, px, py, nta, iena, jeea, ihen, jhen, iad, jstack)
            Call SUBDIV(npb, nsrb, px, py, ntb, ienb, jeeb, ihen, jhen, iad, jstack)

            If iv <> nta + ntb Then
                'ERRORメッセージ＆実行中止
                MessageBox.Show("ERROR IN SUB POLY", "Error終了")
                Me.Close()
            End If

            For i = 1 To iv
                For j = 1 To 3
                    jac(kv(i), j) = 0
                Next j
            Next i
            For i = 1 To nta
                ielm = kv(i)
                For j = 1 To 3
                    mtj(ielm, j) = iena(i, j)
                    If jeea(i, j) <> 0 Then
                        jac(ielm, j) = kv(jeea(i, j))
                    End If
                Next j
            Next i
            For i = 1 To ntb
                ielm = kv(nta + i)
                For j = 1 To 3
                    mtj(ielm, j) = ienb(i, j)
                    If jeeb(i, j) <> 0 Then
                        jac(ielm, j) = kv(nta + jeeb(i, j))
                    End If
                Next j
            Next i
            For i = 1 To iv
                ielm = kv(i)
                For j = 1 To 3
                    ivx = mtj(ielm, j)
                    Call INCR(ivx, ielm, jnb, nei)
                Next j
            Next i

            For i = 1 To iv
                ielm = kv(i)
                For j = 1 To 3
                    jelm = jac(ielm, j)
                    If jelm <> 0 Then GoTo poly130
                    ips = mtj(ielm, j)
                    ipg = mtj(ielm, (j Mod 3) + 1)
                    For k = 1 To jnb(ipg)
                        kelm = nei(ipg, k)
                        For l = 1 To 3
                            iva = mtj(kelm, l)
                            ivb = mtj(kelm, (l Mod 3) + 1)
                            If (iva = ipg) And (ivb = ips) Then
                                jac(ielm, j) = kelm
                                jac(kelm, l) = ielm
                                GoTo poly130
                            End If
                        Next l
                    Next k
poly130:
                Next j
            Next i
        End If

    End Sub


    Private Sub PICK(ByRef iq As Integer, ByRef ip As Integer, ByRef iv As Integer, ByRef kv() As Integer, _
                     ByRef mtj(,) As Integer, ByRef jac(,) As Integer, ByRef map() As Integer, _
                     ByRef npa As Integer, ByRef npb As Integer, ByRef nsra() As Integer, ByRef nsrb() As Integer)

        Dim i As Integer
        Dim ivx As Integer
        Dim jvx As Integer
        Dim jelm As Integer

        npa = 1
        nsra(npa) = ip
        ivx = IVERT(kv(1), ip, mtj)
        npa = npa + 1
        nsra(npa) = mtj(kv(1), (ivx Mod 3) + 1)
        For i = 2 To iv - 1
            jvx = IVERT(kv(i), nsra(npa), mtj)
            jelm = jac(kv(i), jvx)
            If jelm = 0 Then
                npa = npa + 1
                nsra(npa) = mtj(kv(i), (jvx Mod 3) + 1)
            ElseIf map(jelm) = 0 Then
                npa = npa + 1
                nsra(npa) = mtj(kv(i), (jvx Mod 3) + 1)
            End If
        Next i
        npa = npa + 1
        nsra(npa) = iq

        npb = 1
        nsrb(npb) = iq
        ivx = IVERT(kv(iv), iq, mtj)
        npb = npb + 1
        nsrb(npb) = mtj(kv(iv), (ivx Mod 3) + 1)
        For i = iv - 1 To 2 Step -1
            jvx = IVERT(kv(i), nsrb(npb), mtj)
            jelm = jac(kv(i), jvx)
            If jelm = 0 Then
                npb = npb + 1
                nsrb(npb) = mtj(kv(i), (jvx Mod 3) + 1)
            ElseIf map(jelm) = 0 Then
                npb = npb + 1
                nsrb(npb) = mtj(kv(i), (jvx Mod 3) + 1)
            End If
        Next i
        npb = npb + 1
        nsrb(npb) = ip

        If iv <> npa + npb - 4 Then
            'ERRORメッセージ＆実行中止
            MessageBox.Show("ERROR IN SUB PICK", "Error終了")
            Me.Close()
        End If

    End Sub


    Private Sub SUBDIV(ByRef npl As Integer, ByRef nsr() As Integer, ByRef px() As Double, ByRef py() As Double, _
                       ByRef nte As Integer, ByRef ien(,) As Integer, ByRef jee(,) As Integer, _
                       ByRef ihen(,) As Integer, ByRef jhen() As Integer, ByRef iad() As Integer, ByRef jstack() As Integer)
        Dim i As Integer
        Dim j As Integer
        Dim k As Integer
        Dim nnpl As Integer
        Dim nbs As Integer
        Dim ia As Integer
        Dim ib As Integer
        Dim ic As Integer
        Dim ix As Integer
        Dim xa As Double
        Dim ya As Double
        Dim xb As Double
        Dim yb As Double
        Dim see As Double

        nte = 0
        nnpl = npl
        For i = 1 To LTE
            For j = 1 To 3
                ien(i, j) = 0
                jee(i, j) = 0
            Next j
        Next i

subdiv20:
        If 3 <= nnpl Then
            nbs = 1
subdiv30:
            If nbs <= nnpl - 1 Then
                ia = nsr(nbs)
                ib = nsr(nbs + 1)
                ic = nsr(((nbs + 1) Mod nnpl) + 1)
                xa = px(ib) - px(ia)
                ya = py(ib) - py(ia)
                xb = px(ic) - px(ia)
                yb = py(ic) - py(ia)
                see = xa * yb - xb * ya
                If 0.000000000000001 < see Then
                    nte = nte + 1
                    ien(nte, 1) = ia
                    ien(nte, 2) = ib
                    ien(nte, 3) = ic
                    nnpl = nnpl - 1
                    For i = nbs + 1 To nnpl
                        nsr(i) = nsr(i + 1)
                    Next i
                End If
                nbs = nbs + 1
                GoTo subdiv30
            End If
            GoTo subdiv20
        End If

        ix = 0
        For i = 1 To 2 * LTE + 1
            ihen(i, 1) = 0
            ihen(i, 2) = 0
            jhen(i) = 0
            iad(i) = 0
        Next i

        For i = 1 To nte
            For j = 1 To 3
                ia = ien(i, j)
                ib = ien(i, (j Mod 3) + 1)
                For k = 1 To ix
                    If (ihen(k, 1) = ib) And (ihen(k, 2) = ia) Then
                        jee(i, j) = jhen(k)
                        jee(jhen(k), iad(k)) = i
                        GoTo subdiv70
                    End If
                Next k
                ix = ix + 1
                jhen(ix) = i
                iad(ix) = j
                ihen(ix, 1) = ia
                ihen(ix, 2) = ib
subdiv70:
            Next j
        Next i

        Call LAWSON(nte, ien, jee, npl, px, py, jstack)

    End Sub


    Private Sub LAWSON(ByRef nte As Integer, ByRef ien(,) As Integer, ByRef jee(,) As Integer, ByRef npl As Integer, _
                       ByRef px() As Double, ByRef py() As Double, ByRef jstack() As Integer)
        Dim i As Integer
        Dim j As Integer
        Dim itop As Integer
        Dim maxstk As Integer
        Dim ncount As Integer
        Dim ielm As Integer
        Dim il As Integer
        Dim ir As Integer
        Dim jl1 As Integer
        Dim jl2 As Integer
        Dim jl3 As Integer
        Dim jr1 As Integer
        Dim jr2 As Integer
        Dim jr3 As Integer
        Dim iv1 As Integer
        Dim iv2 As Integer
        Dim iv3 As Integer
        Dim iv4 As Integer
        Dim ia As Integer
        Dim ib As Integer
        Dim iswap As Integer
        Dim iedge As Integer
        Dim xx As Double
        Dim yy As Double

        itop = 0
        maxstk = npl
        ncount = 0
        For i = 1 To nte
            ielm = i
            itop = itop + 1
            jstack(itop) = IPUSH(ielm, maxstk, itop, jstack)
        Next i
lawson20:
        If 0 < itop Then
            ncount = ncount + 1
            If LTE < ncount Then
                'ERRORメッセージ＆実行中止
                MessageBox.Show("ERROR IN SUB LAWSON", "Error終了")
                Me.Close()
            End If
            il = jstack(itop)
            itop = itop - 1
            For j = 1 To 3
                jl1 = j
                jl2 = (jl1 Mod 3) + 1
                jl3 = (jl2 Mod 3) + 1
                ir = jee(il, jl1)
                If ir = 0 Then GoTo lawson30
                iv1 = ien(il, jl1)
                iv2 = ien(il, jl2)
                iv3 = ien(il, jl3)
                xx = px(ien(il, jl3))
                yy = py(ien(il, jl3))
                Call EDGE(ir, il, jee, jr1)
                jr2 = (jr1 Mod 3) + 1
                jr3 = (jr2 Mod 3) + 1
                iv4 = ien(ir, jr3)
                Call SWAP(px(iv2), py(iv2), px(iv1), py(iv1), px(iv4), py(iv4), xx, yy, iswap)
                If iswap = 1 Then
                    ia = jee(il, jl2)
                    ib = jee(ir, jr2)
                    ien(il, jl2) = iv4
                    jee(il, jl1) = ib
                    jee(il, jl2) = ir
                    ien(ir, jr2) = iv3
                    jee(ir, jr1) = ia
                    jee(ir, jr2) = il

                    If ia <> 0 Then
                        Call EDGE(ia, il, jee, iedge)
                        jee(ia, iedge) = ir
                    End If
                    If ib <> 0 Then
                        Call EDGE(ib, ir, jee, iedge)
                        jee(ib, iedge) = il
                    End If
                    itop = itop + 1
                    jstack(itop) = IPUSH(il, maxstk, itop, jstack)
                    GoTo lawson20
                End If
lawson30:
            Next j
            GoTo lawson20
        End If

    End Sub


    Private Sub TRFINE(ByRef NEX As Integer, ByRef xmin As Double, ByRef ymin As Double, ByRef dmax As Double, _
                       ByRef NODE As Integer, ByRef px() As Double, ByRef py() As Double, ByRef jnb() As Integer, _
                       ByRef nei(,) As Integer, ByRef NELM As Integer, ByRef mtj(,) As Integer, ByRef jac(,) As Integer, _
                       ByRef idm() As Integer, ByRef iadres() As Integer, ByRef istack() As Integer, ByRef delx() As Double)

        Dim alpha As Double = 3.0
        Dim i As Integer
        Dim j As Integer
        Dim k As Integer
        Dim jz As Integer
        Dim js As Integer
        Dim ndiv As Integer

        Dim range As Double
        Dim px1 As Double
        Dim py1 As Double
        Dim px2 As Double
        Dim py2 As Double
        Dim px3 As Double
        Dim py3 As Double
        Dim px4 As Double
        Dim py4 As Double
        Dim dltd As Double
        Dim delta As Double
        Dim xg As Double
        Dim yg As Double
        Dim xs As Double
        Dim ys As Double
        Dim xp As Double
        Dim yp As Double

        range = dmax
        px1 = 0.0
        py1 = 0.0
        px2 = 1.0
        py2 = 0.0
        px3 = 1.0
        py3 = 1.0
        px4 = 0.0
        py4 = 1.0
        For i = 1 To NEX
            jz = i
            js = NODE + 1
            If dmax < delx(i) Then GoTo trfine60
            ndiv = CInt(range / delx(i))
            delta = 1.0 / CDbl(ndiv)
            dltd = delta / alpha
            dltd = dltd * dltd
            For k = 1 To ndiv + 1
                xg = px2 + (px3 - px2) / CDbl(ndiv) * CDbl(k - 1)
                xs = px1 + (px4 - px1) / CDbl(ndiv) * CDbl(k - 1)
                yg = py2 + (py3 - py2) / CDbl(ndiv) * CDbl(k - 1)
                ys = py1 + (py4 - py1) / CDbl(ndiv) * CDbl(k - 1)
                For j = 1 To ndiv + 1
                    xp = xs + (xg - xs) / CDbl(ndiv) * CDbl(j - 1)
                    yp = ys + (yg - ys) / CDbl(ndiv) * CDbl(j - 1)
                    If TRPLACE(jz, xp, yp, NELM, mtj, jac, idm, px, py, dltd) Then
                        NODE = NODE + 1
                        px(NODE) = xp
                        py(NODE) = yp
                    End If
                Next j
            Next k

            Call DELAUN(jz, js, NODE, NODE, px, py, jnb, nei, NELM, mtj, jac, idm, iadres, istack)
trfine60:
        Next i

    End Sub


    Private Function TRPLACE(ByRef iz As Integer, ByRef xp As Double, ByRef yp As Double, ByRef NELM As Integer, _
                             ByRef mtj(,) As Integer, ByRef jac(,) As Integer, ByRef idm() As Integer, _
                             ByRef px() As Double, ByRef py() As Double, ByRef dltd As Double) As Boolean

        Dim j As Integer
        Dim loc As Integer
        Dim ia As Integer
        Dim ib As Integer
        Dim jr As Integer
        Dim xa As Double
        Dim ya As Double
        Dim xb As Double
        Dim yb As Double
        Dim dpp As Double
        Dim eps As Double
        Dim eps0 As Double
        Dim epl As Double

        Call LOCATE(xp, yp, px, py, mtj, jac, NELM, loc)

        If iz = idm(loc) Then
            For j = 1 To 3
                ia = mtj(loc, j)
                xa = px(ia)
                ya = py(ia)
                dpp = (xa - xp) * (xa - xp) + (ya - yp) * (ya - yp)
                If dpp < dltd Then
                    TRPLACE = False
                    GoTo trplaceend
                End If
                jr = jac(loc, j)
                If iz <> idm(jr) Then
                    ib = mtj(loc, (j Mod 3) + 1)
                    xb = px(ib)
                    yb = py(ib)
                    eps = (yb - ya) * xp - (xb - xa) * yp - xa * (yb - ya) + ya * (xb - xa)
                    eps = eps * eps
                    eps0 = (yb - ya) * (yb - ya) + (xb - xa) * (xb - xa)
                    epl = eps / eps0
                    If epl < dltd Then
                        TRPLACE = False
                        GoTo trplaceend
                    End If
                End If
            Next j
            TRPLACE = True
        Else
            TRPLACE = False
        End If

trplaceend:
    End Function


    Private Sub REMOVE(ByRef NEX As Integer, ByRef index() As Integer, ByRef NELM As Integer, ByRef mtj(,) As Integer, ByRef jac(,) As Integer, ByRef idm() As Integer)
        Dim i As Integer
        Dim j As Integer
        Dim k As Integer
        Dim l As Integer
        Dim iz As Integer
        Dim inelm As Integer
        Dim ielm As Integer
        Dim jelm As Integer
        Dim kelm As Integer
        Dim lelm As Integer
        Dim mkp(3) As Integer
        Dim jkp(3) As Integer
        Dim ikp As Integer
        Dim kedg As Integer
        Dim ledg As Integer

        ielm = 0
        index(1) = 1
        For i = 1 To NEX
            iz = i
            inelm = 0
            For j = index(i) To NELM
                If idm(j) = iz Then
                    ielm = ielm + 1
                    jelm = j
                    inelm = inelm + 1
                    If ielm <> jelm Then
                        For k = 1 To 3
                            mkp(k) = mtj(ielm, k)
                            jkp(k) = jac(ielm, k)
                        Next k
                        ikp = idm(ielm)
                        For k = 1 To 3
                            kelm = jac(jelm, k)
                            If kelm <> 0 Then
                                Call EDGE(kelm, jelm, jac, kedg)
                                jac(kelm, kedg) = ielm + KTE + 1
                            End If
                        Next k
                        For l = 1 To 3
                            lelm = jkp(l)
                            If lelm <> 0 Then
                                Call EDGE(lelm, ielm, jac, ledg)
                                jac(lelm, ledg) = jelm + KTE + 1
                            End If
                        Next l
                        For k = 1 To 3
                            jac(ielm, k) = jac(ielm, k) Mod (KTE + 1)
                            jac(jelm, k) = jac(jelm, k) Mod (KTE + 1)
                            kelm = jac(ielm, k)
                            lelm = jac(jelm, k)
                            For l = 1 To 3
                                jac(kelm, l) = jac(kelm, l) Mod (KTE + 1)
                                jac(lelm, l) = jac(lelm, l) Mod (KTE + 1)
                            Next l
                        Next k
                        For k = 1 To 3
                            jkp(k) = jac(ielm, k)
                        Next k
                        For k = 1 To 3
                            mtj(ielm, k) = mtj(jelm, k)
                            jac(ielm, k) = jac(jelm, k)
                            mtj(jelm, k) = mkp(k)
                            jac(jelm, k) = jkp(k)
                        Next k
                        idm(ielm) = idm(jelm)
                        idm(jelm) = ikp
                    End If
                End If
            Next j
            index(i + 1) = index(i) + inelm
        Next i

        For i = 1 To ielm
            For j = 1 To 3
                If ielm < jac(i, j) Then
                    jac(i, j) = 0
                End If
            Next j
        Next i
        For i = ielm + 1 To NELM
            For j = 1 To 3
                mtj(i, j) = 0
                jac(i, j) = 0
            Next j
        Next i

        NELM = ielm
    End Sub


    Private Sub CHECK(ByRef NELM As Integer, ByRef mtj(,) As Integer, ByRef jac(,) As Integer)
        Dim i As Integer
        Dim j As Integer
        Dim k As Integer
        Dim ielm As Integer
        Dim jelm As Integer
        Dim kelm As Integer
        Dim ia As Integer
        Dim ib As Integer
        Dim ja As Integer
        Dim jb As Integer

        For i = 1 To NELM
            For j = 1 To 3
                ielm = i
                ia = mtj(i, j)
                ib = mtj(i, (j Mod 3) + 1)
                jelm = jac(i, j)
                If jelm = 0 Then GoTo check10
                For k = 1 To 3
                    kelm = jac(jelm, k)
                    If kelm = ielm Then
                        ja = mtj(jelm, k)
                        jb = mtj(jelm, (k Mod 3) + 1)
                        If (ia = jb) And (ib = ja) Then GoTo check10
                        'ERRORメッセージ＆実行中止
                        MessageBox.Show("ERROR IN SUB CHECK", "Error終了")
                        Me.Close()
                    End If
                Next k
                'ERRORメッセージ＆実行中止
                MessageBox.Show("ERROR IN SUB CHECK", "Error終了")
                Me.Close()
check10:
            Next j
        Next i
    End Sub


    Private Sub INCR(ByRef n As Integer, ByRef l As Integer, ByRef jnb() As Integer, ByRef nei(,) As Integer)

        If n <= KTJ Then
            jnb(n) = jnb(n) + 1
            If KCM < jnb(n) Then
                'ERRORメッセージ＆実行中止
                MessageBox.Show("ERROR IN SUB INCR", "Error終了")
                Me.Close()
            End If
            nei(n, jnb(n)) = l
        End If
    End Sub


    Private Sub DECR(ByRef n As Integer, ByRef l As Integer, ByRef jnb() As Integer, ByRef nei(,) As Integer)
        Dim i As Integer
        Dim inb As Integer

        If n <= KTJ Then
            inb = NEIBOR(n, l, jnb, nei)
            For i = inb To jnb(n) - 1
                nei(n, i) = nei(n, i + 1)
            Next i
            nei(n, jnb(n)) = 0
            jnb(n) = jnb(n) - 1
            If jnb(n) = 0 Then
                'ERRORメッセージ＆実行中止
                MessageBox.Show("ERROR IN SUB DECR", "Error終了")
                Me.Close()
            End If
        End If
    End Sub


    Private Sub LOCATE(ByRef xp As Double, ByRef yp As Double, ByRef px() As Double, ByRef py() As Double, _
                       ByRef mtj(,) As Integer, ByRef jac(,) As Integer, ByRef NELM As Integer, ByRef itri As Integer)
        Dim i As Integer
        Dim ia As Integer
        Dim ib As Integer

        itri = NELM
locate10:
        For i = 1 To 3
            ia = mtj(itri, i)
            ib = mtj(itri, (i Mod 3) + 1)
            If (px(ia) - xp) * (py(ib) - yp) < (py(ia) - yp) * (px(ib) - xp) Then
                itri = jac(itri, i)
                GoTo locate10
            End If
        Next i
    End Sub


    Private Sub SWAP(ByRef x1 As Double, ByRef y1 As Double, ByRef x2 As Double, ByRef y2 As Double, ByRef x3 As Double, ByRef y3 As Double, _
                     ByRef xp As Double, ByRef yp As Double, ByRef iswap As Integer)
        Dim x13 As Double
        Dim y13 As Double
        Dim x23 As Double
        Dim y23 As Double
        Dim x1p As Double
        Dim y1p As Double
        Dim x2p As Double
        Dim y2p As Double
        Dim cosa As Double
        Dim cosb As Double
        Dim sina As Double
        Dim sinb As Double

        x13 = x1 - x3
        y13 = y1 - y3
        x23 = x2 - x3
        y23 = y2 - y3
        x1p = x1 - xp
        y1p = y1 - yp
        x2p = x2 - xp
        y2p = y2 - yp
        cosa = x13 * x23 + y13 * y23
        cosb = x2p * x1p + y1p * y2p
        If (0.0 <= cosa) And (0.0 <= cosb) Then
            iswap = 0
        ElseIf (cosa < 0.0) And (cosb < 0.0) Then
            iswap = 1
        Else
            sina = x13 * y23 - x23 * y13
            sinb = x2p * y1p - x1p * y2p
            If (sina * cosb + sinb * cosa) < 0 Then
                iswap = 1
            Else
                iswap = 0
            End If
        End If
    End Sub


    Private Sub EDGE(ByRef l As Integer, ByRef k As Integer, ByRef jac(,) As Integer, ByRef iedge As Integer)
        Dim i As Integer

        iedge = -1
        For i = 1 To 3
            If jac(l, i) = k Then
                iedge = i
                Exit For
            End If
        Next i
        'ERRORメッセージ＆実行中止
        If iedge = -1 Then
            MessageBox.Show("ERROR IN SUB EDGE", "Error終了")
            Me.Close()
        End If
    End Sub


    Private Function IPUSH(ByRef item As Integer, ByRef maxstk As Integer, ByRef itop As Integer, ByRef istack() As Integer) As Integer

        If maxstk < itop Then
            'ERRORメッセージ＆実行中止
            MessageBox.Show("ERROR IN FUNCTION IPUSH", "Error終了")
            Me.Close()
        Else
            IPUSH = item
        End If
    End Function


    Private Function NEIBOR(ByRef n As Integer, ByRef l As Integer, ByRef jnb() As Integer, ByRef nei(,) As Integer) As Integer
        Dim i As Integer

        NEIBOR = -1
        For i = 1 To jnb(n)
            If nei(n, i) = l Then
                NEIBOR = i
                Exit For
            End If
        Next i
        'ERRORメッセージ＆実行中止
        If NEIBOR = -1 Then
            MessageBox.Show("ERROR IN FUNCTION NEIBOR", "Error終了")
            Me.Close()
        End If
    End Function


    Private Function IVERT(ByRef l As Integer, ByRef k As Integer, ByRef mtj(,) As Integer) As Integer
        Dim i As Integer

        IVERT = -1
        For i = 1 To 3
            If mtj(l, i) = k Then
                IVERT = i
                Exit For
            End If
        Next i
        'ERRORメッセージ＆実行中止
        If IVERT = -1 Then
            MessageBox.Show("ERROR IN FUNCTION IVERT", "Error終了")
            Me.Close()
        End If
    End Function


    Private Sub TRDATA(ByRef NEX As Integer, ByRef NODE As Integer, ByRef NELM As Integer, ByRef mtj(,) As Integer, _
                      ByRef jac(,) As Integer, ByRef idm() As Integer, _
                      ByRef px() As Double, ByRef py() As Double, ByRef ifix() As Integer)
        Dim sw As System.IO.StreamWriter
        Dim dat As String
        Dim fname2 As String
        Dim i As Integer

        fname2 = My.Forms.Form1.Label1.Text
        fname2 = System.IO.Path.GetDirectoryName(fname2)
        fname2 = fname2 & "\workf2.csv"
        sw = New System.IO.StreamWriter(fname2, False, System.Text.Encoding.Default)
        dat = NODE.ToString & "," & NELM.ToString : sw.WriteLine(dat)
        For i = 1 To NELM
            dat = mtj(i, 1).ToString & "," & mtj(i, 2).ToString & "," & mtj(i, 3).ToString
            dat = dat & "," & jac(i, 1).ToString & "," & jac(i, 2).ToString & "," & jac(i, 3).ToString
            dat = dat & "," & idm(i).ToString
            sw.WriteLine(dat)
        Next i
        For i = 1 To NODE
            dat = px(i).ToString("0.000") & "," & py(i).ToString("0.000") & "," & ifix(i).ToString
            sw.WriteLine(dat)
        Next i
        sw.Close()

    End Sub

End Class