Discrete time PID control

[Andrew Holme]

I've implemented a DSSS receiver in an FPGA. The code rate NCO is
controlled by an early/late detector. The carrier NCO is controlled by a
Costas Loop. I chose PI gain constants for the control loops by trial and
error but would like to do it more scientifically by analysing stability
using a Bode plot of open loop gain. I always do it this way for analogue
PLLs; but I'm new to z-domain stability analysis.

TimW - I bought your book; but, to be honest, found it quite heavy going.
Nevertheless, adapting some code from the accompanying CD, I came up with
this Scilab script:

kPD = 2^20;
kP = 2^(41-64);
kI = 2^(15-64);
kNCO = 10000;
G = -kPD * (kP + kI/(%z-1)) * kNCO / (%z-1);
G.dt = 1e-3;
scf(0);
clf;
bode(G);

The resulting Bode plot never even passes through unity gain! So is the
system unconditionally stable? Here's a simplified version of the FPGA
code:

`define KP 41
`define KI 15

reg [63:0] phase, rate; // 64 places after the binary point
reg signed [31:0] err; // = Early - Late

always @ (posedge clk) // 10 MHz
if (millisecond) begin
phase <= phase + rate + (err << `KP);
rate <= rate + (err << `KI);
end else
phase <= phase + rate;

It clocks at 10 MHz and applies corrections every millisecond; hence dt=1e-3
and kNCO = 10000.

Is the Bode plot correct?

TIA

 

[Andrew Holme]

I've implemented a DSSS receiver in an FPGA. The code rate NCO is
controlled by an early/late detector. The carrier NCO is controlled by a
Costas Loop. I chose PI gain constants for the control loops by trial and
error but would like to do it more scientifically by analysing stability
using a Bode plot of open loop gain. I always do it this way for analogue
PLLs; but I'm new to z-domain stability analysis.

TimW - I bought your book; but, to be honest, found it quite heavy going.
Nevertheless, adapting some code from the accompanying CD, I came up with
this Scilab script:

kPD = 2^20;
kP = 2^(41-64);
kI = 2^(15-64);
kNCO = 10000;
G = -kPD * (kP + kI/(%z-1)) * kNCO / (%z-1);
G.dt = 1e-3;
scf(0);
clf;
bode(G);

The resulting Bode plot never even passes through unity gain! So is the
system unconditionally stable? Here's a simplified version of the FPGA
code:

`define KP 41
`define KI 15

reg [63:0] phase, rate; // 64 places after the binary point
reg signed [31:0] err; // = Early - Late

always @ (posedge clk) // 10 MHz
if (millisecond) begin
phase <= phase + rate + (err << `KP);
rate <= rate + (err << `KI);
end else
phase <= phase + rate;

It clocks at 10 MHz and applies corrections every millisecond; hence
dt=1e-3
and kNCO = 10000.

Is the Bode plot correct?

TIA

Is 'err' quantized to two values, early or late bang-bang style, or is
it some finer-grain time error signal?

John

err is fine-grained: approx +/- 2^20

Actually, kPD slope should be 2*(2^20) per chip.

 

[Andrew Holme]

I've implemented a DSSS receiver in an FPGA. The code rate NCO is
controlled by an early/late detector. The carrier NCO is controlled by a
Costas Loop. I chose PI gain constants for the control loops by trial
and
error but would like to do it more scientifically by analysing stability
using a Bode plot of open loop gain. I always do it this way for
analogue
PLLs; but I'm new to z-domain stability analysis.

TimW - I bought your book; but, to be honest, found it quite heavy going.
Nevertheless, adapting some code from the accompanying CD, I came up with
this Scilab script:

kPD = 2^20;
kP = 2^(41-64);
kI = 2^(15-64);
kNCO = 10000;
G = -kPD * (kP + kI/(%z-1)) * kNCO / (%z-1);
G.dt = 1e-3;
scf(0);
clf;
bode(G);

The resulting Bode plot never even passes through unity gain! So is the
system unconditionally stable? Here's a simplified version of the FPGA
code:

There's several things wrong with this -- does the loop actually work?

Your gain is HUGE -- I'm seeing nearly 60dB of gain at Fs/2, which should
guarantee wild Fs/2 oscillations.

kP / kI is 67 million some odd, which means your integrator doesn't start
doing anything until Fs/67000000 or so -- that's a pretty low frequency
for an integrator to kick in.

As a minor point, you're including the subtraction in the phase
comparator, meaning that your phase shift criterion is 0 rather than
180 -- that's not how most folks do it, but there's a significant enough
minority that do that I can't quibble _too_ much.

As another minor point, Scilab (and anything that uses polynomials for
this sort of analysis) is vulnerable to numerical errors when the order of
the system gets high. Using

-->G = kPD * tf2ss(kP + kI / (%z - 1)) * tf2ss(kNCO / (%z - 1));

dodges this (note that I'm using my preferred sign convention). It's not
the issue here.

`define KP 41
`define KI 15

reg [63:0] phase, rate; // 64 places after the binary point
reg signed [31:0] err; // = Early - Late

always @ (posedge clk) // 10 MHz
if (millisecond) begin
phase<= phase + rate + (err<< `KP);
rate<= rate + (err<< `KI);
end else
phase<= phase + rate;

It clocks at 10 MHz and applies corrections every millisecond; hence
dt=1e-3
and kNCO = 10000.

Is the Bode plot correct?

I'm not going to pick through that Verilog code for a complete model
unless you pay me, but I know it's not matching the transfer function that
you give. Here's some points that I can pick out:

* If KP were 0, then phase <= phase + rate + (err << 'KP) would
give an effective kP of 1, or if you're viewing your variables
as fractional, then kP of 2^-32. Not 2^-64.

* You are only applying your proportional term, and updating your
integrator, on the millisecond intervals. But you are applying
your integrator state to the phase every NCO sampling interval.
This gives your integrator an extra gain of 10MHz * 1ms.

It would arguably be much better style (and it would certainly
give you smoother operation) if you calculated a frequency
increment command separately from your controller update. You
may not want the slowness that comes with smoothness, and you
may not want to carry the extra storage -- which is where the
arguing comes in.

So: kP is off by 2^-32, kI is off by 10000 * 2^-32, and the gains in your
model are already way too high. If things work, this makes me suspect
that your kPD _must_ be wrong somehow (maybe it's really 2^(20-64)?)

If I make the factor of 10000 correction to kI I get an integrator zero at
around 0.04Hz, which is not terribly unreasonable, if a bit slow. I'm
still not getting a gain of 1 in the region of 90 degrees phase shift, but
I'll leave that part to you.

--

Tim Wescott
Wescott Design Services
http://www.wescottdesign.com

Do you need to implement control loops in software?
"Applied Control Theory for Embedded Systems" was written for you.
See details at http://www.wescottdesign.com/actfes/actfes.html

Tim, thank you. I had applied that kNCO=10000 multiplier to both kI and kP.
It should only be applied to kI for exactly the reason you have pointed out.

The following code is giving me a sensible looking Bode plot:

kPD = 2^20;
kP = 2^(41-64);
kI = 10000 * 2^(18-64);
G = -kPD * (kP + kI/(%z-1)) / (%z-1);
G.dt = 1e-3;
scf(0);
clf;
bode(G);

I'm using `KI=18 in the latest FPGA. The rate was taking several seconds to
settle.

Thanks for the free consultancy I will study the book.