-
Notifications
You must be signed in to change notification settings - Fork 0
/
teste_MBS_rodeiro.py
98 lines (68 loc) · 2.19 KB
/
teste_MBS_rodeiro.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Thu Mar 17 06:48:46 2022
@author: leonardo
"""
import numpy as np
import MBS.MultibodySystem as mbs
import MBS.BodyConnections.Forces
import MBS.marker
import profiles
from bodiesc import rigidBody
from assimulo.solvers import IDA
system = mbs.MultibodySystem('Single wheelset')
'''
Wheelset definition
'''
wst_mass = 650.0
wst_Iy = 70.
wst = rigidBody('Wheelset')
wst.setMass(wst_mass)
wst.setInertiaTensor(np.diag([wst_Iy, wst_Iy, wst_Iy]))
profileMarker1 = wst.addMarker(MBS.marker.marker('Profile 1',np.array([0.56,-0.47,0.0])))
profileMarker2 = wst.addMarker(MBS.marker.marker('Profile 2',np.array([-0.56,-0.47,0.0])))
railMarker1 = system.ground.addMarker(MBS.marker.marker('Profile rail 1',np.array([0.464,-0.66,0.0])))
railMarker2 = system.ground.addMarker(MBS.marker.marker('Profile rail 2',np.array([-0.464,-0.66,0.0])))
'''
Profiles
'''
wst.addProfile(profiles.planarProfile('Wheel 1','./design2.pro', convPar = 1),profileMarker1)
wst.addProfile(profiles.planarProfile('Wheel 2','./design2.pro', convPar = 1),profileMarker2)
wst.profiles[0].rotatePoints(np.pi)
wst.profiles[1].rotatePoints(np.pi)
wst.profiles[1].mirrorVert()
system.ground.addProfile(profiles.planarProfile('Rail 1','./tr68.pro', convPar = 1),railMarker1)
system.ground.addProfile(profiles.planarProfile('Rail 2','./tr68.pro', convPar = 1),railMarker2)
'''
Forces
'''
f = MBS.BodyConnections.Forces.linearSpring_PtP('Mola',1e4)
f.connect(wst,system.ground)
'''
Multibody system definition
'''
system.addBody(wst)
system.addForce(f)
system.gravity = np.array([0,-9.8,0],dtype=np.float64)
system.setupSystem()
problem = system.generate_problem('ind3')
DAE = IDA(problem)
DAE.report_continuously = True
DAE.inith = 1e-5
DAE.num_threads = 12
DAE.suppress_alg = True
outFreq = 10e3 # Hz
finalTime = 1.2
problem.res(0,problem.y0,problem.yd0)
t,p,v=DAE.simulate(finalTime, finalTime * outFreq)
system.postProcess(t, p, v)
q = p[:,:system.n_p]
u = p[:,system.n_p:2*system.n_p]
lam = p[:,2*system.n_p:]
import matplotlib.pyplot as plt
ax = plt.axes()
wst.profiles[0].plotMe(ax)
wst.profiles[1].plotMe(ax)
system.ground.profiles[0].plotMe(ax)
system.ground.profiles[1].plotMe(ax)