Archive for November 15th, 2009

Chain model – simple simulation

Sunday, November 15th, 2009

Wrote a simple script that simulates chain models in Microstation. A series of chain segments of a specific length is defined. Some knots are fix others can move.

 

‘Simulation of a chain model
‘Sebastian Gmelin
’15.11.2009

Option Explicit

Const gravity As Double = 1 / 5
Const springforce As Double = 1 / 8
Const attrForce As Double = 3

Public Type sphere
    number As Long
    radius As Double
    center As Point3d
    newCenter As Point3d
    sElements(3) As Element
    numNeighbours As Long
    Neighbours(10) As Long
    connections(10) As Element
    moveable As Boolean
End Type

Dim spheres() As sphere

Function createSphere(pos As Point3d, rad As Double, number As Long, moveable As Boolean) As sphere
   
    Dim elem As Element
    Dim mat As Matrix3d
   
    createSphere.center = pos
    createSphere.newCenter = pos
    createSphere.radius = rad
    createSphere.moveable = moveable
    createSphere.number = number
    createSphere.numNeighbours = 0
   
    Set elem = CreateLineElement2(Nothing, pos, pos)
    elem.Color = 3
    If moveable Then elem.Color = 1
    elem.LineWeight = 8
    elem.LineStyle = ActiveDesignFile.LineStyles(1)
    ActiveModelReference.AddElement elem
    elem.Redraw
    Set createSphere.sElements(0) = elem
   
    mat = Matrix3dFromAxisAndRotationAngle(2, 0)
    Set elem = CreateEllipseElement2(Nothing, pos, rad, rad, mat, msdFillModeNotFilled)
    elem.Color = 6
    elem.LineWeight = 0
    elem.LineStyle = ActiveDesignFile.LineStyles(2)
    ActiveModelReference.AddElement elem
    elem.Redraw
    Set createSphere.sElements(1) = elem
   
    mat = Matrix3dFromAxisAndRotationAngle(1, Pi / 2)
    Set elem = CreateEllipseElement2(Nothing, pos, rad, rad, mat, msdFillModeNotFilled)
    elem.Color = 6
    elem.LineWeight = 0
    elem.LineStyle = ActiveDesignFile.LineStyles(2)
    ActiveModelReference.AddElement elem
    elem.Redraw
    Set createSphere.sElements(2) = elem
   
    mat = Matrix3dFromAxisAndRotationAngle(0, Pi / 2)
    Set elem = CreateEllipseElement2(Nothing, pos, rad, rad, mat, msdFillModeNotFilled)
    elem.Color = 6
    elem.LineWeight = 0
    elem.LineStyle = ActiveDesignFile.LineStyles(2)
    ActiveModelReference.AddElement elem
    elem.Redraw
    Set createSphere.sElements(3) = elem
   
End Function

Sub moveSphere(s As sphere, vec As Point3d)

    Dim i As Long
   
    s.center = Point3dAdd(s.center, vec)
   
    For i = 0 To UBound(s.sElements)
        Call s.sElements(i).Move(vec)
        s.sElements(i).Rewrite
        s.sElements(i).Redraw
        ‘ActiveDesignFile.Views(1).Redraw
    Next i
   
End Sub

Sub moveConnections()

    Dim s As Long
    Dim c As Long
   
    For s = 0 To UBound(spheres)
        For c = 0 To spheres(s).numNeighbours – 1
            Call spheres(s).connections(c).AsLineElement.ModifyVertex(0, spheres(s).center)
            Call spheres(s).connections(c).AsLineElement.ModifyVertex(1, spheres(spheres(s).Neighbours(c)).center)
            spheres(s).connections(c).Rewrite
            spheres(s).connections(c).Redraw
        Next c
    Next s

End Sub

Sub makeSphereFix(s As Long)

    spheres(s).moveable = False
    spheres(s).sElements(0).Color = 3
    spheres(s).sElements(0).Rewrite
    spheres(s).sElements(0).Redraw

End Sub

Sub makeSphereMoveAble(s As Long)

    spheres(s).moveable = True
    spheres(s).sElements(0).Color = 1
    spheres(s).sElements(0).Rewrite
    spheres(s).sElements(0).Redraw

End Sub

Sub drawConnections(s As Long)

    Dim elem As Element
    Dim points(1) As Point3d
    Dim i As Long
   
    For i = 0 To spheres(s).numNeighbours – 1
        points(0) = spheres(s).center
        points(1) = spheres(spheres(s).Neighbours(i)).center
        Set elem = CreateLineElement1(Nothing, points)
        elem.Color = 2
        elem.LineWeight = 2
        elem.LineStyle = ActiveDesignFile.LineStyles(1)
        ActiveModelReference.AddElement elem
        elem.Redraw
        Set spheres(s).connections(i) = elem
    Next i

End Sub

Sub calculateSphereMove(s As Long)

    Dim i As Long
    Dim k As Long
    Dim direction As Point3d
    Dim dist As Double
    Dim length As Double
    Dim grav As Boolean
    Dim vec As Point3d
    Dim tension As Double
   
    If spheres(s).moveable Then
    tension = 0
        For i = 0 To spheres(s).numNeighbours – 1
            direction = Point3dSubtract(spheres(spheres(s).Neighbours(i)).center, spheres(s).newCenter)
            dist = Abs(Point3dMagnitude(direction))
            ‘Debug.Print (“distance ” + CStr(dist))
            length = spheres(s).radius + spheres(spheres(s).Neighbours(i)).radius
            ‘Debug.Print (“length ” + CStr(length))
            tension = tension + dist – length
        Next i
        tension = tension / (spheres(s).numNeighbours)
        grav = True
        If tension > 0 Then
            spheres(s).newCenter = Point3dAdd(spheres(s).newCenter, Point3dFromXYZ(0, 0, -0.1 * gravity))
        Else
            spheres(s).newCenter = Point3dAdd(spheres(s).newCenter, Point3dFromXYZ(0, 0, -1 * gravity))
        End If
        vec = Point3dFromXYZ(0, 0, 0)
        For i = 0 To spheres(s).numNeighbours – 1
            direction = Point3dSubtract(spheres(spheres(s).Neighbours(i)).center, spheres(s).newCenter)
            dist = Abs(Point3dMagnitude(direction))
            ‘Debug.Print (“distance ” + CStr(dist))
            length = spheres(s).radius + spheres(spheres(s).Neighbours(i)).radius
            ‘Debug.Print (“length ” + CStr(length))
            If (dist – length) > 0 Then
                grav = False
                direction = Point3dScale(direction, attrForce * springforce * (dist – length) / (dist))
            Else
                direction = Point3dScale(direction, springforce * (dist – length) / (dist))
            End If
            vec = Point3dAdd(vec, direction)
        Next i
        spheres(s).newCenter = Point3dAdd(spheres(s).newCenter, vec)
    End If

End Sub

Sub moveSpheres()

    Dim i As Long
    Dim s As Long
    Dim vec As Point3d
   
    For s = 0 To UBound(spheres)
        vec = Point3dSubtract(spheres(s).newCenter, spheres(s).center)
        spheres(s).center = spheres(s).newCenter
   
        For i = 0 To UBound(spheres(s).sElements)
            Call spheres(s).sElements(i).Move(vec)
            spheres(s).sElements(i).Rewrite
            spheres(s).sElements(i).Redraw
            ‘ActiveDesignFile.Views(1).Redraw
        Next i
    Next s

End Sub

 

Sub Main()

    Dim i As Long
    Dim k As Long
    Dim s As Long
   
    Randomize

    ReDim spheres(99)
    For i = 0 To 9
        For k = 0 To 9
            spheres(i * 10 + k) = createSphere(Point3dFromXYZ(i, k, 0), 0.65, i * 10 + k, True)
        Next k
    Next i
   
    For i = 0 To 9
        For k = 0 To 9
            ‘spheres(i * 10 + k) = createSphere(Point3dFromXYZ(i, k, 0), rnd * 0.2 + 0.5, i * 10 + k, True)
            If i < 9 Then
                spheres(i * 10 + k).Neighbours(spheres(i * 10 + k).numNeighbours) = i * 10 + k + 10
                spheres(i * 10 + k).numNeighbours = spheres(i * 10 + k).numNeighbours + 1
            End If
            If k < 9 Then
                spheres(i * 10 + k).Neighbours(spheres(i * 10 + k).numNeighbours) = i * 10 + k + 1
                spheres(i * 10 + k).numNeighbours = spheres(i * 10 + k).numNeighbours + 1
            End If
            If i > 0 Then
                spheres(i * 10 + k).Neighbours(spheres(i * 10 + k).numNeighbours) = i * 10 + k – 10
                spheres(i * 10 + k).numNeighbours = spheres(i * 10 + k).numNeighbours + 1
            End If
            If k > 0 Then
                spheres(i * 10 + k).Neighbours(spheres(i * 10 + k).numNeighbours) = i * 10 + k – 1
                spheres(i * 10 + k).numNeighbours = spheres(i * 10 + k).numNeighbours + 1
            End If
            Call drawConnections(i * 10 + k)
        Next k
    Next i
   
    ‘Call makeSphereFix(0)
    Call makeSphereFix(9)
    ‘Call makeSphereFix(55)
    Call makeSphereFix(90)
    Call makeSphereFix(99)
   
    For i = 0 To 75
‘        Call moveSphere(spheres(Round(Rnd * 99)), Point3dFromXYZ(0, 0, Rnd * 0.5 – 0.25))
‘        Call moveConnections
        For s = 0 To UBound(spheres)
            Call calculateSphereMove(s)
        Next s
        Call moveSpheres
        Call moveConnections
        ActiveDesignFile.Views(1).Redraw
    Next i
   
End Sub

Share