Problem when utilizing the bode function to plot the multiplication of several transfer functions

Hi all,

I'm right now trying to plot the bode plot of a system that contains multiple transfer fucntions, for example: G1*G2*G3.

(all the tranfer functions are declared by using tf() )

However, it just occurred to me that the bode plot of G1*G2*G3 is different from G2*G1*G3. And I felt quite confused about the results since mathematically they should be the same?

==> This is the results of the code below, as you can see, the bode plots are different.

The example code is just similar to:

%% Initialization
omega_D_Val = 50;
M_Val = 4;
Kf_Val = 40;
Mn_Val = 4;
Kfn_Val = 40;
Mn_head_Val = 4/1.5;
Kfn_head_Val = 40;
Ts_Val = 0.0005;
wR_Larger = 3000;
Cf_Val = 1;
Pure_gain = (M_Val/70) * wR_Larger * (1/M_Val);

D_env_Val = 1.5;
K_env_Val = 2500;

Integrator = tf(Ts_Val,[1 -1],Ts_Val);

%% Plotting Part
figure;
set(gcf,'Units','centimeters','position',[0 0 25 30]);
[~, Tenv_None, Nenv_None, ~, OSenv_None]...
= OuterLoop_TransferFcn_withEnv(omega_D_Val, wR_Larger, Ts_Val, (M_Val/70), Mn_Val, M_Val, Kfn_head_Val, Kfn_Val, Kf_Val, D_env_Val, K_env_Val, Cf_Val);

[Gain, ~, Frequency] = bode(Pure_gain * Integrator * Tenv_None * OSenv_None, {1e-4,6283.19});
semilogx(reshape(Frequency,[],1), 20*log10(reshape(Gain,[],1)), "-", 'Color', 'k', 'LineWidth',1.5);
hold on

[Gain, ~, Frequency] = bode(Integrator * Tenv_None * OSenv_None * Pure_gain, {1e-4,6283.19});
semilogx(reshape(Frequency,[],1), 20*log10(reshape(Gain,[],1)), "--", 'Color', 'k', 'LineWidth',1.5);

%% Function for tranfer functions build up
function [G_env, T_env, N_env, Compensator_env, OuterLoop_Sensitivity_Function_env]...
= OuterLoop_TransferFcn_withEnv(omega_D_Val, omega_R_Val, Ts_Val, Mn_head_Val, Mn_Val, M_Val, Kfn_head_Val, Kfn_Val, Kf_Val, D_env_Val, K_env_Val, Cf_Val)
    syms omega_R omega_D Ts Mn_head Mn M Kfn_head Kfn Kf D_env K_env
    syms z
    
    Cf = Cf_Val;
    
    % Alpha and Beta Define
    alpha = (Mn*Kf)/(M*Kfn);
    beta = (Mn*Kfn_head)/(Mn_head*Kfn);

    Assigned_Parameters = [omega_D, omega_R, Ts, Mn_head, Mn, M, Kfn_head, Kfn, Kf, D_env, K_env];
    Assigned_Values = [omega_D_Val, omega_R_Val, Ts_Val, Mn_head_Val, Mn_Val, M_Val, Kfn_head_Val, Kfn_Val, Kf_Val, D_env_Val, K_env_Val];
    
    Denominator = (2*(z - 1)*((2 + omega_R*Ts)*z + (-2 + omega_R*Ts))*(M*z^2 + (M*alpha*omega_D*Ts + D_env*Ts - 2*M)*z + (M + K_env*Ts^2 - D_env*Ts - M*alpha*omega_D*Ts)))...
                + (Mn_head*Cf*omega_R*Ts*((2 + omega_D*Ts)*z + (-2 + omega_D*Ts))*(M*beta*z^3 + (D_env*Ts*beta - 2*M*alpha - M*beta)*z^2 + (K_env*beta*Ts^2 + 4*M*alpha -M*beta)*z + (K_env*beta*Ts^2 - D_env*beta*Ts - 2*M*alpha + M*beta)));
    Denominator = subs(Denominator, Assigned_Parameters, Assigned_Values);
    
    G_tf_Numerator = (Mn_head*omega_R*Ts*((2 + omega_D*Ts)*z + (-2 + omega_D*Ts))*(M*beta*z^3 + (D_env*Ts*beta - 2*M*alpha - M*beta)*z^2 + (K_env*beta*Ts^2 + 4*M*alpha -M*beta)*z + (K_env*beta*Ts^2 - D_env*beta*Ts - 2*M*alpha + M*beta)));
    G_tf_Denominator = (2*(z - 1)*((2 + omega_R*Ts)*z + (-2 + omega_R*Ts))*(M*z^2 + (M*alpha*omega_D*Ts + D_env*Ts - 2*M)*z + (M + K_env*Ts^2 - D_env*Ts - M*alpha*omega_D*Ts)));
    G_tf_Numerator = subs(G_tf_Numerator, Assigned_Parameters, Assigned_Values);
    G_tf_Denominator = subs(G_tf_Denominator, Assigned_Parameters, Assigned_Values);
    [G_tf_Numerator, G_tf_Denominator] = numden(G_tf_Numerator/G_tf_Denominator);
    [bG,~] = coeffs(G_tf_Numerator, z, 'All');
    [bG_D,~] = coeffs(G_tf_Denominator, z, 'All');
    
    T_env_Numerator = M*z^2 + (alpha*omega_D*Ts*M - 2*M)*z + (M - alpha*omega_D*Ts*M);
    T_env_Denominator = M*z^2 + (D_env*Ts - 2*M + alpha*omega_D*Ts*M)*z + (K_env*Ts^2 - D_env*Ts + M - alpha*omega_D*Ts*M);
    T_env_Numerator = subs(T_env_Numerator, Assigned_Parameters, Assigned_Values);
    T_env_Denominator = subs(T_env_Denominator, Assigned_Parameters, Assigned_Values);
    [T_env_Numerator, T_env_Denominator] = numden(T_env_Numerator/T_env_Denominator);
    [bT,~] = coeffs(T_env_Numerator, z, 'All');
    [bT_D,~] = coeffs(T_env_Denominator, z, 'All');
    
    N_env_Numerator = M*z^2 + (D_env*Ts - 2*M)*z + (M + K_env*Ts^2 - D_env*Ts);
    N_env_Denominator = M*z^2 + (D_env*Ts - 2*M + M*alpha*omega_D*Ts)*z + (M + K_env*Ts^2 - D_env*Ts - M*alpha*omega_D*Ts);
    N_env_Numerator = subs(N_env_Numerator, Assigned_Parameters, Assigned_Values);
    N_env_Denominator = subs(N_env_Denominator, Assigned_Parameters, Assigned_Values);
    [N_env_Numerator, N_env_Denominator] = numden(N_env_Numerator/N_env_Denominator);
    [bN,~] = coeffs(N_env_Numerator, z, 'All');
    [bN_D,~] = coeffs(N_env_Denominator, z, 'All');
    
    C_tf_Numerator = (Mn_head*Cf*omega_R*Ts*((2 + omega_D*Ts)*z + (-2 + omega_D*Ts))*(M*beta*z^3 + (D_env*Ts*beta - 2*M*alpha - M*beta)*z^2 + (K_env*beta*Ts^2 + 4*M*alpha -M*beta)*z + (K_env*beta*Ts^2 - D_env*beta*Ts - 2*M*alpha + M*beta)));
    C_tf_Numerator = subs(C_tf_Numerator, Assigned_Parameters, Assigned_Values);
    [C_tf_Numerator, C_tf_Denominator] = numden(C_tf_Numerator/Denominator);
    [bC,~] = coeffs(C_tf_Numerator, z, 'All');
    [bC_D,~] = coeffs(C_tf_Denominator, z, 'All');
    
    OS_tf_Numerator = (2*(z - 1)*((2 + omega_R*Ts)*z + (-2 + omega_R*Ts))*(M*z^2 + (M*alpha*omega_D*Ts + D_env*Ts - 2*M)*z + (M + K_env*Ts^2 - D_env*Ts - M*alpha*omega_D*Ts)));
    OS_tf_Numerator = subs(OS_tf_Numerator, Assigned_Parameters, Assigned_Values);
    [OS_tf_Numerator, OS_tf_Denominator] = numden(OS_tf_Numerator/Denominator);
    [bOS,~] = coeffs(OS_tf_Numerator, z, 'All');
    [bOS_D,~] = coeffs(OS_tf_Denominator, z, 'All');
    
    % Transfer Function Declarations
    G_env = tf(double(bG)/double(bG_D(1)), double(bG_D)/double(bG_D(1)), Ts_Val);
    T_env = tf(double(bT)/double(bT_D(1)), double(bT_D)/double(bT_D(1)), Ts_Val);
    N_env = tf(double(bN)/double(bN_D(1)), double(bN_D)/double(bN_D(1)), Ts_Val);

    Compensator_env = tf(double(bC)/double(bC_D(1)), double(bC_D)/double(bC_D(1)), Ts_Val); % (G*Cf)/(1 + G*Cf)
    OuterLoop_Sensitivity_Function_env = tf(double(bOS)/double(bOS_D(1)), double(bOS_D)/double(bOS_D(1)), Ts_Val); % 1/(1 + G*Cf)
    
end

Answer 1

[Tsai Han Hao]

Hi all,

I've found out that the tf data type is the main problem that causes the difference between the Bode plot of bode(G1*G2*G3) and bode(G2*G3*G1).

If you look through the Algorithms on how Matlab do the bode calculations, you can find out that zpk is the data type that Matlad do the Bode plot.

https://www.mathworks.com/help/ident/ref/lti.bode.html

Therefore, by utilizing bode(zpk(G1)*zpk(G2)*zpk(G3)), the results of the Bode plot will be correct no matter how you arrange the tranfer functions.