Numerical inacurracy when calculating closed-loop system's zeros

Hi there! I would like to start by saying that I admire the great work you're doing! I'm new to Python and open source projects in general, and I've been trying to move from Matlab to Python in my control/signal processing projects. I've been looking for a way to work with internal time delays in control systems, and found out about this amazing project!

I'm having a problem when trying to calculate the zeros from a closed-loop transfer function. The zeros that I'm obtaining have significant numerical inaccuracy when compared to the PID-controller's zeros (they should be the same in this case). The following sample code generates the problem that I'm talking about. Also, there is another (minor) problem: when asking the zeros of a MIMO system, it returns the system's poles and zeros (at the end of the code, there is a part where this error occurred for me):


import numpy as np
import robustcontrol as rc

######################### second order system plus dead time (SOPDT) ##########################
K = 0.9
w_n = 2
csi = 0.3
theta = 0.1

# System's transfer function
G = K*rc.utils.tf([w_n**2], [1, 2*csi*w_n, w_n**2], deadtime = theta)

##################################### Ziegler-Nichols PID Tuning ################################

# critical gain Kcr 
K_cr = rc.utils.margins(G)[0][0] 

# critical oscillation period 
P_cr = np.abs(2*np.pi/rc.utils.margins(G)[3][0])

# PID gains 
Kp = 0.6*K_cr
Ki = Kp/(0.5*P_cr)*rc.utils.tf([1], [1, 0])
Kd = Kp*0.125*P_cr* rc.utils.tf([1, 0],[1])

# PID controller transfer function
C_ZN = Kp + Ki + Kd

################################# Calculating system zeros ###########################################

# Delay free model (zero-order padé approx)
G_delay_free = K*rc.utils.tf([w_n**2], [1, 2*csi*w_n, w_n**2])

# closed-loop considering a delay-free model
H = G_delay_free*C_ZN/(1 + G_delay_free*C_ZN)

# Zeros should be equal in both cases, but there is a significant numerical error (the same does not occur in matlab or using python-control package)
H_zeros = H.zeros()
# returns  array([-2.49594595+0.0365365j, -2.49594595-0.0365365j])  

C_ZN_zeros = C_ZN.zeros()
# returns array([-2.49652882+1.43703181e-05j, -2.49652882-1.43703181e-05j])

# When asking the zeros of a mimo transfer function, it return poles and zeros
H_mimo_style = rc.utils.mimotf(H)
H_mimo_style.zeros()

Thanks for your attention and sorry for the bad coding or bad English, I've been trying to exercise both of those things recently. I'm still learning how to use GitHub too! :)

alchemyst / Skogestad-Python

Numerical inacurracy when calculating closed-loop system's zeros #387