connectorRigidBodySpringDamperTest.py

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

 1#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 2# This is an EXUDYN example
 3#
 4# Details:  Test for RigidBodySpringDamper with springForceTorqueUserFunction;
 5#           the RigidBodySpringDamper can be used to model complicance effects in joints where 
 6#           one axis undergoes large rotations, and the other rotations are small
 7#
 8# Author:   Johannes Gerstmayr
 9# Date:     2021-01-07
10#
11# 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.
12#
13#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
14
15import sys
16sys.path.append('../TestModels')            #for modelUnitTest as this example may be used also as a unit test
17
18import exudyn as exu
19from exudyn.utilities import * #includes itemInterface and rigidBodyUtilities
20import exudyn.graphics as graphics
21import numpy as np
22
23testIsActive = exu.sys.get('testIsActive', False)
24
25SC = exu.SystemContainer()
26mbs = SC.AddSystem()
27
28oGround=mbs.AddObject(ObjectGround(referencePosition= [0,0,0]))
29nGround = mbs.AddNode(NodePointGround(referenceCoordinates=[0,0,0]))
30
31#example with rigid body at [0,0,0], 1kg under initial velocity
32graphicsBody = graphics.Brick(centerPoint=[0,0,0],size=[0.09,0.09,0.2], color=graphics.color.lightred)
33nBody = mbs.AddNode(RigidRxyz(initialVelocities=[0,10,0, 2*pi*4,0,0]))
34oBody = mbs.AddObject(RigidBody(mass=1, inertia=[1,1,1,0,0,0], nodeNumber=nBody, 
35                                visualization=VRigidBody(graphicsData=[graphicsBody])))
36
37mBody = mbs.AddMarker(MarkerNodeRigid(nodeNumber=nBody))
38mGround = mbs.AddMarker(MarkerBodyRigid(bodyNumber=oGround, 
39                                        localPosition = [0,0,0]))
40
41def UFforce(mbs, t, itemIndex, displacement, rotation, velocity, angularVelocity, stiffness, damping, rotJ0, rotJ1, offset): 
42    k = stiffness
43    u = displacement
44    v = velocity
45    w = angularVelocity
46    rot = rotation
47    return [u[0]*k[0][0],u[1]*k[1][1]+v[1]*0.01*k[1][1],u[2]*k[2][2],
48            rot[0]*k[0][0]+w[0]*0.001*k[0][0],rot[1]*k[0][0],rot[2]*k[0][0]]
49
50#markerNumbers and parameters taken from mini example
51k=5000
52mbs.AddObject(RigidBodySpringDamper(markerNumbers = [mGround, mBody], 
53                                    stiffness = np.diag([k,k,k, 0,0,0]), 
54                                    damping = np.diag([0,k*0.01,0, 0,0,0]), 
55                                    offset = [0,1,0, 0,0,0],
56                                    springForceTorqueUserFunction = UFforce))
57
58mbs.Assemble()
59
60tEnd = 0.1
61h=1e-3
62if not testIsActive:
63    tEnd = 1 #parameters sucht that we can see some motion
64    h=1e-5
65
66
67simulationSettings = exu.SimulationSettings()
68simulationSettings.timeIntegration.newton.useModifiedNewton = False #Just for the test; modified Newton is usually faster
69simulationSettings.timeIntegration.numberOfSteps = int(tEnd/h)
70simulationSettings.timeIntegration.endTime = tEnd
71simulationSettings.solution.file.write = False
72
73simulationSettings.timeIntegration.generalizedAlpha.spectralRadius = 1 #no numerical damping
74
75simulationSettings.show.statistics = True
76simulationSettings.timeIntegration.verboseMode = 1
77
78if not testIsActive:
79    SC.renderer.Start()              #start graphics visualization
80    SC.renderer.DoIdleTasks()    #wait for pressing SPACE bar to continue
81
82#start solver:
83mbs.SolveDynamic(simulationSettings)
84
85if not testIsActive:
86    SC.renderer.DoIdleTasks()#wait for pressing 'Q' to quit
87    SC.renderer.Stop()               #safely close rendering window!
88
89p0=mbs.GetObjectOutputBody(oBody, localPosition=[0.1,0.1,0.1], 
90                           variableType = exu.OutputVariableType.Position)
91result = p0[0]+p0[1]
92exu.Print('solution of connectorRigidBodySpringDamperTest=',result) #use x-coordinate
93
94exu.sys['testResult'] = result