springMassFriction.py

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

  1#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  2# This is an EXUDYN example
  3#
  4# Details:  A simple mass-spring system with friction
  5#
  6# Author:   Johannes Gerstmayr
  7# Date:     2024-06-07
  8#
  9# 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.
 10#
 11#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 12
 13import exudyn as exu
 14from exudyn.itemInterface import *
 15from exudyn.utilities import * #includes itemInterface and rigidBodyUtilities
 16import exudyn.graphics as graphics
 17from exudyn.physics import StribeckFunction
 18import matplotlib.pyplot as plt
 19
 20import numpy as np
 21from math import sin, cos, pi, sqrt
 22
 23import time #for sleep()
 24SC = exu.SystemContainer()
 25mbs = SC.AddSystem()
 26
 27#%%+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 28#
 29L = 1
 30spring = 1600;               #stiffness [800]
 31mass = 1;                   #mass
 32force = 200;                  #force amplitude
 33
 34stepSize = 1e-4
 35tEnd = 1
 36
 37#%%+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 38#create model for linear and nonlinear oscillator:
 39# L=0.5
 40# load0 = 80
 41
 42# omega0=np.sqrt(spring/mass)
 43# f0 = 0.*omega0/(2*np.pi)
 44# f1 = 1.*omega0/(2*np.pi)
 45
 46# #user function for spring force
 47# def springForce(mbs, t, itemIndex, u, v, k, d, offset, mu, muPropZone):
 48#     k=mbs.variables['stiffness']
 49#     d=mbs.variables['damping']
 50#     if mbs.variables['mode'] == 0:
 51#         return k*u + v*d
 52#     else:
 53#         #return 0.1*k*u+k*u**3+v*d
 54#         return k*u+1000*k*u**3+v*d #breaks down at 13.40Hz
 55
 56# mode = 0 #0...forward, 1...backward
 57
 58#user function for load
 59def UserLoad3D(mbs, t, load):
 60    f = force
 61    if t >= 0.5: f = 0.2*force
 62    return [f, 0, 0]
 63
 64def UserSpringForce(mbs, t, itemIndex, u, v, k, d, f0):
 65    muFn = 20
 66    ff = muFn*StribeckFunction(v, 1, 0.1)
 67    return u*k+ff
 68
 69gSphere=graphics.Sphere(radius=0.1*L, nTiles=32, color=graphics.color.orange)
 70gGroundList = [graphics.Brick(centerPoint=[-0.1,0,0],size=[0.2,0.4,0.4],color=graphics.color.grey)]
 71
 72oGround = mbs.CreateGround(graphicsDataList=gGroundList)
 73oMass = mbs.CreateMassPoint(referencePosition=[L,0,0],
 74                            mass=mass,
 75                            graphicsDataList=[gSphere])
 76
 77mMass = mbs.AddMarker(MarkerBodyPosition(bodyNumber = oMass))
 78oSpringDamper = mbs.CreateSpringDamper(itemNumbers=[oGround, oMass],
 79                                       referenceLength = L,
 80                                       stiffness = spring,
 81                                       springForceUserFunction=UserSpringForce,
 82                                       drawSize = 0.1*L,color=graphics.color.dodgerblue)
 83
 84mbs.AddLoad(Force(markerNumber=mMass, loadVector=[0,0,0],
 85                      loadVectorUserFunction=UserLoad3D))
 86
 87sPos = mbs.AddSensor(SensorBody(bodyNumber=oMass, storeInternal=True,
 88                                outputVariableType=exu.OutputVariableType.Position))
 89
 90#%%+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 91mbs.Assemble()
 92
 93
 94SC.visualizationSettings.view0.window.globalFontSize = 12
 95SC.visualizationSettings.openGL.lineWidth = 4
 96SC.visualizationSettings.openGL.multiSampling = 4
 97SC.visualizationSettings.general.graphicsUpdateInterval = 0.005
 98#SC.visualizationSettings.view0.window.renderWindowSize=[1024,900]
 99SC.visualizationSettings.view0.window.renderWindowSize=[1600,1000]
100SC.visualizationSettings.general.showSolverInformation = False
101SC.visualizationSettings.view0.scene.drawCoordinateSystem = False
102
103SC.visualizationSettings.loads.fixedLoadSize=0
104SC.visualizationSettings.loads.loadSizeFactor=0.5
105SC.visualizationSettings.loads.drawSimplified=False
106SC.visualizationSettings.loads.defaultSize=1
107SC.visualizationSettings.loads.defaultRadius=0.01
108
109SC.visualizationSettings.general.autoFitScene = True #otherwise, renderState not accepted for zoom
110
111#++++++++++++++++++++++++++++++++++++++++
112#setup simulation settings and run interactive dialog:
113simulationSettings = exu.SimulationSettings()
114simulationSettings.timeIntegration.generalizedAlpha.spectralRadius = 1
115simulationSettings.solution.file.write = True
116simulationSettings.solution.file.writePeriod = 0.005 #data not used
117simulationSettings.solution.sensors.writePeriod = 0.002 #data not used
118simulationSettings.solution.file.information = 'mass-spring-friction-oscillatior'
119simulationSettings.timeIntegration.verboseMode = 1 #turn off, because of lots of output
120
121simulationSettings.timeIntegration.numberOfSteps = int(tEnd/stepSize)
122simulationSettings.timeIntegration.endTime = tEnd
123simulationSettings.timeIntegration.newton.useModifiedNewton = True
124
125simulationSettings.show.computationTime = True
126# simulationSettings.numberOfThreads = 1
127
128#%%+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
129#set up interactive window
130
131
132mbs.SolveDynamic(simulationSettings=simulationSettings,
133                 solverType=exu.DynamicSolverType.ExplicitMidpoint)
134
135mbs.SolutionViewer()