.. _examples-springdampermasspointsystem: ****************************** SpringDamperMasspointSystem.py ****************************** You can view and download this file on Github: `SpringDamperMasspointSystem.py `_ .. code-block:: python :linenos: #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ # This is a EXUDYN example # # Details: The 3D movement of a point mass system is simulated. # The point masses are connected by spring-damper elements. # The system ist modelled in accordance with the Bachelor thesis of Thomas # Pfurtscheller (SS/WS19) # # Author: Holzinger Stefan # Date: 2019-10-14 # # 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. # #+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ import exudyn as exu from exudyn.itemInterface import * from exudyn.utilities import * #includes itemInterface and rigidBodyUtilities import exudyn.graphics as graphics #only import if it does not conflict import numpy as np print('EXUDYN version='+exu.GetVersionString()) #++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ # Modelling of MBS starts here # define number of nodes in mesh numberOfNodesX = 6 numberOfNodesY = 4 numberOfNodes = numberOfNodesX*numberOfNodesY # define undeformed spring length springDamperElementLength = 1 springDamperDiagnoalElementLength = np.sqrt( springDamperElementLength*springDamperElementLength + springDamperElementLength*springDamperElementLength ) # define mass point mass and spring-damper properties massMassPoint = 1/2 # each node has half mass, since elements have mass elementStiffness = 1000 elementDamping = 5 activateElement = 1 # Flags clambRightSide = True # clamb right side of mesh. If false, only left side of system is clambed useSpringDamperElementsX = True useSpringDamperElementsY = True useSpringDamperElementsDiagonal = True useSpringDamperElementsOffDiagonal = True createSolutionFile = False # this part of the code is used for convergence analysis # if hend < h0, step size will be reduced until this condition is false h0 = 1.0e-5 hend = 1.0e-5 stepSizeReductionFactor = 2 # factor which is used to reduce time step size (increase number of steps) stepSizeList = list() stepSizeList.append(h0) ctr1 = 0 while stepSizeList[ctr1] > hend: stepSizeList.append( stepSizeList[ctr1]/stepSizeReductionFactor ) ctr1 = ctr1 + 1 # simulation settings tEnd = 5.0 # s # simulate system until tend # actual system definition starts here for i in range( len(stepSizeList) ): SC = exu.SystemContainer() mbs = SC.AddSystem() simulationSettings = exu.SimulationSettings() #takes currently set values or default values simulationSettings.solutionSettings.coordinatesSolutionFileName = 'BASpringDamperSystem h=' + str(stepSizeList[i]) + '.txt' print(int( tEnd/stepSizeList[i] ) ) writeStep = 0.01 simulationSettings.timeIntegration.numberOfSteps = int( tEnd/stepSizeList[i] ) simulationSettings.timeIntegration.endTime = tEnd simulationSettings.solutionSettings.writeSolutionToFile = createSolutionFile simulationSettings.solutionSettings.solutionWritePeriod = writeStep simulationSettings.solutionSettings.outputPrecision = 16 simulationSettings.displayComputationTime = True simulationSettings.timeIntegration.verboseMode = 1 #simulationSettings.linearSolverType = exu.LinearSolverType.EigenSparse #simulationSettings.timeIntegration.generalizedAlpha.spectralRadius = 1 #simulationSettings.displayStatistics = True SC.visualizationSettings.bodies.showNumbers = True SC.visualizationSettings.nodes.defaultSize = 0.1 SC.visualizationSettings.markers.defaultSize = 0.05 SC.visualizationSettings.connectors.defaultSize = 0.1 ####################### define nodes ######################################### nodeCounter = 0 # zero based for i in range(numberOfNodesY): for j in range(numberOfNodesX): nodeName = "node " + str(nodeCounter) # increase node counter nodeCounter = nodeCounter + 1 if j == 0: nodeDict = {"nodeType": "PointGround", "referenceCoordinates": [0.0, i*springDamperElementLength, 0.0], "name": nodeName} elif nodeCounter == (i+1)*numberOfNodesX and clambRightSide == True: nodeDict = {"nodeType": "PointGround", "referenceCoordinates": [j*springDamperElementLength, i*springDamperElementLength, 0.0], "name": nodeName} else: nodeDict = {"nodeType": "Point", "referenceCoordinates": [j*springDamperElementLength, i*springDamperElementLength, 0.0], "initialCoordinates": [0.0, 0.0, 0.0], "name": nodeName} # add node to mbs mbs.AddNode(nodeDict) ####################### create objects ####################################### for i in range(numberOfNodes): nodeDict = mbs.GetNode(i) nodeName = nodeDict["name"] nodeNumber = mbs.GetNodeNumber(nodeName) massPointName = "mass point - " + nodeName objectDict = {"objectType": "MassPoint", "physicsMass": massMassPoint, "nodeNumber": nodeNumber, "name": massPointName} mbs.AddObject(objectDict) ####################### add markers ########################################## bodyMassMarkerName = list() nodePositionMarkerName = list() for i in range(numberOfNodes): nodeDict = mbs.GetNode(i) nodeName = nodeDict["name"] nodeNumber = mbs.GetNodeNumber(nodeName) bodyMassMarkerName.append("body mass marker - " + nodeName) nodePositionMarkerName.append("node position marker - " + nodeName) # body mass - used for garvitational load bodyMassMarkerDict = {"markerType": "BodyMass", "bodyNumber": exu.ObjectIndex(nodeNumber), #nodeNumber=bodyNumber "name": bodyMassMarkerName[i]} nodePositionMarkerDict = {"markerType": "NodePosition", "nodeNumber": nodeNumber, "name": nodePositionMarkerName[i]} mbs.AddMarker(bodyMassMarkerDict) mbs.AddMarker(nodePositionMarkerDict) ####################### add loads ############################################ for i in range( len(bodyMassMarkerName) ): markerNumberOfBodyMassMarker = mbs.GetMarkerNumber(bodyMassMarkerName[i]) loadName = "gravity " + str(bodyMassMarkerName[i]) loadDict = {"loadType": "MassProportional", "markerNumber": markerNumberOfBodyMassMarker, "loadVector": [0.0, 0.0, -9.81], "name": loadName} mbs.AddLoad(loadDict) ####################### add connectors ####################################### ctr1 = 0 elementCtr = 0; # spring-damper elements in x-direction if useSpringDamperElementsX == True: for i in range(numberOfNodesY): for j in range(numberOfNodesX-1): markerNumberOfPositionMarkerL = mbs.GetMarkerNumber(nodePositionMarkerName[ctr1]) markerNumberOfPositionMarkerR = mbs.GetMarkerNumber(nodePositionMarkerName[ctr1+1]) ctr1 = ctr1 + 1 springDamperElementName = "spring damper" + str(elementCtr) objectDict = {"objectType": "ConnectorSpringDamper", "markerNumbers": [markerNumberOfPositionMarkerL, markerNumberOfPositionMarkerR], "stiffness": elementStiffness, "damping": elementDamping, "force": 0, "referenceLength": springDamperElementLength, "activeConnector": activateElement, "name": springDamperElementName, "VdrawSize": 0.0} mbs.AddObject(objectDict) elementCtr = elementCtr + 1 ctr1 = ctr1 + 1 # spring-damper elements in y-direction if useSpringDamperElementsY == True: ctr1 = 0 for i in range(numberOfNodesY-1): for j in range(numberOfNodesX): markerNumberOfPositionMarkerL = mbs.GetMarkerNumber(nodePositionMarkerName[ctr1]) markerNumberOfPositionMarkerR = mbs.GetMarkerNumber(nodePositionMarkerName[ctr1+numberOfNodesX]) ctr1 = ctr1 + 1 springDamperElementName = "spring damper" + str(elementCtr) objectDict = {"objectType": "ConnectorSpringDamper", "markerNumbers": [markerNumberOfPositionMarkerL, markerNumberOfPositionMarkerR], "stiffness": elementStiffness, "damping": elementDamping, "force": 0, "referenceLength": springDamperElementLength, "activeConnector": activateElement, "name": springDamperElementName, "VdrawSize": 0.0} mbs.AddObject(objectDict) elementCtr = elementCtr + 1 # spring-damper elements in off-diagnoal-direction if useSpringDamperElementsOffDiagonal == True: ctr1 = 0 for i in range(numberOfNodesY-1): #range(numberOfNodes): for j in range(numberOfNodesX-1): markerNumberOfPositionMarkerL = mbs.GetMarkerNumber(nodePositionMarkerName[ctr1]) markerNumberOfPositionMarkerR = mbs.GetMarkerNumber(nodePositionMarkerName[ctr1+numberOfNodesX+1]) ctr1 = ctr1 + 1 springDamperElementName = "spring damper" + str(elementCtr) objectDict = {"objectType": "ConnectorSpringDamper", "markerNumbers": [markerNumberOfPositionMarkerL, markerNumberOfPositionMarkerR], "stiffness": elementStiffness, "damping": elementDamping, "force": 0, "referenceLength": springDamperDiagnoalElementLength, "activeConnector": activateElement, "name": springDamperElementName, "VdrawSize": 0.0} mbs.AddObject(objectDict) elementCtr = elementCtr + 1 ctr1 = ctr1 + 1 # spring-damper elements in diagnoal-direction if useSpringDamperElementsDiagonal == True: ctr1 = 0 for i in range(numberOfNodesY-1): for j in range(numberOfNodesX-1): markerNumberOfPositionMarkerL = mbs.GetMarkerNumber(nodePositionMarkerName[ctr1+1]) markerNumberOfPositionMarkerR = mbs.GetMarkerNumber(nodePositionMarkerName[ctr1+numberOfNodesX]) ctr1 = ctr1 + 1 springDamperElementName = "spring damper" + str(elementCtr) objectDict = {"objectType": "ConnectorSpringDamper", "markerNumbers": [markerNumberOfPositionMarkerL, markerNumberOfPositionMarkerR], "stiffness": elementStiffness, "damping": elementDamping, "force": 0, "referenceLength": springDamperDiagnoalElementLength, "activeConnector": activateElement, "name": springDamperElementName, "VdrawSize": 0.0} mbs.AddObject(objectDict) elementCtr = elementCtr + 1 ctr1 = ctr1 + 1 ####################### assemble and solve mbs ############################### mbs.Assemble() #print(mbs) exu.StartRenderer() mbs.SolveDynamic(simulationSettings, solverType = exudyn.DynamicSolverType.ODE23) SC.WaitForRenderEngineStopFlag() exu.StopRenderer() #safely close rendering window!