Problem Set-1 c Week 11Challenge : Transient Handling
the following steps are followed for both FWD and F1:
- first input provieded
- A matrix is inputted using A=[-(cf+cr)/(m*u),(((b*cr)-(a*cf))/(m*u))-u;((b*cr)-(a*cf))/(jz*u),-((((a^2)*cf)+((b^2)*cr))/(jz*u))]
- b matrix is inpurted usihg B=[cf/m;(a*cf)/jz]
- output for the lateral speed specified using cl and dz matrix
cl=[10]
di=0 - output for yaw rate is specified usingc2 and d2 matrix:
C2=[0,1]
d2=0 - output for lateral acceleration is specified using C3 and d3 matrix: c3=A(1;)+u*[0,1]
d3=B(1)
where A(1,:) means the First row in A matrix with all the columns and b(1) means the first element in the b matrıx - time array is form 0 and 2 with 0.01 steps in between
- sinusoidal steering input function is declared using st=(pi/180)*5*Sin(2*pi*7*t)
- then all the outputs are combined in and C and A column array
- a state-space system is formed from the A, B, C.D matrix
- linear simulation is performed on the system using the Lsim command to simulate the system dynamic behavior with time with the steering function and time as an input, this
the linear simulation will be performed on all the outputs which will simulate them concerning time with the steering function as input - the inputs and outputs of the state space system are named using input and output channels naming commands. sys. input name and sys. Output name
- bode command is used to store the gain magnitude phase and Frequency In an array
- to turn the 3-dimensional matrix Of the magnitude intO one dimension, a for loop Is used to store the magnitude value from the 3-dimensional matrix into one
dimensional array - it can be plotted using a log scale with frequency and the frequency is changed to Hz
- subplot command with 2 rows and 2 columns with the first element is the sinusoidal Output for steering with time
- subplot command with 2 rows and 2 columns with the second element is the sinusoidal Output for lateral speed which is plotted using the first column of y array with all the rows
- subplot command with 2 rows and 2 columns with the third element is the sinusoidal output for lateral speed which is plotted using the second column of y array with all the rows
- subplot command with 2 rows and 2 columns with the fourth element is the sinusoidal output for lateral speed which is plotted using the third column of y array with all the rows
- step response of the system was plotted using step command
- all the step response parameters are stored in an array s using the step info command, this will be used for correlation with analytical calculations
- using the damp command Matlab automatically calculates the natural frequency and the damping ratio of the system and stores them in an array
- calculating the values was done using
- the following coefficient is calculated using
YB=Cf+Cr
YR=(a*Cf-b*Cr)/u
YS=-Cf
NB=a*Cf-b*Cr
NR=(((a^2)*Cf}+((b^2)^Cr))/u
NS=-a*Cf - wn is calculated using
Wn=sqrt(((YB*NR)+(NB*m*u)-(YR*NB))/(m*u*jz)) - damping ratio theta is calculated using
Theta=abs(-((NR/jz)+(YB/(m*u))*(1/2*Wn))) - damped natural frequency Is Calculated using
Wd=Wn*sqrt(1-(Theta)^2) - Then rising time is calculated using
tr=((m*u*NS)/(YS*NB)-(YB*NS)))
M script for FWD passenger car:
clear all
close all
clc
l = 2.538984 ; % Wheelbase (m)
a = 0.890016; % cg to front axle (m)
b = l-a; % cg to rear axle (m)
m = 907.18474; % Mass (kg)
jz = 1128 ; % moment of inertia (kg-m^2)
cf = 61171.9437; % front cornering stiffness (N/rad)
cr = 50976.6191; % rear cornering stiffness (N/rad)
u = 100/3.6; % m/s
A=[-(cf+cr)/(m*u),(((b*cr)-(a*cf))/(m*u))-u;((b*cr)-(a*cf))/(jz*u),-((((a^2)*cf)+((b^2)*cr))/(jz*u))]
B=[cf/m;(a*cf)/jz]
c1=[1,0]; %lateral speed
d1=0;
c2=[0,1]; % yaw rate
d2=0;
c3=A(1,:)+u*[0,1]; %lateral acceleration
d3=B(1);
t=[0:0.01:2]
st=(pi/180)*5*sin(2*pi*7*t)
C=[c1;c2;c3]
D=[d1;d2;d3]
sys=ss(A,B,C,D)
Y=lsim(sys,st,t)
sys.Inputname={'st'}
sys.Outputname={'lateral speed','yaw rate','lateral acceleration'}
[mag1,phase1,omega1] = bode(sys({'yaw rate'},{'st'}))
for i=1:length(mag1)
magnitude1(i)=mag1(1,1,i)
i=i+1
end
figure(1)
loglog(omega1/(2*pi),magnitude1,'linewidth',2)
title('yaw gain frequency response')
xlabel('frequency(Hz)')
ylabel('gain')
grid on
grid on
figure(2)
subplot(2,2,1)
plot(t,st,'linewidth',2)
title('steering')
xlabel('time(s)')
ylabel('steering angle(rad)')
grid on
subplot(2,2,2)
plot(t,Y(:,1),'linewidth',2)
title('lateral speed')
xlabel('time(s)')
ylabel('lateral speed(m/s)')
grid on
subplot(2,2,3)
plot(t,Y(:,2),'linewidth',2)
title('yawrate')
xlabel('time(s)')
ylabel('yaw rate(rad/sec)')
grid on
subplot(2,2,4)
plot(t,Y(:,3),'linewidth',2)
title('lateral accelaration')
xlabel('time(s)')
ylabel('lateral acceleration(m/s^2)')
grid on
figure(3)
step(sys)
grid on
S=stepinfo(sys)
[wn theta]=damp(sys)
Cf = -61171.9437; % front cornering stiffness (N/rad)
Cr = -50976.6191
YB=Cf+Cr
YR=(a*Cf-b*Cr)/u
YS=-Cf
NB=a*Cf-b*Cr
NR=(((a^2)*Cf)+((b^2)*Cr))/u
NS=-a*Cf
Wn=sqrt(((YB*NR)+(NB*m*u)-(YR*NB))/(m*u*jz))
Theta=abs(-((NR/jz)+(YB/(m*u)))*(1/(2*Wn)))
Wd=Wn*sqrt(1-(Theta)^2)
tr=((m*u*NS)/((YS*NB)-(YB*NS)))
Correlation at 7hz:
Gain (from yaw rate and steering rate) =0.0935/0.0855 = 1.093
From frequency response, it can be seen that gain is approx. 1.1
Steady-state yaw gain refers to the gain of the system when t approaches infinity or at extremely low frequency, for example, in the step response, it achieves a steady-state output value of 5.4 for a unity step input, therefore the steady-state gain is 5.4.
rise time correction;
There is a minor difference between the analytical calculation of rising time and the automated calculation using the step info function in MATLAB.
Resut:
YB=-1.1258*10^5
YR=6.1172*10^4
YS=6.1172*10^4
NB=2.9615*10^4
NR=-6.7344*10^3
NS=5.4444*10^4
Wn=7.1912 rad/s
Theta=0.72245
Wd=4.9564 rad/s
tr=0.1755 S
M script for F1 car:
clear all
close all
clc
l = 2.901696 ; % Wheelbase (m)
a = 1.740408; % cg to front axle (m)
b = l-a; % cg to rear axle (m)
m = 779.725284; % Mass (kg)
jz = 865 ; % moment of inertia (kg-m^2)
cf = 276803.935; % front cornering stiffness (N/rad)
cr = 429733.794; % rear cornering stiffness (N/rad)
u = 100/3.6; % m/s
A=[-(cf+cr)/(m*u),(((b*cr)-(a*cf))/(m*u))-u;((b*cr)-(a*cf))/(jz*u),-((((a^2)*cf)+((b^2)*cr))/(jz*u))]
B=[cf/m;(a*cf)/jz]
c1=[1,0] %lateral speed
d1=0
c2=[0,1] % yaw rate
d2=0
c3=A(1,:)+u*[0,1] %lateral acceleration
d3=B(1)
t=[0:0.01:2]
st=(pi/180)*5*sin(2*pi*7*t)
C=[c1;c2;c3]
D=[d1;d2;d3]
sys=ss(A,B,C,D)
Y=lsim(sys,st,t)
sys.Inputname={'st'}
sys.Outputname={'lateral speed','yaw rate','lateral acceleration'}
[mag1,phase1,omega1] = bode(sys({'yaw rate'},{'st'}))
for i=1:length(mag1)
magnitude1(i)=mag1(1,1,i)
i=i+1
end
figure(1)
loglog(omega1/(2*pi),magnitude1,'linewidth',2)
title('yaw gain frequency response')
xlabel('frequency')
ylabel('yaw gain')
grid on
grid on
figure(2)
subplot(2,2,1)
plot(t,st,'linewidth',2)
title('steering')
xlabel('time(s)')
ylabel('steering angle(rad)')
grid on
subplot(2,2,2)
plot(t,Y(:,1),'linewidth',2)
title('lateral speed')
xlabel('time(s)')
ylabel('lateral speed(m/s)')
grid on
subplot(2,2,3)
plot(t,Y(:,2),'linewidth',2)
title('yawrate')
xlabel('time(s)')
ylabel('yaw rate(rad/sec)')
grid on
subplot(2,2,4)
plot(t,Y(:,3),'linewidth',2)
title('lateral accelaration')
xlabel('time(s)')
ylabel('lateral acceleration(m/s^2)')
grid on
figure(3)
step(sys)
grid on
S=stepinfo(sys)
[wn theta]=damp(sys)
Cf = -61171.9437; % front cornering stiffness (N/rad)
Cr = -50976.6191; %rear
YB=Cf+Cr
YR=(a*Cf-b*Cr)/u
YS=-Cf
NB=a*Cf-b*Cr
NR=(((a^2)*Cf)+((b^2)*Cr))/u
NS=-a*Cf
Wn=sqrt(abs(((YB*NR)+(NB*m*u)-(YR*NB))/(m*u*jz)))
Theta=abs(-((NR/jz)+(YB/(m*u)))*(1/(2*Wn)))
Wd=Wn*sqrt(1-(Theta)^2)
tr=abs((m*u*NS)/((YS*NB)-(YB*NS)))
Correlation at 7 Hz:
Gain (from the yaw rate and steering rate) = 1.03/0.0855 12.04
From frequency response, it can be seen that gain is approx. 12.05
The steady-state gain of 30 is found in the time plot as well as the frequency plot
There is a minor difference between the analytical calculation of rising time and the automated calculation using the step info function in MATLAB.
Results:
YB=-7.0654*10^5
YR=177.87
YS=2.7680*10^5
NB=1.7293810^4
NR=-1.4585*10^4
NS=-4.8175*10^5
Wn=13.3077 rad/s
Theta=0.9837
Wd=2.3934 rad/s
tr=0.1058 s