ANCFslidingJoint.py
You can view and download this file on Github: ANCFslidingJoint.py
1#+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2# This is an EXUDYN example
3#
4# Details: 3D ANCF Cable element with sliding joint test
5#
6# Author: Johannes Gerstmayr
7# Date: 2026-03-02
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.utilities import ObjectANCFCable, VObjectANCFCable, InertiaCuboid, MarkerBodyRigid,\
15 NodeGenericData, ObjectJointSliding, MarkerBodyBeamShape
16import exudyn.graphics as graphics #only import if it does not conflict
17from exudyn.beams import GenerateStraightLineANCFCable
18
19import numpy as np
20
21#create an environment for mini example
22SC = exu.SystemContainer()
23mbs = SC.AddSystem()
24
25L = 2
26hg = 1 #ground offset y
27rCable = 0.01
28oGround=mbs.CreateGround(referencePosition= [0,0,0],
29 graphicsDataList=[graphics.CheckerBoard([0.5*L,-hg,0], normal=[0,1,0],size=3),
30 graphics.Cylinder([0,0,0],[0,-hg,0],radius=rCable*2, color=graphics.color.orange),
31 graphics.Cylinder([L,0,0],[0,-hg,0],radius=rCable*2, color=graphics.color.orange)
32 ])
33
34rhoA = 7800*rCable**2*np.pi
35EA = 0.25*100000.
36EI = 1
37
38nCables = 1
39for i in range(nCables):
40 p0 = np.array([0,0,i*0.1])
41 p1 = p0 + [L,0,0]
42
43 cable = ObjectANCFCable(physicsMassPerLength=rhoA,
44 physicsBendingStiffness = EI,
45 physicsBendingDamping = EI*0.02,
46 physicsAxialStiffness=EA,
47 physicsAxialDamping=EA*0.02,
48 visualization=VObjectANCFCable(radius = rCable),
49 )
50
51 ancf=GenerateStraightLineANCFCable(mbs=mbs,
52 positionOfNode0=p0, positionOfNode1=p1,
53 numberOfElements=48, #converged to 4 digits
54 cableTemplate=cable, #this defines the beam element properties
55 massProportionalLoad = [0,-9.81,0],
56 fixedConstraintsNode0 = [1,1,1, 0,1,1], #add constraints for pos and rot (r'_y,r'_z)
57 fixedConstraintsNode1 = [1,1,1, 0,1,1], #add constraints for pos and rot (r'_y,r'_z)
58 )
59 #ancf=[cableNodeList, cableObjectList, loadList, cableNodePositionList, cableCoordinateConstraintList]
60
61lElem = mbs.GetObject(ancf[1][0])['physicsLength']
62
63slidingCoordinateInit = 0.1*L
64initialLocalMarker = int(slidingCoordinateInit/lElem) #second element
65
66hy = 0.4 #height of rigid body
67sz = 5*rCable #z-displacement of rigid body
68cubeLengths = [5*rCable,hy-2*rCable,10*rCable]
69oRigid = mbs.CreateRigidBody(referencePosition=[slidingCoordinateInit, -0.5*hy,-sz],
70 inertia=InertiaCuboid(2*7800, cubeLengths),
71 gravity=[0,-9.81,0],
72 graphicsDataList=[graphics.Brick(size=cubeLengths, color=graphics.color.dodgerblue),
73 graphics.Sphere([0,0.5*hy,sz], radius=2*rCable, color=graphics.color.red)]
74 )
75mRigidTop = mbs.AddMarker(MarkerBodyRigid(bodyNumber=oRigid, localPosition=[0,0.5*hy,sz]))
76
77addSlidingJoint = True
78if addSlidingJoint:
79
80
81 cableMarkerList = []#list of MarkerBodyBeamShape
82 offsetList = [] #list of offsets counted from first cable element; needed in sliding joint
83 offset = 0 #first cable element has offset 0
84 for item in ancf[1]: #create markers for cable elements
85 m = mbs.AddMarker(MarkerBodyBeamShape(bodyNumber = item))
86 cableMarkerList += [m]
87 offsetList += [offset]
88 offset += lElem
89
90 nodeDataSJ = mbs.AddNode(NodeGenericData(initialCoordinates=[initialLocalMarker,slidingCoordinateInit],numberOfDataCoordinates=2)) #initial index in cable list
91 slidingJoint = mbs.AddObject(ObjectJointSliding(markerNumbers=[mRigidTop,cableMarkerList[initialLocalMarker]],
92 constrainRotations=[0,0,0],
93 slidingMarkerNumbers=cableMarkerList, slidingMarkerOffsets=offsetList,
94 nodeNumber=nodeDataSJ))
95
96
97#assemble and solve system for default parameters
98mbs.Assemble()
99
100endTime=10
101stepSize = 0.5e-3
102
103simulationSettings = exu.SimulationSettings()
104
105#simulationSettings.solutionSettings.writeSolutionToFile = False
106simulationSettings.solutionSettings.solutionWritePeriod = 0.02 #data not used
107simulationSettings.solutionSettings.binarySolutionFile = True
108simulationSettings.solutionSettings.outputPrecision = 6 #float
109simulationSettings.solutionSettings.sensorsWritePeriod = 0.002 #data not used
110simulationSettings.timeIntegration.verboseMode = 1 #turn off, because of lots of output
111simulationSettings.linearSolverType = exu.LinearSolverType.EigenSparse
112simulationSettings.parallel.numberOfThreads = 1
113simulationSettings.displayComputationTime = True
114simulationSettings.displayStatistics = True
115
116simulationSettings.timeIntegration.numberOfSteps = int(endTime/stepSize)
117simulationSettings.timeIntegration.endTime = endTime
118simulationSettings.timeIntegration.newton.useModifiedNewton = True
119
120#simulationSettings.timeIntegration.simulateInRealtime = True
121#simulationSettings.timeIntegration.realtimeFactor = 0.5
122
123SC.visualizationSettings.openGL.multiSampling = 4
124
125SC.visualizationSettings.general.graphicsUpdateInterval = 0.02
126SC.visualizationSettings.view0.window.renderWindowSize=[1200,1024]
127SC.visualizationSettings.view0.camera.perspective = 0.5
128SC.visualizationSettings.openGL.light0.shadow = 0.3
129SC.visualizationSettings.openGL.light0.position = [2,10,2,1]
130SC.visualizationSettings.nodes.show = False
131SC.visualizationSettings.loads.show = False
132SC.visualizationSettings.connectors.show = False
133#+++++++++++++++++++++++++++++++++++++++++++++++++++++++
134
135SC.renderer.Start()
136SC.renderer.SetModelView(zoom=1.335652,rotationVector=[0.7933899,0.5903195,0.2917689],centerPoint=[1.012,-0.349,0])
137SC.renderer.DoIdleTasks()
138
139mbs.SolveDynamic(simulationSettings)
140
141SC.renderer.DoIdleTasks()
142SC.renderer.Stop() #safely close rendering window!
143
144mbs.SolutionViewer()