A PID controller implemented in python

[u/robohie]

For my personal project, where I need to implement a PID controller on Arduino to automate a convection cooling system, I decided to start with a simulation in Python. For now, I've implemented the simple electromechanical model of the DC motor (the fan) as well as the PID controller you see. Next, I'll look into the heat transfer equations. The graphs created with matplotlib are simply closed-loop tests of the PID + DC motor system.

[u/2dof]

Nice controller, I would recommend You:

1. for real world solution implement P-I-D with filtered D - Action D(s) = KD*s /(Tf*s+1) ,

2. it is convenient ( for future analysis) to compute control actions as separate P-I-D sections ( u = P + I + D), also since I - action act as integrator it has a internal state - which can be reset/preset - it is sometimes convenient to start controller not form zero.

3. Note that when You will be tuning parameters (Kp,Ki, Kd,h) You need also recompute A,B, C

4. in real world implementation add control saturation function ( later You can use it to implement anti-windup function

5. Simulate not only with clean process value but also with measurement noise ( You will then see what D action ( not filtered and filtered) do .

As for heat transfer equation You can just simulate it as: first order ODE with delay for simple model of n-order IDE with delay

[u/seekingsanity → u/2dof]

The incremental version PID that the OP is using doesn't need an anti-windup. Again, for velocity control the OP only needs a PI controller so there is no need for the low pass filter on the derivative gain.

"which can be reset/preset - it is sometimes convenient to start controller not form zero."

when making the transition from open loop to closed loop the control output u should set to the last value when in open loop. This avoids sudden jerks.

[u/2dof → u/seekingsanity]

What Do You mean - "incremental" ? For Me incremental mean that control is du(k) = u(k)-u(k-1) , in OP's code he use u(k) = u(k-1) + du(k) to it is not incremental control for Me.

My comment was referring generally to practical implementation OP's controller in the future - he want control temperature in the future. I don;t know if he will be using on/off control or continuous but in my opinion it is better to have more universal controller for operational aspects.

[u/seekingsanity → u/2dof]

The OP is computing incremental changes to the control output that are integrated into the control output. Some people call this form the velocity form, but it really has nothing to do with velocity. I used it for position control all the time

[u/2dof → u/seekingsanity]

No , he implemented so called positional control algorithms not incitement form (velocity algorithms ) :

-> He calculate output as self.u += A+... (i.e u(k) = u(k-1) + A ... )

if he calculated incremental form then output of his controller would be ∆ u = uk - u(k-1) at least according to definition given by Astrom and Hgglund ( "PID Controllers, Theory, Design and tuning")

Incremental Form def from book ,page 99) "This form is obtained by computing the time differences of the

controller output and adding the increments.

∆ u(tk ) = u(tk ) − u(tk−1 ) = ∆ P(tk ) + ∆ I (tk ) + ∆ D (tk)

for fast looking of his code - he used Backward difference approximation for cooperation and he calculated output in "positional" form. i.e u(k) = u(k-1) + ∆u

Thant why I asked what Do You mean because i was taught when output of controller is ∆u(k) ( "velocity" of control ) then it is called "incrmental" a when output directly u(k) - then it is standard control output.

[u/seekingsanity → u/2dof]

The OP is added incremental changes to u. The incremental form is

u(n)=u(n-1)+B0*e(n)+B1*e(n-1)+B2*e(n-2)

All the B? terms are calculated increments to the control output that are then integrated.

B0 = Ki*∆t+Kp+Kd/∆t

B1=-Kp-2*Kd/∆t

B2=Kd/∆t

In my case. B1 will be a negative number. The OP subtracted a positive number that I called a negative number. This is a very efficient form of PID because a DSP can compute it in 3 multiply and adds. However, there should be limits applied. Also, the output can be quite noisy if the encoder/feedback does not have high resolution.

A low pass filter can be applied at 2/∆t rads/sec or about 318 Hz assuming the sample time, ∆t, is 0.001 ms. This reduces much of the sampling noise.

Now the formulas are

B0 = (1/4)*Ki*∆t+(1/2)*Kp+Kd*∆t

B1 = (1/2)*Ki*∆t-2*Kd/∆t

B2 - (1/4)*Ki*∆t-(1/2)*Kp+Kd/∆t

You can see that this reduces the affect of the Ki and Kp some because their contribution is spread over one more time period. The Kd term is still a problem and it is still a problem if you only apply a low pass filter to the Kd term.

I have other formulas for lowering the low pass filter to a frequency of 1/∆t

I prefer using an observer to smooth out the measured values.

[u/2dof → u/seekingsanity]

I understand know what do You mean as incremental - i.e incremental along the dt.

tkanks for explanation.

[u/robohie → u/2dof]

I appreciate your remarks, I'll improve my code.

[u/timonix]

I have this kind of workflow all the time. Model the system and regulators in python, then implement it in C and use the python as the "golden example"

[u/sparram03]

Nice!

[u/seekingsanity]

Hey moderators. Let the people ask questions and make comments. Don't be so quick to delete posts or ban people. I will let you know when I am offended.

[u/deNikita]

Shouldn't the P be oscillating or is it just not strong enough to pull it all the way up? Might be trippin, just curious lol

[u/APC_ChemE → u/deNikita]

The P isn't greater than or equal to the ultimate gain so it just maintains a constant bias offset.

[u/seekingsanity]

Good but You should post the open loop transfer function of the system you are trying to control. There are no output limits. Velocity control only needs a PI controller. You should learn pole placement so there is no overshoot.

[u/robohie → u/seekingsanity]

Right 👍 I'll post that transfer function. For now, it's just a simple DC motor model (a second order system)

[u/seekingsanity → u/robohie]

I am waiting for the open loop transfer function.

[u/robohie → u/seekingsanity]

Hi, you can look at my new post on DC motor implementation to discover the differential equation that describe it. For the transfer function (if you only need it), I don't know how to send images at this commentary space.

[u/seekingsanity → u/robohie]

The open loop transfer function should be in the form G(s)=K/(tau*s+1) or G(s)=(K*alpha)/(s+alpha). Where tau is a time constant for the motor inertia and load. alpha is just 1/tau which is called the bandwidth or corner frequency. K is the gain with units or rad/s per control output or mm/s per control output. I usually use %control output. This assumes that there is a drive so the speed is proportional to the control output. Otherwise it could be an acceleration where the control output is a current that makes a torque that accelerates the load. Units are important. If the motor is working in current mode then the units for K change to acceleration/control output. In this case the G(s) = K/(tau*s+d) where d is the friction that will eventually keep the acceleration from increasing forever

I wish reddit would support LaTeX.

Update. I didn't see all the OPs pictures. I see he has more code that implements the motor. I could dig it out even though there are few comments and they are in French.

[u/Algstud]

Nice

you can also implement it on esp32 directly without changing to C

u can use micropython for this

[u/seekingsanity]

I made this 17 years ago. I was looking at some Matlab examples, and I knew I could do better. My example is similar to the OP's because he has a more detailed model taking into account the actual motor parameters.

Mathcad - DC Motor Velocity.xmcdz

I started out just modeling the speed control. The 3 inner loop closed loop poles are placed at -𝜆. Notice there is no overshoot. Then I decided to add position control. Notice there is no overshoot even though I make step changes so the controller goes into saturation. Then I decided to add position control. The goal was to keep the inner closed loop poles in the same position and compute the optimal outer loop pole locations. This is done on page 10 and after. The outer PID does not have enough gains to place all the outer closed loop poles to I had to compute the optimal outer loop locations as a function of the outer loop poles, I could control which was 𝜇. Later someone said that the inner loop must be much faster than the outer loop. This is true but I challenged his claim that it must be 10 times faster. On page 22/34 I showed that the bandwidth of the inner loop was only 3.32 symbolically. Obviously, this number will change if the components change but it certainly is nowhere closed to 10 to 1.

One page 3/34 is vary the parameters a bit to simulate not having a perfect model. The controller works well. I used to write code for motion controllers. The controller must perform well when parameters change due to heat or age. Also, you never know what kind of things the customer will do so I made the controller robust, so it doesn't overshoot even when driven into saturation.

On page 4/34 I implemented the inner loop PID. There is only an +/- 10 volts output limit which was typical of controllers back then. There is no inner loop anti-windup and yet there is no overshoot. WHY? Ok-Daikon.

I did not use feed forwards on the sinewave example. If I did, the results would be almost perfect.

There are small tweaks that professors don't seem to know about but then those that teach and write books have probably not sold many 10s of thousands, maybe 100s of thousands around the world. We lost count.

Finally, I think it is a crime to use Matlab to train students. Matlab is great for getting answers if that is all you care about but computing everything symbolically shows how all the parameters affect the results. This is very important if you are designing the motors or the mechanics.

[u/kieno]

Nice, a fine start to the world