# Discretization of a PID controller

I've been given a plant transfer function to control trough a PID. The transfer function is:

$$G_{\small{p}}(s) = \frac{1}{(1+s)(1+2s)}$$

Doing a little work with matlab I got this PID controller transfer function:

$$G_{\small{c}}(s) = 3 + \frac{1}{s} + \frac{2s}{1+\frac{s}{60}}$$

So $$K_{\small{p}} = 3, K_{\small{i}} = 1, K_{\small{d}} = 2$$. Please note that I needed to put an high frequency pole into derivative term to keep it real. That pole has been put reasonably away from the bandwidth limit of the plant (cut-off frequency about $$2 rad/s$$).

Checking the step response on matlab everything seems to work really fine.

Now I need to discretize such controller function and provide anti-windup functionality, I'm not an expert but since I need an antiwindup I think I need to divide the integral action from the proportional and derivative action.

Using Tustin method and a sampling time of $$10^{-3}s$$ I got these:

uP[k] = 3*e[k]
uI[k] = 0.0005*(e[k] + e[k-1]) + uI[k-1]
uD[k] = 116.5*(e[k] - e[k-1]) + 0.9417 * uD[k-1]


of course the resulting command will be u[k] = uP[k] + uI[k] + uD[k]

The C code become:

double control(const double y, const double r)
{
static double e[2] = { 0, 0 };  // e[0] -> e_k-1, e[1] -> e_k
static double uD[2] = { 0, 0 }; // uD[0] -> uD_k-1, uD[1] -> uD_k
static double uI[2] = { 0, 0 }; // uI[0] -> uI_k-1, uI[1] -> uI_k
double uP, u;

e[1] = r - y;

uP = Kp * e[1];
uI[1] = Ki * (e[1] + e[0]) + uI[0];
uD[1] = Kd1 * (e[1] - e[0]) + Kd2 * uD[0];

u = uP + uI[1] + uD[1];

// Cut-Off
if(u>MAX_CMD){
u = MAX_CMD;
}
else if(u<MIN_CMD) {
u = MIN_CMD;
}
//

e[0] = e[1];
uD[0] = uD[1];
uI[0] = uI[1];

return u;
}


I think each step is correct but I don't feel very comfortable with such high derivative terms and low integral. Moreover I cannot figure out an effective antiwindup feature (I always used it on PI controllers never on PID). Something like that should work?

    // Cut-Off
if(u>MAX_CMD){
u = MAX_CMD;
uI[1] = 0; // Antiwindup
}
else if(u<MIN_CMD) {
u = MIN_CMD;
uI[1] = 0; // Antiwindup
}
//

• 1 Why did you make uI[1]=0 (effectively uI[0] is also made zero immediately afterwards). Do you have a text book reference for this type of anti wind-up implementation? 2 "I don't feel very comfortable with such high derivative terms and low integral" Is that part of the question? It seems to be un related to the anti-wind-up part of the question. 3 why is a condition on u checked for anti-wind-up rather than uI it self?
– AJN
Feb 13 at 15:42
• Re: the high D coefficient - I didn't follow your calculation, but yes I would question it. Likely to amplify noise unnecessarily. Feb 13 at 16:17
• @PeteW I've not been asked to take noise into account. Feb 13 at 18:05
• @PeteW it looks like the OP just asked Matlab, and it coughed up answers where the $1 / \Delta t$ and $\Delta t$ terms were implicit -- if you account for that, then the gains should look a lot more like the time-domain gains. Feb 14 at 4:24
• @weirdgyn - See here: dsp.stackexchange.com/questions/58533/… ... derivative, and even the identity function $y(t)=u(t)$ are both causal and anti-causal. In continuous-time, it's not an issue. In my view it only becomes an issue when you do discretize, since something like $\frac{u[n]-u[n-1]}{T_s}$ most closely approximates $du/dt$ at a point in time a half-sample ago, rather than "now". If I do the analysis in continuous time, and my discrete implementation has sufficient sample rate, I can worry about something else. Feb 29 at 19:42

Re: anti-windup

Almost there. The essential part of the logic that's worked for me:

• If the output is saturated to the maximum, then do not increase the integrator variable (but let it fall if it wants)

• If the output is saturated to the minimum, then do not decrease the integrator variable (but let it rise if it wants)

It can be refined further, but that's 90% of it

• Interesting .. I'll try this! Feb 14 at 9:15

Your gains look perfectly fine. The Tustin approximation* works by noting that $$s \simeq \frac 2 T \frac{z - 1}{z + 1} \tag 1$$ then substituting the expression in $$z$$ for $$s$$ everywhere that it occurs. So your integrator gain gets multiplied by about $$T / 2$$, and your derivative gain gets multiplied by about $$2 / T$$ -- except your derivative gain is modified by the pole at $$z = 0.9417$$, so that nice tidy multiplication isn't quite so obvious.

As to integrator anti-windup, there is no best way to implement it. It's a nonlinear solution to a nonlinear problem**. There are a few popular ways which I suggest you research. It looks like you're limiting the integrator state to the minimum or maximum actuator drive -- this is the easiest and most basic form of integrator anti-windup.

* Not method. It is the Tustin approximation. You are approximating a transfer function in $$s$$ with one in $$z$$, and it's best to never forget that, especially when your sampling rate gets within a factor of ten of your fastest poles or zeros.

** Even if your plant described by a transfer function and is thus presumed to be linear, the use of integrator anti-windup implies the existence of some limit on the actuator strength -- since, from the point of view of your controller, the actuator is part of the plant, so have to take the nonlinearity into account.

• The sampling time has been selected as $10^{-3}s$ accordingly with the fastest pole (60 rad/s) . Feb 14 at 8:41
• I know all the discretization theory Tustin, Euler and other substitution methos... but the typical problems I work on always lead to PI controller and never full PID .. this time have been a little bit suprised by the differences between gain factor of the two parts. Feb 14 at 8:44
• in my control theory books they refer to Tustin method (or bilinear transform) ... of course Tustin method use the first order Padè approximant in the form you reported but they are called methods o transforms. Feb 29 at 8:11