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