taskmanagerTest.py

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

  1#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  2# This is an EXUDYN example
  3#
  4# Details:  Systematic tests with several types of models, object types, and constraints;
  5#           test for implicit, explicit and static solvers for number of threads in task manager
  6#
  7# Author:   Johannes Gerstmayr
  8# Date:     2024-10-08
  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 * #includes itemInterface and rigidBodyUtilities
 16import exudyn.graphics as graphics
 17
 18import numpy as np
 19
 20testIsActive = exu.sys.get('testIsActive', False)
 21
 22testIsActive = True
 23
 24SC = exu.SystemContainer()
 25mbs = SC.AddSystem()
 26
 27resultTotal = 0
 28#general dimentions
 29a = 1
 30b = 0.5
 31g = [0,-9.81,0]
 32mass = 0.2
 33stiffness = 1000
 34damping = 2
 35load = 2
 36totalCases = 0
 37
 38#user function for spring force
 39def springForce(mbs, t, itemIndex, u, v, k, d, offset): #changed 2023-01-21:, mu, muPropZone):
 40    return k*u+k*u**3+v*d
 41
 42#user function for load
 43def userLoadVector(mbs, t, loadVector):
 44    return np.sin(4*pi*t)*np.array(loadVector)
 45
 46nMasses = 1
 47caseConstraints = 0
 48caseODE1 = 0
 49caseUser = 0
 50caseStatic = 0
 51nThreads = 1
 52
 53for nMasses in [0, 1, 2, 4]:
 54    for caseConstraints in range(2):
 55        for caseStatic in range(2):
 56            for caseUser in range(2):
 57                for caseODE1 in range(2):
 58                    for nThreads in [1,2,5]:
 59                        mbs.Reset()
 60                        
 61                        if caseStatic and caseODE1: continue
 62                        if nMasses == 0 and not caseODE1: continue
 63                        
 64                        useExplicitSolver = (nMasses == 2 and caseConstraints == 0 and caseStatic == 0)
 65
 66                        if not caseUser:
 67                            springForce = 0
 68                            
 69                        case = [nThreads, caseODE1, nMasses, caseConstraints, caseUser, 
 70                                caseStatic, int(useExplicitSolver)]
 71                        
 72                        ground0 = mbs.CreateGround()
 73                        ground1 = mbs.CreateGround(referencePosition=[0,b,0])
 74                        bLast0 = ground0
 75                        bLast1 = ground1
 76                        
 77                        for i in range(nMasses):
 78                            xOff = i*a
 79                            b0=mbs.CreateMassPoint(referencePosition=[xOff+a, 0, 0],
 80                                                mass=mass,
 81                                                drawSize = 0.1*a, color=graphics.color.red,
 82                                                create2D=True,
 83                                                gravity = g)
 84                            b1=mbs.CreateMassPoint(referencePosition=[xOff+a, b, 0],
 85                                                mass=mass,
 86                                                drawSize = 0.1*a, color=graphics.color.red,
 87                                                create2D=True,
 88                                                gravity = g)
 89                        
 90                            mbs.CreateSpringDamper(itemNumbers=[bLast1, b1], stiffness=stiffness, damping=damping,
 91                                                   springForceUserFunction=springForce)
 92                            if caseConstraints:
 93                                mbs.CreateDistanceConstraint(itemNumbers=[bLast0, b0])
 94                                mbs.CreateDistanceConstraint(itemNumbers=[bLast1, b0])
 95                                mbs.CreateDistanceConstraint(itemNumbers=[b0, b1])
 96                            else:
 97                                mbs.CreateSpringDamper(itemNumbers=[bLast0, b0], stiffness=stiffness, damping=damping)
 98                                mbs.CreateSpringDamper(itemNumbers=[bLast1, b0], stiffness=stiffness, damping=damping)
 99                                mbs.CreateSpringDamper(itemNumbers=[b0, b1], stiffness=stiffness, damping=damping)
100                            
101                            bLast0 = b0
102                            bLast1 = b1
103                        
104                        if caseUser:
105                            mbs.CreateForce(itemNumber=bLast1, 
106                                            loadVector=[10,20,0],
107                                            loadVectorUserFunction=userLoadVector)
108                        
109                        sLast = mbs.AddSensor(SensorBody(bodyNumber=bLast0, storeInternal=True,
110                                              outputVariableType=exu.OutputVariableType.Displacement))
111                        
112                        sODE1 = None
113                        if caseODE1:
114                            #set up a 2-DOF system
115                            nODE1 = mbs.AddNode(NodeGenericODE1(referenceCoordinates=[0,0,1],
116                                                                initialCoordinates=[1,0,0],
117                                                                numberOfODE1Coordinates=3))
118                            
119                            #build system matrix and force vector
120                            #undamped mechanical system with m=1, K=100, f=1
121                            #additionally 1 ODE1 coordinate
122                            Amat = np.array([[0,1,0],
123                                          [-100,0,0],
124                                          [0,0,-10]])
125                            bVec = np.array([0,1,50])
126                            
127                            oGenericODE1 = mbs.AddObject(ObjectGenericODE1(nodeNumbers=[nODE1],
128                                                                           systemMatrix=Amat,
129                                                                           rhsVector=bVec))
130                            sODE1 = mbs.AddSensor(SensorNode(nodeNumber=nODE1, storeInternal=True,
131                                                  outputVariableType=exu.OutputVariableType.Coordinates))
132                        
133                        
134                        mbs.Assemble()
135                        
136                        tEnd = 0.1*20
137                        h = 0.01
138                        simulationSettings = exu.SimulationSettings()
139                        simulationSettings.timeIntegration.newton.useModifiedNewton = False #Just for the test; modified Newton is usually faster
140                        simulationSettings.timeIntegration.numberOfSteps = int(tEnd/h)
141                        simulationSettings.timeIntegration.endTime = tEnd
142                        # simulationSettings.show.statistics = True
143                        # simulationSettings.show.computationTime = True
144                        simulationSettings.timeIntegration.verboseMode = 0
145                        simulationSettings.staticSolver.verboseMode = 0
146                        
147                        simulationSettings.parallel.numberOfThreads = nThreads
148                        
149                        SC.visualizationSettings.nodes.drawNodesAsPoint = False
150                        SC.visualizationSettings.nodes.tiling = 32
151                        SC.visualizationSettings.view0.window.alwaysOnTop = True
152                        
153                        #start solver:
154                        if not testIsActive:
155                            SC.renderer.Start()
156                            SC.renderer.DoIdleTasks()
157                        
158                        if caseStatic:
159                            mbs.SolveStatic(simulationSettings)
160                        else:
161                            if not useExplicitSolver:
162                                mbs.SolveDynamic(simulationSettings)
163                            else:
164                                solverType = exu.DynamicSolverType.RK44
165                                if caseUser and not caseODE1:
166                                    solverType = exu.DynamicSolverType.VelocityVerlet
167
168                                mbs.SolveDynamic(simulationSettings,
169                                                 solverType=solverType)
170                        
171                        if not testIsActive:
172                            SC.renderer.Stop()
173                        
174                        #+++++++++++++++++++++++++++++++++++++
175                        #evaluate results
176                        x = mbs.GetSensorValues(sLast)
177                        resultTotal += sum(x)
178                        output = 'tip-disp='+str(list(x))
179                        
180                        if caseODE1:
181                            x = mbs.GetSensorValues(sODE1)
182                            output += ', ODE1='+str(list(x))
183                            resultTotal += sum(x)
184                        
185                        if True:
186                            exu.Print('** CASE '+str(totalCases)+' **:',case,
187                                      '; RESULTS:',output)
188
189                        totalCases += 1
190
191#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
192#evaluate final (=current) output values
193resultTotal *= 0.001 #to avoid too large non-deterministic round-off errors due to multithreading
194exu.Print('\ntotal cases:', totalCases)
195exu.Print('result taskmanagerTest=', resultTotal)
196
197exu.sys['testResult'] = resultTotal #-0.17210771618100057 
198
199# mbs.SolutionViewer()
200
201#%%++++++++++++++++
202if False:
203    mbs.PlotSensor(sLast, components=[1])
204    if sODE1 != None:
205        mbs.PlotSensor(sODE1, components=[0,1,2])