ConvexContactTest.py

You can view and download this file on Github: ConvexContactTest.py

  1#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  2# This is an EXUDYN example
  3#
  4# Details:  Test with ObjectContactConvexRoll, which models a roll of a mechanum wheel or any other roll
  5#           which is described by a polynomial profile
  6#
  7# Author:   Peter Manzl
  8# Date:     2021-12-21
  9#
 10# Copyright:This file is part of Exudyn. Exudyn is free software. You can redistribute it and/or modify it under the terms of the Exudyn license. See 'LICENSE.txt' for more details.
 11#
 12#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 13
 14import exudyn as exu
 15from exudyn.utilities import * #includes itemInterface and rigidBodyUtilities
 16import exudyn.graphics as graphics #only import if it does not conflict
 17
 18useGraphics = True #without test
 19#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 20#you can erase the following lines and all exudynTestGlobals related operations if this is not intended to be used as TestModel:
 21try: #only if called from test suite
 22    from modelUnitTests import exudynTestGlobals #for globally storing test results
 23    useGraphics = exudynTestGlobals.useGraphics
 24except:
 25    class ExudynTestGlobals:
 26        pass
 27    exudynTestGlobals = ExudynTestGlobals()
 28#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 29
 30SC = exu.SystemContainer()
 31mbs = SC.AddSystem()
 32
 33
 34# create Ground and graphics:
 35scb, g1, g2 = 0.5, 0.92, 0.72
 36graphGround = graphics.CheckerBoard(point= [0,0,0], normal= [0,0,1], size= scb, color= [g1, g1, g1, 1.0],
 37                                        alternatingColor= [g2, g2, g2, 1.0], nTiles= 12)
 38oGround = mbs.CreateGround(referencePosition = [0,0,0],
 39                           graphicsDataList=[graphGround])
 40mGround = mbs.AddMarker(MarkerBodyRigid(bodyNumber=oGround,localPosition = [0,0,0], name='Ground'))
 41nGround0=mbs.AddNode(NodePointGround(referenceCoordinates = [0,0,0]))
 42
 43
 44poly = [-3.6, 0,  1.65e-02] # coefficients of the polynomial creating the rolling body
 45length = 0.1                # length of the roller
 46contour =  [[-length/2, 0]]
 47x = np.linspace(start = - length/2, stop =length/2, num=51)
 48for i in range(np.size(x)):
 49    contour+= [[x[i], np.polyval(poly, x[i])]]
 50contour += [ [length/2, 0]] # to create a closed contour
 51graphRoll = [graphics.SolidOfRevolution([0,0,0], [1,0,0], contour, graphics.color.lightred[0:3]+[1],
 52                                          alternatingColor=graphics.color.blue, nTiles = 32)]
 53
 54InertiaRoll = InertiaCylinder(density=7800, length=length, outerRadius=3e-3, axis=0)
 55bRoll = mbs.CreateRigidBody(inertia = InertiaRoll,
 56                            referencePosition = [0,0,poly[-1]*1.2],
 57                            referenceRotationMatrix =RotationMatrixY(np.pi/16),
 58                            initialAngularVelocity = RotationMatrixY(np.pi/16) @ np.array([-1000,0,0]),  # in Global coordinates
 59                            initialVelocity= [0,0,0],
 60                            gravity = [0,0,-9.81],
 61                            graphicsDataList = graphRoll)
 62
 63nRoll = mbs.GetObject(bRoll)['nodeNumber']
 64
 65nData = mbs.AddNode(NodeGenericData(initialCoordinates=[0,0,0], numberOfDataCoordinates=3))
 66mRoll = mbs.AddMarker(MarkerNodeRigid(nodeNumber=nRoll))
 67CConvexRoll = mbs.AddObject(ObjectContactConvexRoll(markerNumbers=[mGround, mRoll],
 68                        nodeNumber=nData, contactStiffness=1e3, contactDamping=1, dynamicFriction = 0.9,
 69                       staticFrictionOffset = 0, viscousFriction=0, exponentialDecayStatic=1e-3,
 70                       frictionProportionalZone=1e-4, rollLength=length, coefficientsHull=poly,
 71                       visualization={'show': True, 'color': graphics.color.lightgreen}))
 72
 73sBody = mbs.AddSensor(SensorBody(bodyNumber=bRoll, #fileName='PosRoller.txt',
 74                                 storeInternal=True,
 75                                 outputVariableType=exu.OutputVariableType.Position, visualization={'show': False}))
 76
 77mbs.Assemble()
 78
 79h = 5e-4   #test
 80tEnd = 0.1 #test
 81#tEnd = 0.1*20 #for simulation
 82sims=exu.SimulationSettings()
 83sims.timeIntegration.generalizedAlpha.spectralRadius=0.7
 84sims.timeIntegration.endTime = tEnd
 85sims.timeIntegration.numberOfSteps = int(tEnd/h) #original: 1e-3, fails now in Newton
 86sims.timeIntegration.verboseMode = 0
 87sims.timeIntegration.stepInformation = 3 #do not show step reduction
 88sims.solutionSettings.coordinatesSolutionFileName = 'solution/coordinatesSolution.txt'
 89# sims.timeIntegration.newton.absoluteTolerance = 1e-8
 90# sims.timeIntegration.newton.relativeTolerance = 1e-6
 91
 92if useGraphics:
 93    sims.timeIntegration.verboseMode = 1
 94    sims.timeIntegration.stepInformation = 3+128+256
 95    exu.StartRenderer()
 96    mbs.WaitForUserToContinue()
 97mbs.SolveDynamic(sims)
 98if useGraphics:
 99    SC.WaitForRenderEngineStopFlag()
100    exu.StopRenderer() #safely close rendering window!
101
102sol = mbs.systemData.GetODE2Coordinates()
103exudynTestGlobals.testResult = np.sum(sol[:2])
104exu.Print('result of ConvexContactTest=',exudynTestGlobals.testResult)
105# %%
106if useGraphics:
107    #pos = np.loadtxt('PosRoller.txt', delimiter=',', comments='#')
108    pos = mbs.GetSensorStoredData(sBody)
109    exu.Print('End Pos: {}'.format(pos[-1,:]))
110
111
112    mbs.PlotSensor(sBody,[0,1,2])
113
114
115if useGraphics and False:
116    SC.visualizationSettings.general.autoFitScene = False
117    SC.visualizationSettings.general.graphicsUpdateInterval=0.02
118
119    sol = LoadSolutionFile('solution/coordinatesSolution.txt', safeMode=True)#, maxRows=100)
120    print('start SolutionViewer')
121    mbs.SolutionViewer(sol)