pickleCopyMbs.py

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

  1#%%+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  2# This is an EXUDYN example
  3#
  4# Details:  Two rigid bodies, where the first body moves on a straight rail (fixed on ground) and the second is mounted with a revolute joint on the first body;
  5#           the particular reason for this test example is the copying of the model as well as load save;
  6#           shows reading and writing mbs with dictionaries and export by HDF5 files
  7#
  8# Author:   Johannes Gerstmayr
  9# Date:     2024-10-13
 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 exudyn as exu
 16from exudyn.utilities import * #includes itemInterface and rigidBodyUtilities
 17import exudyn.graphics as graphics #only import if it does not conflict
 18from math import sin, cos
 19import copy
 20
 21useGraphics = True #without test
 22#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 23#you can erase the following lines and all exudynTestGlobals related operations if this is not intended to be used as TestModel:
 24try: #only if called from test suite
 25    from modelUnitTests import exudynTestGlobals #for globally storing test results
 26    useGraphics = exudynTestGlobals.useGraphics
 27except:
 28    class ExudynTestGlobals:
 29        pass
 30    exudynTestGlobals = ExudynTestGlobals()
 31#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 32
 33
 34SC = exu.SystemContainer()
 35mbs = SC.AddSystem()
 36
 37
 38#%%+++++++++++++++++++++++++++++++++++
 39L = 0.5
 40w = 0.05
 41ground = mbs.CreateGround(graphicsDataList=[graphics.CheckerBoard(point=[0,0,-0.1*L], size=10)])
 42
 43b0size = [0.5*L,L,0.2*L]
 44p0 = np.array([0.5,1,0])
 45b0 = mbs.CreateRigidBody(referencePosition=p0,
 46                         initialVelocity=[0,1,0],
 47                         inertia=InertiaCuboid(500, b0size),
 48                         graphicsDataList=[graphics.Brick(size=b0size,
 49                                                          color=graphics.color.dodgerblue)])
 50b1size=[w,L,w]
 51p1 = p0 + [0,0.5*L,0.5*b0size[2]+w]
 52b1 = mbs.CreateRigidBody(referencePosition=p1,
 53                         initialVelocity=[0,1,0],
 54                         initialAngularVelocity=[0,0,4*pi],
 55                         inertia=InertiaCuboid(1000, b1size),
 56                         graphicsDataList=[graphics.Brick(size=b1size,
 57                                                          color=graphics.color.red)])
 58
 59mbs.CreateRevoluteJoint(bodyNumbers=[b0,b1],
 60                        position=p1-[0,0.5*L,0.5*w],
 61                        axis = [0,0,1],
 62                        axisRadius=0.5*w, axisLength=2.2*w)
 63
 64mbs.CreatePrismaticJoint(bodyNumbers=[ground,b0],
 65                         position=p0,
 66                         axis = [0,1,0],
 67                         axisRadius=0.5*w, axisLength=5*L)
 68
 69# UFtorque = 0
 70def UFtorque(mbs, t, loadVector):
 71    f = np.heaviside(t - 1,1)
 72    return [0,0,-f*2*cos(pi*2*t)]
 73
 74load = mbs.CreateTorque(bodyNumber=b1, loadVector=[0,0,0],
 75                        loadVectorUserFunction=UFtorque)
 76
 77sPos = mbs.AddSensor(SensorBody(bodyNumber=b0, storeInternal=True,
 78                                outputVariableType=exu.OutputVariableType.Position))
 79
 80
 81#del mbs
 82
 83#%%+++++++++++++++++++++++++++++++++++
 84#use dictionaries and pickle to load/save MBS
 85import pickle
 86
 87mbsDict = mbs.GetDictionary()
 88
 89#save MUST be done before solving, as solver structures and settings cannot be saved!
 90#erase variables and system variables as they may contain some solver stuff (if solved once)
 91#  that should not / can not be pickled
 92mbsDict['variables'] = {}
 93mbsDict['systemVariables'] = {}
 94
 95#save mbs data (could also use dill)
 96with open('solution/mbs.pkl', 'wb') as f:
 97    pickle.dump(mbsDict, f) #, pickle.HIGHEST_PROTOCOL)
 98
 99#load mbs data
100with open('solution/mbs.pkl', 'rb') as f:
101    mbsCopy = pickle.load(f)
102
103SC = exu.SystemContainer()
104mbs2 = SC.AddSystem()
105mbs2.SetDictionary(mbsDict)
106
107#+++++++++++++++++++++++
108
109hasH5py = False
110try:
111    import h5py
112    hasH5py = True
113except ImportError:
114    exu.Print('h5py not available; skipping HDF5 test')
115
116if hasH5py:
117    from exudyn.advancedUtilities import SaveDictToHDF5, LoadDictFromHDF5
118    #use HDF5 load / save may not work for all cases, as some types are not implemented
119    #save MUST be done before solving, as solver structures and settings cannot be saved!
120    mbsDict = mbs.GetDictionary()
121    mbsDict['variables'] = {}
122    mbsDict['systemVariables'] = {}
123    SaveDictToHDF5('solution/mbs.h5', mbsDict) #problems with Python functions as sub-dicts
124    mbsCopy2 = LoadDictFromHDF5('solution/mbs.h5', globals())
125    # print('loaded_data:\n', mbsCopy2)
126
127    SC = exu.SystemContainer()
128    mbs2 = SC.AddSystem()
129    mbs2.SetDictionary(mbsCopy)
130
131#%%+++++++++++++++++++++++++++++++++++
132#ALTERNATIVE: work with copy of mbs:
133#we have to create a new Systemcontainer:
134#mbs2 = copy.copy(mbs) #alternative way to directly copy
135# SC = exu.SystemContainer()
136# SC.AppendSystem(mbs2)
137
138mbs2.Assemble()
139
140#check mbs2:
141#print(mbs2)
142
143simulationSettings = exu.SimulationSettings() #takes currently set values or default values
144
145tEnd = 1
146h = 1e-3
147simulationSettings.timeIntegration.numberOfSteps = int(tEnd/h)
148simulationSettings.timeIntegration.endTime = tEnd
149simulationSettings.timeIntegration.verboseMode = 1
150# simulationSettings.timeIntegration.simulateInRealtime = True
151
152simulationSettings.timeIntegration.newton.useModifiedNewton = True
153simulationSettings.timeIntegration.generalizedAlpha.spectralRadius = 0.6 #0.6 works well
154# simulationSettings.displayComputationTime = True
155
156simulationSettings.solutionSettings.solutionInformation = "rigid body tests"
157SC.visualizationSettings.nodes.defaultSize = 0.05
158SC.visualizationSettings.openGL.multiSampling = 4
159SC.visualizationSettings.openGL.lineWidth = 2
160SC.visualizationSettings.window.renderWindowSize = [2000,1600]
161
162if useGraphics:
163    SC.renderer.Start()
164    mbs2.WaitForUserToContinue()
165
166mbs2.SolveDynamic(simulationSettings)
167
168if useGraphics:
169    SC.renderer.DoIdleTasks()
170    SC.renderer.Stop() #safely close rendering window!
171
172
173#+++++++++++++++++++++++++++++++++++++++++++++
174uTotal = 0.1 * sum(mbs2.GetSensorValues(sPos)) #sPos is still the correct index
175exu.Print('uTotal=',uTotal)
176
177exudynTestGlobals.testResult = uTotal
178#+++++++++++++++++++++++++++++++++++++++++++++