Both the tractor and semi-trailer are assumed to be rigid bodies connected at a hitching point (node 30). This hitching point imposes constraints on the motion of the rigid bodies (nodes 10 for the tractor and 20 for the trailer). As a result, the vehicle model exhibits eight independent degrees of freedom: vertical bouncing (yT) and pitching (θT) motions of the tractor’s center of gravity (node 10), pitching motion (θS) of the semi-trailer’s center of gravity (node 20), and vertical hop motions (y1, y2, y3, y4, y5) of the five axle assemblies at nodes 1, 101, 201, 301, and 401.
The constraint is yS (vertical bouncing of the trailer) = yT +h1θT +h3θS, where h1 (from node 10 to node 30) and h2 (from node 20 to node 30) are the distances from center of gravity to the hitching point for the tractor and the trailer.
https://imgur.com/vrAfAwfIn the provided code, I attempted to model the hitching poind as a hinge connection between the two systems using two zero-length elements with the vertical degree of freedom released at the support (node 32). The general principle works, but the issue arises when changing the stiffness properties of the spring elements used for the hinge. This leads to variations in the reaction forces.
Below are the results from a static analysis, showing the obtained reaction forces at the nodes compared to their respective reference values:
Attached is the code I used for this analysis.Node 0: Reaction = 44559.36 N (Reference = 39090 N)
Node 100: Reaction = 85718.90 N (Reference = 78750 N)
Node 200: Reaction = 61634.03 N (Reference = 65113 N)
Node 300: Reaction = 47806.91 N (Reference = 51953 N)
Node 400: Reaction = 33979.79 N (Reference = 38793 N)
Node 32: Reaction = 0.00 N
Code: Select all
import openseespy.opensees as ops
import opsvis as opsv
import matplotlib.pyplot as plt
import numpy as np
# =============================================================================
# Material Definitions
# =============================================================================
g = 9.81 # Gravity
# Material properties for tires and suspension
ku1 = 175.e4 # Tire stiffness
ku2 = 175.e4
ku3 = 350.e4
ku4 = 350.e4
ku5 = 350.e4
ks1, cs1 = 6.e6, 1.e4 # Suspension stiffness and damping
ks2, cs2 = 6.e6, 1.e4
ks3, cs3 = 10.e6, 2.e4
ks4, cs4 = 10.e6, 2.e4
ks5, cs5 = 10.e6, 2.e4
# Mass properties
mb1 = 3100
mb2 = 20000
ib1 = 5000
ib2 = 120000
mu1 = 750
mu2 = 750
mu3 = 1100
mu4 = 1100
mu5 = 1100
# =============================================================================
# # Model
# =============================================================================
ops.wipe()
ops.model('basic', '-ndm', 2, '-ndf', 3)
# Material
ops.uniaxialMaterial('Elastic', 1, ku1)
ops.uniaxialMaterial('Elastic', 2, ks1, cs1)
ops.uniaxialMaterial('Elastic', 3, ku2)
ops.uniaxialMaterial('Elastic', 4, ks2, cs2)
ops.uniaxialMaterial('Elastic', 5, ku3)
ops.uniaxialMaterial('Elastic', 6, ks3, cs3)
ops.uniaxialMaterial('Elastic', 7, ku4)
ops.uniaxialMaterial('Elastic', 8, ks4, cs4)
ops.uniaxialMaterial('Elastic', 9, ku5)
ops.uniaxialMaterial('Elastic', 10, ks5, cs5)
# =============================================================================
# Model Generation
# =============================================================================
x_ax1 = 0 # Axle 1
x_ax2 = 4.5
x_ax3 = 11
x_ax4 = 12.2
x_ax5 = 13.4 # Axle 5
x_mb1 = 1 # Center of mass 1 (tractor body)
x_mb2 = 9 # Centar of mass 2 (trailer)
x_hg = 4 # Hinge location
### Tractor Body
# Axle 1 nodes and elements
ops.node(0, x_ax1, 0)
ops.fix(0, 1, 1, 1)
ops.node(1, x_ax1, 0)
ops.mass(1, 0, mu1, 0)
ops.node(2, x_ax1, 0)
ops.element('zeroLength', 1, 0, 1, '-mat', 1, 1, '-dir', 1, 2)
ops.element('zeroLength', 2, 1, 2, '-mat', 2, 2, '-dir', 1, 2, '-doRayleigh')
# Axle 2 nodes and elements
ops.node(100, -x_ax2, 0)
ops.fix(100, 1, 1, 1)
ops.node(101, -x_ax2, 0)
ops.mass(101, 0, mu2, 0)
ops.node(102, -x_ax2, 0)
ops.element('zeroLength', 3, 100, 101, '-mat', 3, 3, '-dir', 1, 2)
ops.element('zeroLength', 4, 101, 102, '-mat', 4, 4, '-dir', 1, 2, '-doRayleigh')
# Center node and rigid links (tractor body)
ops.node(10, -x_mb1, 0)
ops.mass(10, 0, mb1, ib1)
ops.node(11, -x_mb1, 0)
ops.rigidLink('-beam', 10, 2)
ops.rigidLink('-beam', 10, 102)
ops.rigidLink('-beam', 11, 1)
ops.rigidLink('-beam', 11, 101)
### Trailer
# Axle 3 nodes and elements
ops.node(200, -x_ax3, 0)
ops.fix(200, 1, 1, 1)
ops.node(201, -x_ax3, 0)
ops.mass(201, 0, mu3, 0)
ops.node(202, -x_ax3, 0)
ops.element('zeroLength', 5, 200, 201, '-mat', 5, 5, '-dir', 1, 2)
ops.element('zeroLength', 6, 201, 202, '-mat', 6, 6, '-dir', 1, 2, '-doRayleigh')
# Axle 4 nodes and elements
ops.node(300, -x_ax4, 0)
ops.fix(300, 1, 1, 1)
ops.node(301, -x_ax4, 0)
ops.mass(301, 0, mu4, 0)
ops.node(302, -x_ax4, 0)
ops.element('zeroLength', 7, 300, 301, '-mat', 7, 7, '-dir', 1, 2)
ops.element('zeroLength', 8, 301, 302, '-mat', 8, 8, '-dir', 1, 2, '-doRayleigh')
# Axle 5 nodes and elements
ops.node(400, -x_ax5, 0)
ops.fix(400, 1, 1, 1)
ops.node(401, -x_ax5, 0)
ops.mass(401, 0, mu5, 0)
ops.node(402, -x_ax5, 0)
ops.element('zeroLength', 9, 400, 401, '-mat', 9, 9, '-dir', 1, 2)
ops.element('zeroLength', 10, 401, 402, '-mat', 10, 10, '-dir', 1, 2, '-doRayleigh')
# Center node and rigid links (trailer)
ops.node(20, -x_mb2, 0)
ops.mass(20, 0, mb2, ib2)
ops.node(21, -x_mb2, 0)
ops.rigidLink('-beam', 20, 202)
ops.rigidLink('-beam', 20, 302)
ops.rigidLink('-beam', 20, 402)
ops.rigidLink('-beam', 21, 201)
ops.rigidLink('-beam', 21, 301)
ops.rigidLink('-beam', 21, 401)
# Rigid link to connect the tractor and trailer
ops.node(30, -x_hg, 0)
ops.node(31, -x_hg, 0)
ops.node(32, -x_hg, 0)
ops.fix(32, 1, 0, 1)
ops.element('zeroLength', 11, 32, 30, '-mat', 6, 6, '-dir', 1, 2)
ops.element('zeroLength', 12, 31, 30, '-mat', 6, 6, '-dir', 1, 2)
ops.rigidLink('-beam', 11, 31)
ops.rigidLink('-beam', 20, 30)
ops.rigidLink('-beam', 21, 31)
# Visualize model
opsv.plot_model()
plt.show()
# Define full load applied to the center node
load1_x = 0
load1_y = -1 * (mu1+mu2+mb1) * g # Total weight converted to force
load2_x = 0
load2_y = -1 * (mu3+mu4+mu5+mb2) * g
ops.timeSeries('Constant', 1) # Use constant load over the analysis step
ops.pattern('Plain', 1, 1) # Load pattern
ops.load(10, load1_x, load1_y, 0) # Apply load to center node
ops.load(20, load2_x, load2_y, 0)
# Static analysis settings
ops.test('NormDispIncr', 1e-3, 10, 0)
ops.algorithm('KrylovNewton')
ops.system('UmfPack')
ops.constraints('Transformation')
ops.integrator('LoadControl', 1.0)
ops.analysis('Static')
# Perform static analysis
ok = ops.analyze(1) # Perform analysis in one step
ops.reactions()
# Reaction tracking
fixed_node_tags = [0, 100, 200, 300, 400, 32] # Nodes with fixed supports
reactions = {tag: ops.nodeReaction(tag, 2) for tag in fixed_node_tags}
# Analysis results
if ok < 0:
print("Static analysis failed.")
else:
print("Static analysis successful.")
for tag, reaction in reactions.items():
print(f"Node {tag}: Reaction = {reaction:.2f} N")