Google Groups no longer supports new Usenet posts or subscriptions. Historical content remains viewable.
Dismiss

PID autotuning - not working for heating application

90 views
Skip to first unread message

Frank W.

unread,
Nov 5, 2009, 12:37:21 PM11/5/09
to
Is there a control engineering expert here?

I could us a bit of help on how to implement a PID autotune function for a
heating application (a small boiler).

My current PID autotune function produces no reliable results. I use the
relay feedback method (�strom and H�gglund) but it seems that the ultimate
period Tu - one of the paramters determined with the relay feedback test -
is directly correlated to the relay output step u. u is an arbitrary value
which makes Tu an arbitary value. Since Tu is required to compute Ti and Td
(e.g. Ti = 0.5 Tu = Ziegler-Nichols), autotune is not possible.

I believe this is because the machine heats very quickly (1500W boiler) and
cools down very slowly, so the process value isn't a sinusoid. Here is an
example graphics, green line is setpoint, red plot is process value and the
grey vertical bars represent heat output u:

http://img196.imageshack.us/img196/7409/examplejr.jpg

Regardless of u, the temperature always dips the same, small amount below
setpoint - because the machine cools slowly, the temperature can not fall
far below setpoint before it's reacting to the heat. It then shoots up by an
amount that is proportional to u. As a result, the plot resembles mountains
above setpoint. If u is large, temperature will shoot high above the
setpoint and take very long to cool down. Big mountains, so Tu gets large,
Tu ~ u = my problem.

Since all PID temperature controllers have Autotune, there must be a
solution for this problem. Any ideas?


Vladimir Vassilevsky

unread,
Nov 5, 2009, 2:09:43 PM11/5/09
to

Frank W. wrote:

> Is there a control engineering expert here?

You can contact Tim Wescott who is the expert in control systems,
however the problem is not really a rocket science.

> I could us a bit of help on how to implement a PID autotune function for a
> heating application (a small boiler).

It looks like a huge lag in the control loop. PID approach won't work
very well. Is your heater on-off type or proportionally regulated?

> Since all PID temperature controllers have Autotune, there must be a
> solution for this problem. Any ideas?

To begin with, abandon autotune and see if you can get good control by
manual adjustment of parameters.


Vladimir Vassilevsky
DSP and Mixed Signal Design Consultant
http://www.abvolt.com

Bob

unread,
Nov 5, 2009, 5:43:57 PM11/5/09
to
As someone else said: try manually adjusting parameters, instead. Do you
really need autotune? If the boiler/heater/load/setpoint don't change
much, then twiddle the parameters until you like the result.

That said, there are at least a dozen algorithms for autotuning. They
fall into two basic methods (relay and step), but the parameters they
generate can be quite different due to the different desired responses
(e.g. "fastest to reach a certain error after a step change in setpoint"
vs "little or no overshoot after a disturbance", etc.). You have to be
clear about what you're trying to do before you bother trying to tune it.

Large delay and overshoot are common to heater control loops. Placing
the sensor as close as possible to the heater can help, but then you
might not be measuring the thing you want hot. I just had a (maybe)
similar problem that included a requirement for a smooth ramp of the
temperature. That led me to close the loop around the slope (calculated
with an extremely simple FIR filter) during ramping. The slope-loop
worked far better than trying get the temperature-loop to follow a
ramping setpoint.

YMMV
Bob

Tim Wescott

unread,
Nov 6, 2009, 12:48:04 AM11/6/09
to
On Thu, 05 Nov 2009 13:09:43 -0600, Vladimir Vassilevsky wrote:

> Frank W. wrote:
>
>> Is there a control engineering expert here?
>
> You can contact Tim Wescott who is the expert in control systems,
> however the problem is not really a rocket science.

Or you can wait for him to get around to doing the newsgroup thang.

(thanks for the plug).

I'll bet there's rockets out there with temperature loops.

>> I could us a bit of help on how to implement a PID autotune function
>> for a heating application (a small boiler).
>
> It looks like a huge lag in the control loop. PID approach won't work
> very well. Is your heater on-off type or proportionally regulated?
>
>> Since all PID temperature controllers have Autotune, there must be a
>> solution for this problem. Any ideas?
>
> To begin with, abandon autotune and see if you can get good control by
> manual adjustment of parameters.
>

Ditto -- try it by hand, and go from there.

--
www.wescottdesign.com

Tim Wescott

unread,
Nov 6, 2009, 12:55:05 AM11/6/09
to
On Thu, 05 Nov 2009 18:37:21 +0100, Frank W. wrote:

> Is there a control engineering expert here?
>
> I could us a bit of help on how to implement a PID autotune function for
> a heating application (a small boiler).
>
> My current PID autotune function produces no reliable results. I use the

> relay feedback method (Åstrom and Hägglund) but it seems that the


> ultimate period Tu - one of the paramters determined with the relay
> feedback test - is directly correlated to the relay output step u. u is
> an arbitrary value which makes Tu an arbitary value. Since Tu is
> required to compute Ti and Td (e.g. Ti = 0.5 Tu = Ziegler-Nichols),
> autotune is not possible.

Are you sure you haven't just stopped too soon in turning the output step
down? I see where you're having trouble with this; is there some reason
you can't turn u down even further? Things may start looking better when
you get it down to where the temperature curve is roughly symmetrical
around the average.

> I believe this is because the machine heats very quickly (1500W boiler)
> and cools down very slowly, so the process value isn't a sinusoid. Here
> is an example graphics, green line is setpoint, red plot is process
> value and the grey vertical bars represent heat output u:
>
> http://img196.imageshack.us/img196/7409/examplejr.jpg
>
> Regardless of u, the temperature always dips the same, small amount
> below setpoint - because the machine cools slowly, the temperature can
> not fall far below setpoint before it's reacting to the heat. It then
> shoots up by an amount that is proportional to u. As a result, the plot
> resembles mountains above setpoint. If u is large, temperature will
> shoot high above the setpoint and take very long to cool down. Big
> mountains, so Tu gets large, Tu ~ u = my problem.
>
> Since all PID temperature controllers have Autotune, there must be a
> solution for this problem. Any ideas?

Why do you want to use autotune? If this is a product that you're
working on, and if the boiler design is the same for all of the parts
that you ship, then you should just need one tuning.

Autotune is for when you're selling a shrink-wrapped controller that has
to work for anything -- and autotune is often considered to be a good way
to get the tuning in the ballpark so that a human can get involved and
actually make it work right.

Unless you need to adapt to a variety of different, unexpected boiler
combinations -- that don't change as the boiler is operating -- then
there's not much value in autotune, IMHO.

--
www.wescottdesign.com

Frank W.

unread,
Nov 6, 2009, 9:15:08 AM11/6/09
to
"Tim Wescott" <t...@seemywebsite.com> wrote:

> you can't turn u down even further?

A bit. Sensor resolution is the limit. The temperature hardly drops below
the setpoint during relay test iterations because the system cools so slowly
(it drops 0.3 degree below setpoint in case of my machine's 200ml boiler),
so making it nearly symmetric creates a PV signal with very few steps.

> Things may start looking better when

I didn't say that the results are bad. It's just that I don't know which of
the results to pick since every variation of the relay output step produces
a different ultimate period. The lower u is made, the lower Tu gets (there
is no lower limit).

> Why do you want to use autotune? If this is a product that you're
> working on, and if the boiler design is the same for all of the parts
> that you ship, then you should just need one tuning.

The algorithm is for a semi-commercial PID controller for coffee machines,
so autotune can not be left out since users expect it. Boiler sizes rougly
vary from 200ml to a few thousand ml. Wattage varies not so much (from
1250W-2500W).


Tim Wescott

unread,
Nov 6, 2009, 8:17:15 PM11/6/09
to

It sounds like an ideal application for autotune -- you have a wide and
uncontrolled product range, so you can't just make a 'one size fits all'
tuning. On the other hand your plants are constrained to a certain class
of systems, so your basic problem isn't so wide open as to be impossible.

Any sort of tuning exercise involves making some assumptions about the
system you're controlling, poking it with some excitation, seeing how it
responds, then formulating a controller for the thing. In some cases you
have a step in there where you create a mathematical model of your plant,
in others you don't. The process is the same with autotuning, except
that you hope to have a machine doing the tedious bits for you, instead
of you having to run around the world tuning up coffee makers.

The Aostroem-Haggluend method never develops that explicit model, nor
does it work with just any old excitation -- it just pokes the system in
a specific way, extracts some parameters that aren't really closely
related to what control theorists are used to (but which are practically
quite useful), then gets PID tuning parameters from those measured system
parameters.

Your problem is that you're not really poking the system in quite the
right way, so the Aostroem-Haggluend method is breaking down for you.

Clearly there are two ways you can go: one, poke the system the _right_
way, and two, abandon the A-H method. You _could_ use excitation and
response with more 'university' style adaptive control methods to derive
a system model and develop a controller -- this is not a bad way to go if
you're comfortable wrangling the math. Or you can try to fix your
excitation, or you can use the open-loop A-H method of data collection.

Have you tried turning your drive way down and adding some substantial
hysteresis to your relay function? I.e. instead of the relay being a
static rule that says "output = on if above threshold, output = off if
below", you make the rule say "go from off to on if above a low
threshold, go from on to off if above a high threshold". Make the two
thresholds far enough apart (5 or 10 degrees?) to insure you capture good
data, then find a value for your drive that seems to give good results
for a variety of boiler sizes.

Or you can use the open-loop A-H -- this is actually the original method
that Ziegler and Nichols used (as far as I know). In it you start from
steady state, make a step change in the drive, and look for a ramp
response in the output. You deduce gain and delay from the slope of the
ramp and the delay from that ramp to an ideal 'immediate' ramp. In your
case you probably want to take off from a fairly hot temperature, so you
may want to get up to your launch temperature with a really, really slow
integrator-only loop or some other means that'll get you at a steady
temperature safely. Once you're at a constant temperature and drive, you
can step the drive and start collecting data points on the response.

Or hire me :-).

--
www.wescottdesign.com

Frank W.

unread,
Nov 9, 2009, 2:36:52 AM11/9/09
to
"Tim Wescott" <t...@seemywebsite.com> wrote:

> you're not really poking the system in quite the right way

Possible. The original method is: if the process value crosses sepoint, the
heat is turned on (below setpoint) or off (above). One has to determine
steady state output (us) first and apply steps to the heating element on
this basis: heat = us +- u. This limits heat output if you want to avoid
getting into saturation: if us = 6%, as it is in my case, then you can not
have larger output steps than 6% (to avoid getting below 0%), ie. power
swings from 0% to 12%. Since that's too low to generate peaks that I can
measure with good resolution, I deviate from the inherent step limit of the
original algorithm and excite the system with larger steps (e.g. 0 - 25% =
+-12.5%). I'm now trying to minimize the possible effect of this deviation
my incrementally reducing u until I have swings that are just big enough to
be measured with resaonable resolution.

> Or you can use the open-loop A-H -- this is actually the original method
> that Ziegler and Nichols used (as far as I know). In it you start from
> steady state, make a step change in the drive, and look for a ramp
> response in the output.

I don't think this would work any better because this test requires special
osciallations, too. There are very fast PV rises in my system because the
sensor is on the outside of a thick, small die-cast aluminium boiler and
close to strong heating elements. And there are very slow declines, dictated
by the small, natural heat dissipation of a water-filled boiler (no flow).
That means I can't find the special P that creates oscillations because all
do: even a small P will cause forced oscillatations of PV because the system
is dampend so much in one direction.


Tim Wescott

unread,
Nov 9, 2009, 11:33:08 AM11/9/09
to
On Mon, 09 Nov 2009 08:36:52 +0100, Frank W. wrote:
<begin restored context>

Or you can use the open-loop A-H -- this is actually the original method

that Ziegler and Nichols used (as far as I know). In it you start from
steady state, make a step change in the drive, and look for a ramp

response in the output. You deduce gain and delay from the slope of the
ramp and the delay from that ramp to an ideal 'immediate' ramp. In your
case you probably want to take off from a fairly hot temperature, so you
may want to get up to your launch temperature with a really, really slow
integrator-only loop or some other means that'll get you at a steady
temperature safely. Once you're at a constant temperature and drive, you
can step the drive and start collecting data points on the response.

<end restored context>

I hear an attempt at a logical, 'software engineer's' description of the
behavior of the system. "First event A happens, which causes event B,
then event C happens, which causes event D, etc.". This is sensible for
most problems that one solves with embedded code, but it isn't a way that
works for most control loops -- mostly because in a dynamic system like
this, even if the driving force can be separated into discrete events,
the response to the driving force gets spread out all over everything.

I do not believe that your temperature response to heat input is as
nonlinear as you think it is. This is partially because in my experience
thermal systems themselves are pretty linear (and therefore act the same
in both directions), and partially because looking at your graph I don't
see any true signs of underdamped behavior -- I see a system that is
responding to asymmetrically, yes, but it looks like a linear response to
an asymmetrical input.

To hazard a guess, I'd say that one way or another when you command a 6%
drive to your heater you're getting 1/4 as much heat as when you command
a 12% drive -- almost as if you're commanding a steady voltage to the
heater, or you're commanding the timing on an SCR drive, or something.
The response that you're getting looks very much like the response you'd
expect with a 3/4 'on' drive, and a -1/4 'off' drive, which is consistent
with an ambient cooling at your set point temperature that's 1/4 as
strong as the heating you get when the thing is on.

If you're really exciting a nonlinearity in the heating system you want
to find it and characterize it now, rather than struggling with
mysterious voodoo code for the rest of your career.

My recommended method would be to do it by analysis, but since you're
already set up for measurements, making a chart of heater commands vs.
final temperature, from 0 up through whatever gets to the maximum water
temperature you can tolerate (I'd go for a gentle simmer, unless that
causes some non-reversible safety mechanism to go off). If the chart
shows a curved line you have a nonlinearity in your drive that you need
to deal with (otherwise you have a nonlinearity somewhere else).

--
www.wescottdesign.com

Frank W.

unread,
Nov 9, 2009, 2:32:29 PM11/9/09
to
"Tim Wescott" <t...@seemywebsite.com> wrote:

> To hazard a guess, I'd say that one way or another when you command a 6%
> drive to your heater you're getting 1/4 as much heat as when you command a
> 12% drive

95C steady state with 7% heat and 80C steady state with 3.5% heat. Heating
from room temperature to 95C takes ~1min, cooling takes hours (cooling from
95C to 80C takes 15min).

Ultimate gain can be computed just fine, it doesn't matter much what output
steps are used during autotuning. The �str�m/H�gglund result is always ~25%
Kc gain and that works fine and seems right; manual tuning suggests a
similiar value.

It was only ultimate period that gave me problems. Different excitations
during autotuning result in different periods. However, after the recent
changes (minimizing excitation with the goal of creating autotuning peaks
that are as small as possible but large enough to exceed sensor resolution,
thus turning the curve more into a sinus), it seems I get fairly consistent
results that work: Td is calcluated as 2.5s. Manual tuning suggests 1.6s.
That's good enough for me, so I consider the case closed. Thanks for your
help, Tim.


Oliver Betz

unread,
Nov 10, 2009, 3:49:01 AM11/10/09
to
"Frank W." wrote:

>"Tim Wescott" <t...@seemywebsite.com> wrote:
>
>> To hazard a guess, I'd say that one way or another when you command a 6%
>> drive to your heater you're getting 1/4 as much heat as when you command a
>> 12% drive

>95C steady state with 7% heat and 80C steady state with 3.5% heat. Heating

Did you consider possible long time constants for "steady state"? Slow
warming of enclosures can ruin your model extraction.

And what is "heat", IOW how do you control the heater - phase angle,
full wave? As Tim wrote, there might be a nonlinear function. And if
you control a SCR without mains synchronisation, you might get huge
uncertainties for low power values.

>from room temperature to 95C takes ~1min, cooling takes hours (cooling from
>95C to 80C takes 15min).

1k/min : 70k/min = 1.4%, that's not so far off from 3.5% to assume
nonlinear behaviour. And your cooling experiment might be influenced
by a hot environment.

Oliver
--
Oliver Betz, Munich
despammed.com might be broken, use Reply-To:

Frank W.

unread,
Nov 10, 2009, 4:40:10 AM11/10/09
to
"Oliver Betz" <ob...@despammed.com> wrote:
> And what is "heat", IOW how do you control the heater - phase angle,
> full wave?

Full wave. Phase angle is not permitted here with so much wattage.

> you control a SCR without mains synchronisation, you might get huge
> uncertainties for low power values.

The relay has a ZC circuit. 50hz power and 1s cycle time means 100 zero
crossings/cycle, ie. 1% granularity, equals 0.5% average error.


Datesfat Chicks

unread,
Nov 10, 2009, 11:54:36 AM11/10/09
to
"Frank W." <frankw...@mailinator.com> wrote in message
news:7lggquF...@mid.dfncis.de...

>
> Since all PID temperature controllers have Autotune, there must be a
> solution for this problem. Any ideas?

As you probably know from control theory, the basic theory of a PID
controller is that you have a system described by a set of linear
differential equations that is inherently unstable or has some performance
problems. As a result you strap a PID controller onto it (with said
controller also described by its own linear differential equations), and the
resulting system (now described by linear differential equations which are a
mathematical mix of the underlying system and the PID controller) has better
characteristics.

Did you notice that there is a word that appears many times in my
description above?

Want to guess what the word is?

That word is "linear".

A system with a time delay is not described by linear differential
equations. Strapping a PID controller onto it is bad math.

One of the more classic examples is a shower or an industrial process that
mixes fluids of varying temperature and the sensor is located substantially
downstream from the mixing value. This is a pure time delay. My shower at
home is like that. I turn the water a little hotter. Nothing happens. I
turn it a little more hotter. Nothing happens. Then I turn it a little
more hotter. Then the wave of hot liquid hits me and I scream in agony.

Over time, I've adapted to my shower. I don't burn myself anymore.

I think the control algorithms you want to use for a system like yours fall
outside the range of PID. I'm sure there is a body of theory that covers
it, but I don't know what that is.

I would heat the system full bore for a fixed period of time, then stop and
wait to see how the temperature catches up. And work from there.

The best control strategy for that system isn't going to be PID. That is a
non-linear system.

Datesfat

Datesfat Chicks

unread,
Nov 10, 2009, 12:00:47 PM11/10/09
to
"Datesfat Chicks" <datesfa...@gmail.com> wrote in message
news:qZSdncz-fbqoB2TX...@giganews.com...

From Wikipedia's entry on PID controllers:

http://en.wikipedia.org/wiki/PID_controller

<QUOTE>
Another problem faced with PID controllers is that they are linear. Thus,
performance of PID controllers in non-linear systems (such as HVAC systems)
is variable. Often PID controllers are enhanced through methods such as PID
gain scheduling or fuzzy logic. Further practical application issues can
arise from instrumentation connected to the controller. A high enough
sampling rate, measurement precision, and measurement accuracy are required
to achieve adequate control performance.
</QUOTE>

I don't know what the best control strategy is, but it ain't PID.

Datesfat

Tim Wescott

unread,
Nov 10, 2009, 12:54:21 PM11/10/09
to

I've been resisting forking this over into the control newsgroup: now
it's compelling.

Systems with delay can be perfectly linear, as well as time invariant --
they just can't be described by ordinary differential equations with a
finite number of states.

To be linear, a system only needs to satisfy the superposition property.
A delay element satisfies superposition just fine.

And while a PID controller may not be the theoretically best controller
for a system with delay, in many cases it's not a bad choice at all. PID
controllers can and will give perfectly satisfactory service with plants
that have significant delay. The thousands, if not millions, of PID
controllers in mills and factories around the world that are controlling
plants whose responses are dominated by delay certainly belie any
declaration that the PID controller isn't a good choice to control a
plant with delay.

None of the above is intended to minimize the difficulty in analyzing and
designing a truly optimal controller for a plant with pure delay --
that's an exercise that can make your brain hurt, and fast. And nothing
of the above is intended to chase you away from taking plant delays more
directly into account if a discrete-state controller such as a PID won't
let you eke the performance that you need out of your plant.

But in the absence of significant nonlinearities or time varying behavior
you can use all the analysis tools that are suitable for linear time
invariant systems on a system with delays just fine. You can do good
design work, without ever having to explicitly write out the differential
equations, much less solving them.

So if you don't want to get lost in Mathemagic Land searching for
performance that isn't necessary for your product's success, a good ol'
PID controller may be exactly the optimal controller -- in terms of
adequate performance and reasonable engineering time -- even if it
doesn't satisfy any egghead academic measure of "optimal" for the
particular plant you're trying to control.

--
www.wescottdesign.com

Datesfat Chicks

unread,
Nov 10, 2009, 1:24:00 PM11/10/09
to
"Tim Wescott" <t...@seemywebsite.com> wrote in message
news:rPCdnSIjis5QNWTX...@web-ster.com...

>
> Systems with delay can be perfectly linear, as well as time invariant --
> they just can't be described by ordinary differential equations with a
> finite number of states.

Hi Tim,

I might have missed something significant here.

It is my assumption that a system with a pure time delay is inherently
non-linear.

Let's take my shower example with a pure delay in the pipes ...

With no delay, you can just say that

Temperature(t) = Valve_Position

or perhaps with a little thermal mass thrown in you can say that:

d Temperature / dt = K * (Valve_Position - Temperature)

where of course I'm assuming that valve position and water temperature are
the same thing.

The first is I think a 0'th order linear differential equation and the
second is a 1st-order LDE.

But how would you linearize a system with a pure time delay, exactly?

The shower example with a pure pipe delay between the shower valve and my
skin is fine.

Thanks, Datesfat

Jerry Avins

unread,
Nov 10, 2009, 3:00:43 PM11/10/09
to
Datesfat Chicks wrote:
> "Tim Wescott" <t...@seemywebsite.com> wrote in message
> news:rPCdnSIjis5QNWTX...@web-ster.com...
>>
>> Systems with delay can be perfectly linear, as well as time invariant
>> -- they just can't be described by ordinary differential equations with a
>> finite number of states.
>
> Hi Tim,
>
> I might have missed something significant here.
>
> It is my assumption that a system with a pure time delay is inherently
> non-linear.

Superposition is sufficient proof of linearity. What comes out of a pipe
(assuming that there is no mixing in transit) is almost a delayed linear
superopsition of what is pushed into it, but it is not linear because it
is not a pure delay. When the input velocity increases because both hot
and cold water are flowing, the delay time decreases. Superposition
doesn't strictly apply because the time to look isn't well defined.

Any delay pushes a servo system toward unstable. That's not a linearity
problem.

> Let's take my shower example with a pure delay in the pipes ...
>
> With no delay, you can just say that
>
> Temperature(t) = Valve_Position
>
> or perhaps with a little thermal mass thrown in you can say that:
>
> d Temperature / dt = K * (Valve_Position - Temperature)
>
> where of course I'm assuming that valve position and water temperature
> are the same thing.

There's also the time it takes the valve to move.

> The first is I think a 0'th order linear differential equation and the
> second is a 1st-order LDE.
>
> But how would you linearize a system with a pure time delay, exactly?

It's already linear. Just nasty.

> The shower example with a pure pipe delay between the shower valve and
> my skin is fine.

But, as I wrote above, a pipe is onlt a pure delay as long as the flow
is constant.

Jerry
--
Engineering is the art of making what you want from things you can get.
¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯

pnachtwey

unread,
Nov 10, 2009, 3:45:24 PM11/10/09
to
On Nov 5, 9:37 am, "Frank W." <frankw_use...@mailinator.com> wrote:
> Is there a control engineering expert here?
Yes.

> Since all PID temperature controllers have Autotune, there must be a
> solution for this problem.
Some auto tuning algorithms are pretty crude and I think some
companies implement the simplest auto tuning algorithms just to say
they have one. It is a marketing thing. How are the customers to
know the difference? If they knew the difference they wouldn't need
the auto tuning.

> Any ideas?
Turn down the gain for one. It looks like the control output , gray
bar, amplitude is variable and not a simple on/off SSR.
Because the plant responds so quickly to the output I think the PID
update time should be faster but then the lack of resolution
will be a limitation. Using the derivative gain may be handy to but
again the lack of resolution will be a factor.
A tried and true method is the feed forward. If one know how much
control output is require for every temperature then one can compute
the control output and be with a percent or two without even using the
PID gains. The PID gains are just used to correct for errors in the
feed forward predictions.
You should be able to get very close with a feed forward and a
proportional gain. The integrator will simply remove the last bit of
error.

Auto tune works great when the plant matches the model AND the system
is properly excited. If not then auto tune doesn't work so well. Are
the water levels in the boiler constant? If not can the water level
be detected? Obviously the temperature will rise more slow or not as
fast if the boiler is full of water as opposed to being empty.

Does the tuning need to be better for a coffee machine? Manual tuning
should be easy enough.

When I get serious I prefer to make a model of the system. The
excitation requires a few steps in the control output. The steps
should be big enough so the response is clearly more significant than
the noise. The dead time should be clearly visible and the same goes
for the rise time. One can often get the plant gain, time constant and
dead time by inspection. If not then one must log the control output
and the response. Then one fits a FOPDT or SOPDT model to your coffee
machine. I do this using Scilab's optim or lsqrsolve. Then one goes
to the www.controlguru.com site and uses the IMC, internal model
control, equations there for computing the gains. For a FOPDT
Kc=t1/((tc+dt)*K)
ti=t1
Where:
Kc is the controller gain. It should have units of %output per
degree.
ti is the integrator time constant which is usually in minutes.
t1 is the plant time constant
dt is the dead time
K is the plant gain. This has units of degree per % control output.
tc is the closed loop time constant. tc is made smaller for more
aggressive tuning and longer for more conservative tuning. The value
for tc is dependent on the plant time constant and dead time. The
formula for tc on on the www.controlguru.com site.

You can see from the equation that the dead time reduces the gain.
Also, if the plant gain changes due to the water level then K will be
a function of water and the controller gain will then be modified as a
function of water level.
Kc=t1/((tc+dt)*K(WaterLevel))

I would start here. Once the plant gain, time constant and dead time
are know there are other options. On the www.controlguru.com site
there are formulas for adding the derivative time constant td. One
can also implement a Smith Predictor. A Smith Predictor isn't that
hard to implement. The hard part is finding the plant parameters.

I talked to someone else on the phone about making a high end coffee
or espresso machine a few months back. Our product was gross overkill
but would have made a nice research tool.

Peter Nachtwey

Tim Wescott

unread,
Nov 10, 2009, 5:52:57 PM11/10/09
to

Well, if it's already linear you don't linearize it.

Take the system y = h(x, t) ==> y(t) = x(t - td). Testing this with
superposition we get

y1(t) = x1(t - td),
y2(t) = x2(t - td),

y1(t) + y2(t) = x1(t - td) + x2(t - td)

which is both h(x1, t) + h(x2, t) and h(x1 + x2, t) -- therefore the
system is linear.

Note that as Jerry points out a shower isn't necessarily a linear system,
unless your shower valve insures a constant flow and the pipes don't have
any turbulence. Let a vastly simplified version be

y(t) = x(t - kd * x(t)),

(this doesn't capture the delay behavior in even a perfect pipe)

Then we try superposition:

y1(t) = x1(t - kd * x1(t)),
y2(t) = x2(t - kd * x2(t)).

This does _not_ equal the system output to the sum:

ys(t) = x1(t - kd * (x1(t) + x2(t))) + x2(t - kd * (x1(t) + x2(t))).

so this system isn't linear -- but not for the reason that you thought.

--
www.wescottdesign.com

Tim Wescott

unread,
Nov 10, 2009, 5:59:27 PM11/10/09
to
On Tue, 10 Nov 2009 12:45:24 -0800, pnachtwey wrote:

> On Nov 5, 9:37 am, "Frank W." <frankw_use...@mailinator.com> wrote:
>> Is there a control engineering expert here?
> Yes.
>> Since all PID temperature controllers have Autotune, there must be a
>> solution for this problem.
> Some auto tuning algorithms are pretty crude and I think some companies
> implement the simplest auto tuning algorithms just to say they have one.
> It is a marketing thing. How are the customers to know the difference?
> If they knew the difference they wouldn't need the auto tuning.

-- snip --

> Auto tune works great when the plant matches the model AND the system is
> properly excited. If not then auto tune doesn't work so well. Are the
> water levels in the boiler constant? If not can the water level be
> detected? Obviously the temperature will rise more slow or not as fast
> if the boiler is full of water as opposed to being empty.
>
> Does the tuning need to be better for a coffee machine? Manual tuning
> should be easy enough.

-- snip --

> I talked to someone else on the phone about making a high end coffee or
> espresso machine a few months back. Our product was gross overkill but
> would have made a nice research tool.
>
> Peter Nachtwey

Speaking of overkill, why not just use an on-off thermostat with some
hysteresis? Energy consumption? Lack of cool?

For that matter, how about modeling the whole range of coffeemakers that
you expect (or that you're willing to warrant the device as working).
Then see if you can get reasonable performance with just one tuning
setting*. Then instead of telling the customer he can mess around with
autotune, you can tell him he can just plug & play.

* I've done this on blood pressure monitor; we expected to need four
different tunings for the whole range of cuffs from 'premie arm' to 'fat
guy thigh', and ended up needing just one.

--
www.wescottdesign.com

pnachtwey

unread,
Nov 10, 2009, 6:51:33 PM11/10/09
to
On Nov 10, 2:59 pm, Tim Wescott <t...@seemywebsite.com> wrote:
> Speaking of overkill, why not just use an on-off thermostat with some
> hysteresis?
Good idea, then there are no gains to tweak so coffee can be made
right away and the coffee drinking can begin. How accurate must the
temperature control be on a coffee maker and who would be able to tell
the difference between 85 and 86 degrees C? 85-86 degrees seems to be
a few degrees too low unless Frank W knows something the rest of the
world doesn't.

Peter Nachtwey

Frank W.

unread,
Nov 11, 2009, 3:54:22 AM11/11/09
to
pnachtwey" <pnac...@gmail.com> wrote:

> unless Frank W knows something the rest of the world doesn't

Don't speak for the rest of the world if you are an ignoramus.

- With a thermostat, temperature swings are far more pronounced than you
claim. It's not just two degrees: The PID has reduced temperature swings of
~15 degree Celsius with the original thermostat to nearly nothing. Such
variations had a noticeable influence on taste (cold = sour). Here's a
PID-vs-thermostat temperature plot for a different, slightly higher-end 500$
home coffee machine which has even bigger 23.9C swings:
http://www.espressoparts.com/espressoparts/content/graphics/pagebuilder/pidsd3c_silvia_tempplot.jpg
- An equally important benefit is that with a PID instead of a thermostat,
the PV is adjustable. Different beans/roasts call for different brew
temperature - roughly between 85 and 99 degree Celsius.

RockyG

unread,
Nov 11, 2009, 4:51:50 AM11/11/09
to
> Here's a
>PID-vs-thermostat temperature plot for a different, slightly higher-end
500$
>home coffee machine which has even bigger 23.9C swings:
>http://www.espressoparts.com/espressoparts/content/graphics/pagebuilder/pidsd3c_silvia_tempplot.jpg
>- An equally important benefit is that with a PID instead of a
thermostat,
>the PV is adjustable. Different beans/roasts call for different brew
>temperature - roughly between 85 and 99 degree Celsius
>
Just for interest, is the boiler pressurised? I saw that the graph show
temepratures higher than 100 deg C. and wondered how they were achieved.


---------------------------------------
This message was sent using the comp.arch.embedded web interface on
http://www.EmbeddedRelated.com

pnachtwey

unread,
Nov 11, 2009, 10:33:43 AM11/11/09
to
On Nov 11, 12:54 am, "Frank W." <frankw_use...@mailinator.com> wrote:

> pnachtwey" <pnacht...@gmail.com> wrote:
> > unless Frank W knows something the rest of the world doesn't
>
> Don't speak for the rest of the world if you are an ignoramus.
You have pissed me off. I searched the internet, the rest of the
world, to find the correct brewing temperature and normal limits which
apparently you don't know and know one else has bothered to do. DON'T
YOU THINK THAT IS IMPORTANT????
I see have wasted my time.

Peter Nachtwey

Tim Wescott

unread,
Nov 11, 2009, 11:35:13 AM11/11/09
to
On Wed, 11 Nov 2009 09:54:22 +0100, Frank W. wrote:

> pnachtwey" <pnac...@gmail.com> wrote:
>
>> unless Frank W knows something the rest of the world doesn't
>
> Don't speak for the rest of the world if you are an ignoramus.

Actually Peter is one of those "control experts" that you were looking
for. If anything he can claim more pertinent experience at this stuff
than me -- and if a potential customer came to me with a problem that
seemed to be in his bailiwick I wouldn't hesitate* to send them his way.

Control problems seem to attract overly complicated solutions; so one of
the things a responsible control engineer always needs to do is to ask
"is there some bone-head simple way to do this that will be adequate".

> - With a thermostat, temperature swings are far more pronounced than you
> claim. It's not just two degrees: The PID has reduced temperature swings
> of ~15 degree Celsius with the original thermostat to nearly nothing.
> Such variations had a noticeable influence on taste (cold = sour).
> Here's a PID-vs-thermostat temperature plot for a different, slightly
> higher-end 500$ home coffee machine which has even bigger 23.9C swings:
> http://www.espressoparts.com/espressoparts/content/graphics/pagebuilder/
pidsd3c_silvia_tempplot.jpg
> - An equally important benefit is that with a PID instead of a
> thermostat, the PV is adjustable. Different beans/roasts call for
> different brew temperature - roughly between 85 and 99 degree Celsius.

Not to say that you aren't right in the end, but there are thermostats
and there are thermostats. Design your software for relay control, keep
the hysteresis low and the setpoint adjustable, and you may find that a
swing of a degree or two is quite achievable. I don't know if it would
be adequate, but you can do a heck of a lot better than a $0.1 snap-
action bit of bimetal and a switch.

Of course, if "PID control" has become a sales point, then you'd better
have it in there.

* I may do so wistfully if work were slow, but I'd still do it.

--
www.wescottdesign.com

RRogers

unread,
Nov 11, 2009, 1:18:16 PM11/11/09
to
On Nov 10, 11:54 am, Tim Wescott <t...@seemywebsite.com> wrote:
> On Tue, 10 Nov 2009 11:54:36 -0500, Datesfat Chicks wrote:
> > "Frank W." <frankw_use...@mailinator.com> wrote in message

I have recently done a thermal MIMO PID controller that ended up
preforming adequately despite using very simple controls.
Some comments:
Even the simplest differential description ends up with an infinite
number of state/poles.
Most real thermal systems have little tabs and things that foul up
theoretical analysis.
Therefore: you can start with simple mathematical models to estimate
requirements but you always end up with approximations.
Pole zero analysis in this case is almost worthless except to roughly
get started.
Bode and/or Nichol's chart analysis (I used both) works very well;
but ..
You have to get and use the experimental data. You can use that
directly or find a sufficiently good model for the system.
You should establish a "process" for the tuning and experiments; the
system you take the data on will undoubtedly not be the one that ends
up being manufactured.
Gotcha's: Scilab's system identification processes are unstable
dealing with this type of system. They can be used to attempt
modelling but tread carefully and double check.
When taking the data, the room/environmental temperature will do
everything it can to confound the experiment.
Don't worry about the lower frequencies, go to where the phase starts
to shift significantly.
For the Bode/Nichols derived compensation just redo the experiment
(which you probably will) to clarify the standard compensation region
round the Bode criterion; 180 degrees +- one or two decades.
Try to give at least hints to how the tuning was done for the
"outsourced" maintenance people who have to maintain the tuning after
the mechanical assembly is altered; unless you want to come back and
start over yourself in a year.

Really, really examine the code to make sure you don't "windup". I
was forced to rely on programmers in another group and I had study the
experimental results for a while to realize that the anti-windup code
just clipped the output not the integrator.

Ray

Frank W.

unread,
Nov 11, 2009, 7:16:25 PM11/11/09
to
"RockyG" <Rober...@gmail.com> wrote:

> Just for interest, is the boiler pressurised? I saw that the graph show
> temepratures higher than 100 deg C

The boiler of the Rancilio machine, which uses an internal, immersed heating
element, is pressurized. But not all are. For example the boiler of a Gaggia
Coffee, which uses external heating elements, possesses a steam exhaust tube
for pressure relief which is closed (well, almost) by a ball floating before
a hole during shots only.


pnachtwey

unread,
Nov 12, 2009, 9:43:13 PM11/12/09
to
I agree with the last paragraph.
However, I have had a lot of success with identifying systems poles
and zero. I can then place both where I want with the controller
gains.

I didn't know Scilab has a system identification function, but I have
used the lsqrsolve and optim successfully.

Peter Nachtwey


RRogers

unread,
Nov 13, 2009, 9:44:52 AM11/13/09
to

Interesting, I have thought about going that route but opted for a
more conventional process; System Identification routines. But that
wasn't very satisfactory. I have a problem in that I like to continue
along routes until I really understand why they don't work. Sometimes
I think that half my brain is autistic.
Once I get my system identification code reorganized (with or without
a gui) I plan to test it against my data and some available test cases
from NICONET. Although they don't seem to be MIMO. In biological
testing equipment you are forced into MIMO situations in order get the
required temperature accuracy over large testing areas and
environmental conditions. In addition mammalian reactions are tuned
to constant temperature within a narrow band; 37degC in our case
(presuming no aliens in the group). I was actually looking forward
to doing that; I had never had use MIMO before. Wasn't so enthused
after a while; the design process is a lot more complicated and the
tools were not robust.
Once I resolve (or at least identify) the problems perhaps I will
compare the results with lsqrsolve. If your interested I will post a
link here; but don't expect anything soon. I am just settling into
Mexico, and am not as fast as I used to be.

Ray

pnachtwey

unread,
Nov 13, 2009, 12:45:46 PM11/13/09
to
When you have MIMO test data why don't you share it with us. I would
like to have a crack at too.
It would be helpful to know what I am fitting data too though so I can
get the general form the equations right. I don't know anything about
your field of study.

The trick is how you use optim() and lsqrsolve(). The best system
identification uses Runge-Kutta to integrate the model's system of
differential equations.

For MIMO systems you will need to use optim(). optim() can optimize a
cost function. lsqrsolve() requires two arrays of data, the actual
data and the estimated data. I don't know how you would do this if
you have two sets of actual data and two sets of estimated data.

Peter Nachtwey

RRogers

unread,
Nov 14, 2009, 10:44:34 AM11/14/09
to

I don't quite understand your approach; it seems different from what I
had in mind. I have multiple sets of experimental data consisting of
three stimulus/drive columns and three columns of resulting temerature
data; together with a multitude of other columns of other temperature
readings for thermal design of the overall assembly.
My hypothetical approach to raw curve fitting type of modeling:
Write out the ABCD equations with unknown coefficients and try to find
the coefficients; which are linear (superficially) coefficients
applied to the data. Having an adequate model in hand, then I thought
I would use optim() to find the control gains in the closed loop.
This is not what you are describing. My formulation was just a
passing thought and certainly has a lot of problems I haven't
resolved.
Your comments don't fall in line with this, so why not tell me yours.

Brief technical details follow (of interest only to those who enjoy
these things):
The system consists of three heaters and three sensors; actually
far more sensors for the data, but the others were temporary and
informational for the rest of the machine and not used in control.
The system consists of a disk holding something like 20 test strips
and rotating the strips under a dispenser and then under an optical
head; so each of the test strips rotated to have a drop of sample
deposited and then put under the optical head to monitor the reaction
development. One of the heating systems was a buffer to isolate the
test disk from the room. The other two are more precise and localized
controls that control the sample tray fairly precisely to 37 degC.
The reason for two heaters: one controls most of the circular sample
disk consisting of 20 or so test strips that have been entered; the
other heater brings the incoming test strips up to temperature from
the room temperature when they are inserted. The original specs were
that the samples had to be at 37degC +- .1degC when the reaction was
occuring, warmup in 5 minutes, ambient/room temperature 18degC to
30degC. I designed the control system to be .02 degC accurate at the
tray thermistor, control loop closure at power up inside of two
minutes, PID controls around the principal MIMO directions (the
thermistors were placed reasonably close to the individual heaters).
The last part was to make the programming (done by another group) and
maintenance easier; requiring less skilled people and the end
performance was adequate. The problems involved were:
1) I couldn't put the thermistors where I actually wanted them,
2) I couldn't be hyper conservative and truly insulate the assembly
( the mechanical people had more than enough problems) so I had to
rely on chunks of metal smoothing out the spatial frequencies. Of
course the assembly had variations across it anyway.
3) I never had the final machine available during testing because the
mechanical people needed to know about thermal problems before the
design was finished, and I didn't want to be the person holding up
release after the machine was finished. The was only a problem during
testing since one set of readings would be different from a set taken
later.
4) The sys-id routines were not robust and had to be watched very
carefully. In fact I ended up using the DC gains of the models as the
first quality determiner. Then I would look at various residuals to
determine the real quality. Usually the test data was split in half
(or so) so the model wouldn't be just regurgitating data back to me.
The first half was used to determine a model and then the model was
used to predict the second half; the resulting residual time series
were then examined. I wanted the residuals to be below .1 degC (1
part out of 370) or so but never got there due to inadequacies in the
model, and I had to settle for 2 degC; the slack/error was taken up
when the loops were closed. Apparently the sys-id routines want
random inputs; whereas people are more comfortable with large step
inputs. I have both types of data.
5) All of the heater systems talk to each other and the environment
thermally; the reason for MIMO approach.
6) Severe organizational problems with people who had never done
instrument design before (: That's a different story.

What is driven home is the fact that you are just looking for an
adequate model of reality in thermal situations; not looking for
"truth". The mechanical assembly can not reduced to anything less
than a FEA analysis; which I couldn't get the department to
institute. It's not a trivial thing to incorporate in a design
process. Having done a partial survey I think COMSOL is a pretty good
multiphysics tool and does have the ability to incorporate spice
models between objects like a thermistor (actually a point) and a
heater.

And so on, I have more information. None of this relates to any
proprietary information; except if I come up with a better process I
can answer questions from the engineer who has to redo the system
after they make changes to the mechanical design. The design changes
are inevitable and occasionally people get back to me with questions.

If you really want some data I can post it on an FTP sight. The
project is done and I am retired so there is no hurry. The data is
not clean and has a lot of confounding disturbances; OTOH there is a
lot of it :) I am still interested in determining a better process
for establishing good models; although I am inclined to fix up the sys-
id functions so that higher order approximations don't lead to
(wildly) worse and worse predictions. That is just nonsense.
Be aware that my criteria are DC gain and residuals; and any comment
on the modelling will probably be oriented around that. If your
interested in my code; my SCILAB program does produce a lot of
outputs, BODE and Nichols charts; but is not finished code in the
sense that some parameters are done with I/O, and some parameters are
entries in the code. There are shortcomings, I never did a good Bode
plot of the raw data, just of the models. I kept meaning to but that
requires a lot of filtering to be meaningful.

Hope I haven't bored you to much Peter.

Ray

JCH

unread,
Nov 14, 2009, 12:24:49 PM11/14/09
to

"RRogers" <rero...@plaidheron.com> schrieb im Newsbeitrag
news:d45695c8-7d48-4841...@f20g2000prn.googlegroups.com...


See simple example with differential equation of order 2:

* http://home.arcor.de/janch/janch/_control/20081123-real-system-model/

I try to find the best possible process transfer function (page 1) by using
approximation methods on the basis of some measured values (page 2).

Thereafter I have a benchmark test scheme (page 3) with a program (page 4)
that automatically finds the best PID parameters using the IAE criteria.

This could be done for process identifications up to differential equations
of degree 6.


--
Regards JCH

RRogers

unread,
Nov 14, 2009, 4:25:39 PM11/14/09
to
clip..........
> ...
>
> read more »

Okay I have uploaded the file that corresponds to step inputs. This
one is fairly clean.
http://www.plaidheron.com/ray/temp
static-test.jpg
static-test.xls
Should get you there. If there is a permission problem let me know; I
will resolve.

The .jpg is a graph to get the idea. T-11 is included to verify the
environment didn't change much.
The .xls is: sheet 1 graphs, sheet static-test is the long
experimental run covering about 4 hours
Cols: T-1,2,3 are the three direct thermistors used later for control
Cols: M,N,O are the PWM drives, 0-100%, to the corresponding heaters;
the trailing columns can be ignored
The intermediate columns are various sensors distributed away from the
actively controled points.

Let me know and I (or you ) can cross-verify your model against other
experimental runs.

I have other experimental data sets that are less clear; some are
basically random inputs to try to satisfy the sys-id programs.

Ray

JCH

unread,
Nov 15, 2009, 7:14:59 AM11/15/09
to

"RRogers" <rero...@plaidheron.com> schrieb im Newsbeitrag
news:3d4e61d7-69d7-4431...@x5g2000prf.googlegroups.com...


Basically refering to

* http://home.arcor.de/janch/janch/_control/20081123-real-system-model/

Can you approach the best possible ODE (process transfer function) in a
range of order <= 6?

C6 y'''''' + C5 y''''' + C4 y'''' + C3 y''' + C2 y'' + C1 y' + y = K

Decimal commas!

Example data points: 30

0 0
0,062 0
0,124 0,002
0,187 0,012
0,249 0,04
0,311 0,093
0,373 0,17
0,435 0,266
0,498 0,373
0,56 0,48
0,622 0,581
0,684 0,671
0,746 0,748
0,809 0,811
0,871 0,861
0,933 0,899
0,995 0,929
1,057 0,95
1,12 0,966
1,182 0,977
1,244 0,984
1,306 0,99
1,368 0,993
1,431 0,996
1,493 0,998
1,555 0,999
1,617 1
1,679 1
1,741 1
1,804 1,001


--
Regards JCH

My solution see down here:

Decimal commas!
1,048734E-06 y'''''' + 6,2427E-05 y''''' + 0,001548347 y'''' + 0,02048154
y''' + 0,1523982 y'' + 0,6047773 y' + y = 1,000953

pnachtwey

unread,
Nov 15, 2009, 9:12:08 AM11/15/09
to
On Nov 14, 7:44 am, RRogers <rerog...@plaidheron.com> wrote:
> I don't quite understand your approach; it seems different from what I
> had in mind.  I have multiple sets of experimental data consisting of
> three stimulus/drive columns and three columns of resulting temerature
> data; together with a multitude of other columns of other temperature
> readings for thermal design of the overall assembly.
>      My hypothetical approach  to raw curve fitting type of modeling:
> Write out the ABCD equations with unknown coefficients and try to find
> the coefficients; which are linear (superficially) coefficients
> applied to the data.  Having an adequate model in hand, then I thought
> I  would use optim() to find the control gains in the closed loop.
> This is not what you are describing.  My formulation was just a
> passing thought and certainly has a lot of problems I haven't
> resolved.
> Your comments don't fall in line with this, so why not tell me yours.
Why not use the principle of superimposition. Test each heater with
respect to each sensor
and then find the FOPDT or SOPDT coefficients that work
For the first temperature sensor you have a FOPDT formula that looks
like
t1'=A1*t1+B11*u1(t-dt11)+B12*u2(t-dt12)+B13*u3(t-dt13)+C
Where:
t1 is the temperature a sensor 1
A1 is the system time constant at temperature sensor 1. This is
basically exp(-t/tau1).
B11 is the input coupling of heater 1 to sensor 1.
B12 is the input coupling of heater 2 to sensor 1.
B13 is the input coupling of heater 3 to sensor 1.
u1(t-dt11) is the heater 1 signal for time t.
dt11 is the dead time from heater 1 to sensor 1.
C is the ambient temperature. It had better be the same for all
test unless the ambient temperature is really changing.
It is easy to ID B11 B12 and B13 if they are turned on 1 at
a time but the starting point should be ambient temperature or
some steady state. When done you would have this
t1'=A1*t1+B11*u1(t-dt11)+B12*u2(t-dt12)+B13*u3(t-dt13)+C
t2'=A2*t2+B21*u1(t-dt21)+B22*u2(t-dt22)+B23*u3(t-dt23)+C
t3'=A3*t3+B31*u1(t-dt31)+B32*u2(t-dt32)+B33*u3(t-dt33)+C

All the coefficient could probably be ID at once but then it would
be much harder to get exact values. It is best to do small sections
at a time and rely on superimposition.

The way I ID a system is like this
http://www.deltamotion.com/peter/PDF/Mathcad%20-%20Sysid%20SOPDT.pdf

1. On page 1/10 I define the ideal SOPDT system. I chose different
value to to see how the well the system identification works under
different conditions Notice that there is dead time and I don't assume
all the poles are at the same location like others on this newsgroup.
2. At the bottom of page 2/10 I generate the test data that is later
to be used for system identification. I add noise the to ideal data
just to simulate reality a bit. The CO(t) function is a few steps.
The function can be arbitrary but I have found that the excitation is
critical to the identification. Dead times and time constants are
determined more accurately if the are step or rapid changes. The gain
and ambient coefficients are determined more accurate if the are steps
at different levels.
3. One page 3/10 I plot and save the generated test data. I can post
it on my FTP site for you to practice with. Notice that this data has
dead time and two poles that aren't at the same location. I could
have added more noise but the quasi-Newton method seems to filter it
out well.
4. One page 4/10 the system identification is done. Mathcad's Minerr
function can be like either Scilab's optim() function or lsqrsolve
function depending on the option chosen. I chose the quasi-Newton
optimization which is similar to the optim() function. Runge-Kutta is
used to integrate the differential equation. The differential equation
doesn't need to be linear. I could easily put a none linear term in
there like one that changes the gain as a function of temperature.
This happens with heat exchangers because of the LMTD. Fluid systems
are often of the form
v'=g/m-K*v^2. It is easy to ID non linear system IF you know the
general form of the equation and just need to ID the constants.
Notice that the ID'd poles are closer together than the real poles. I
have notice that system identification tends to ID the poles closer
together than what they really are. Notice that I all ID a dead time
and an ambient temperature. This is something that JCH does not do.
At the bottom of the MSE(), mean squared error function, is where I
calculate the mean squared error between the estimated temperature and
the actual or test data temperature. The Minerr function adjusts
Kp,t1,t2,thetap, and C till the MSE is minimized. You can see the
results are not perfect but that is reality.
5. On page 5/10 the actual or perfect response is compare to the
estimated response. The response looks close, almost identical, even
though the system identification puts the estimated poles closer
together. Also notice that a good system identification routine can
ID systems that are excited by more than just a step change. In fact
they must must be able to do system identification with arbitrary
excitation. Above I said the excitation is the key to doing system
identification. One key is the make multiple steps at different
levels. This is very important in computing the gain and computing
the gain when it isn't linear. Heat exchanger's gain changes because
of LMTD. ( log mean temperature difference ).
6. On page 6/10 PID gains are calculated using the estimated plant
parameters found by system identification. My formula is a little
more complex that the IMC formulas but the response is faster/better
for the same closed loop time constant. I doubt the extra complexity
in the formula is worth the effort for most applications.
7 Page 7/10 simulates the PID control of the original system using the
gains calculated from the system identification. Notice that feed
back noise is simulated as well as the dead time.
8. Page 8/10. The simulation show the response. The response isn't
perfect because there was noise in the original data used to do the
system identification. The system identification is not perfect
because the poles are closer together than they should be and I
simulated noise on the feedback but this is closer to reality.
9 Page 9/10 uses the internal model gain formulas that I got from the
www.controguru.com site. They work well too and are much simpler they
don't work quite as mine. I should have provided a IAE value for my
gains and the IMC gains for comparison.
10 page 10/10 shows the IMC response is a little slower but most would
be please with it.

I would use the above technique one at time with each heater and
temperature sensor. Actually one can excite each of the heaters one
at a time but the data for the 3 temperature sensors at the same time.

I posted a link to a scilab program that does the same thing many
years ago but no one seemed interested.

JCH, you should copy this so your program can handle dead times,
arbitrary inputs, and poles that are not all at the same place. What
you appear to be missing the quasi-Newton code( BFGS) or Levenberg-
Marquardt code that allows you to do proper optimization. I bet you
use a grid search.

Peter Nachtwey

pnachtwey

unread,
Nov 15, 2009, 9:17:52 AM11/15/09
to
On Nov 14, 1:25 pm, RRogers <rerog...@plaidheron.com> wrote:
> clip..........
>
> > ...
>
> > read more »
>
> Okay  I have uploaded the file that corresponds to step inputs.  This
> one is fairly clean.http://www.plaidheron.com/ray/temp

> static-test.jpg
> static-test.xls
> Should get you there.  If there is a permission problem let me know; I
> will resolve.
>
> The .jpg is a graph to get the idea.  T-11 is included to verify the
> environment didn't change much.
> The .xls is: sheet 1 graphs, sheet static-test is the long
> experimental run covering about 4 hours
> Cols: T-1,2,3  are the three direct thermistors used later for control
> Cols: M,N,O are the PWM drives, 0-100%, to the corresponding heaters;
> the trailing columns can be ignored
> The intermediate columns are various sensors distributed away from the
> actively controled points.
>
> Let me know and I (or you ) can cross-verify your model against other
> experimental runs.
>
> I have other experimental data sets that are less clear; some are
> basically random inputs to try to satisfy the sys-id programs.
>
> Ray
When starting the identification process the system must be at steady
state. The three temperature sensors are at different temperatures.
That could be steady state for a combination of heater outputs but it
is hard to know. If all the heaters started at the same ambient
temperature then I know the system was at steady state.

Peter Nachtwey

RRogers

unread,
Nov 15, 2009, 9:53:39 AM11/15/09
to
On Nov 15, 6:14 am, "JCH" <ja...@nospam.arcornews.de> wrote:
> "RRogers" <rerog...@plaidheron.com> schrieb im Newsbeitragnews:3d4e61d7-69d7-4431...@x5g2000prf.googlegroups.com...

>
>
>
> > clip..........
> >> ...
>
> >> read more »
>
> > Okay  I have uploaded the file that corresponds to step inputs.  This
> > one is fairly clean.
> >http://www.plaidheron.com/ray/temp
> > static-test.jpg
> > static-test.xls
> > Should get you there.  If there is a permission problem let me know; I
> > will resolve.
>
> > The .jpg is a graph to get the idea.  T-11 is included to verify the
> > environment didn't change much.
> > The .xls is: sheet 1 graphs, sheet static-test is the long
> > experimental run covering about 4 hours
> > Cols: T-1,2,3  are the three direct thermistors used later for control
> > Cols: M,N,O are the PWM drives, 0-100%, to the corresponding heaters;
> > the trailing columns can be ignored
> > The intermediate columns are various sensors distributed away from the
> > actively controled points.
>
> > Let me know and I (or you ) can cross-verify your model against other
> > experimental runs.
>
> > I have other experimental data sets that are less clear; some are
> > basically random inputs to try to satisfy the sys-id programs.
>
> Basically refering to
>
> *http://home.arcor.de/janch/janch/_control/20081123-real-system-model/

We seem to have a disconnect here.
The system is MIMO which means that a finite model would have a set of
simultaneous differential equations. In my case three independent
variables drives and three dependent variables; leading to three
simultaneous differential equations whose order varies with the number
of state variables needed for an adequate description. Including the
room temperature we actually have four drives. Including the various
components inside the instrument (motors, solenoids, and doors) we
would have more; but for the sake of simplicity I took 3 drives and 3
sensors and treated the other drives as disturbances. A design
assumption that could have been rendered wrong by results; but then I
would have had to add more sensors and possibly more heaters.
The reason for the 3 heaters and sensors is to establish control over
extended mechanical assemblies having basically an infinite numbers of
internal states. Although the higher order states are rapidly
suppressed by the heat equation when the metal thermal time constant
is short.
As an illustration: The simple case of the sun warming a piece of
ground through the seasons. The result is basically that a 20 degC
surface variation causes .5 degC variation 2 meters down with a six
month lag; with the transfer function having an infinite number of
poles and a continuously rolling phase shift going through 180 deg
over and over. This imposes constraints when you are trying to hurry
it up via control systems. These numbers are "representative" since I
am remembering; I do have the book Bell Labs book somewhere that
solves the equation.
Alternately: Writing the Green's function for the internal temperature
of a bar heated at the surfaces requires an infinite degree polynomial
resulting in an infinite number of poles in the Laplace xform. But
the significance of higher poles drops down exponentially, so they
don't matter unless you try to wrap a control loop and close the loop
with time constants that are comprable.

And so on
Ray

RRogers

unread,
Nov 15, 2009, 9:54:00 AM11/15/09
to
On Nov 15, 6:14 am, "JCH" <ja...@nospam.arcornews.de> wrote:
> "RRogers" <rerog...@plaidheron.com> schrieb im Newsbeitragnews:3d4e61d7-69d7-4431...@x5g2000prf.googlegroups.com...

>
>
>
> > clip..........
> >> ...
>
> >> read more »
>
> > Okay  I have uploaded the file that corresponds to step inputs.  This
> > one is fairly clean.
> >http://www.plaidheron.com/ray/temp
> > static-test.jpg
> > static-test.xls
> > Should get you there.  If there is a permission problem let me know; I
> > will resolve.
>
> > The .jpg is a graph to get the idea.  T-11 is included to verify the
> > environment didn't change much.
> > The .xls is: sheet 1 graphs, sheet static-test is the long
> > experimental run covering about 4 hours
> > Cols: T-1,2,3  are the three direct thermistors used later for control
> > Cols: M,N,O are the PWM drives, 0-100%, to the corresponding heaters;
> > the trailing columns can be ignored
> > The intermediate columns are various sensors distributed away from the
> > actively controled points.
>
> > Let me know and I (or you ) can cross-verify your model against other
> > experimental runs.
>
> > I have other experimental data sets that are less clear; some are
> > basically random inputs to try to satisfy the sys-id programs.
>
> Basically refering to
>
> *http://home.arcor.de/janch/janch/_control/20081123-real-system-model/

We seem to have a disconnect here.

RRogers

unread,
Nov 15, 2009, 4:47:23 PM11/15/09
to
> The way I ID a system is like thishttp://www.deltamotion.com/peter/PDF/Mathcad%20-%20Sysid%20SOPDT.pdf
> 9 Page 9/10 uses the internal model gain formulas that I got from thewww.controguru.comsite.  They work well too and are much simpler they

> don't work quite as mine.   I should have provided a IAE value for my
> gains and the IMC gains for comparison.
> 10 page 10/10 shows the IMC response is a little slower but most would
> be please with it.
>
> I would use the above technique one at time with each heater and
> temperature sensor.   Actually one can excite each of the heaters one
> at a time but the data for the 3 temperature sensors at the same time.
>
> I posted a link to a scilab program that does the same thing many
> years ago but no one seemed interested.
>
> JCH, you should copy this so your program can handle dead times,
> arbitrary inputs, and poles that are not all at the same place. What
> you appear to be missing the quasi-Newton code( BFGS) or Levenberg-
> Marquardt code that allows you to do proper optimization.  I bet you
> use a grid search.
>
> Peter Nachtwey

Peter,
Sorry I didn't answer earlier; I was answering JCH and putting
the data up. I considered the route you illustrated, and perhaps I
should have tried harder. But in your example in the above text you
implicitly assumed (by writing the equation that way) that a good
description was accessible through a single state variable per heater/
sensor, and I ran into intellectual problems trying to have the
flexibility for extension.
I did read your earlier posting and will reread it.
The problem I had was this:
Suppose the correct set of terms for sensor one was:
x'=Ax+Bu
y=Cx+Du
Where u is heater power, y is the sensor readings, and x is the
internal state vector larger than u or y . Now a set of individual
SISO readings and using supposition would result in individual state
vectors xi and Ai,Bi,Ci,Di . Some the state vectors might be shared
between the individual responses and some not. How do you determine
which ones are shared or not? I am sure I can make up a circuit with
discrete components that would illustrate this. This interpretation
problem could be resolvable, but I didn't see how to do it. If I am
missing your point please post an example with more state variables.

It does seem as though some generalization of Thevin's equivalence
circuit theorem might be possible and applicable. Maybe using
Telegin's theorem? But these thoughts are just meanderings of my
mind.

Ray

RRogers

unread,
Nov 15, 2009, 6:05:18 PM11/15/09
to

Peter,
Okay, I will post that experiment but it's not as clean. Since
I only had shared access to the prototype I couldn't let the machine
cool down long enough for a real restart, and (of course) the room
temperature changed. These thermal systems have really long "tails";
some sections (plastic) absorb heat and let it out very slowly.

Ray

RRogers

unread,
Nov 15, 2009, 6:05:29 PM11/15/09
to
On Nov 15, 8:17 am, pnachtwey <pnacht...@gmail.com> wrote:

Peter,

RRogers

unread,
Nov 15, 2009, 7:43:57 PM11/15/09
to
On Nov 15, 8:17 am, pnachtwey <pnacht...@gmail.com> wrote:

Peter,

RRogers

unread,
Nov 15, 2009, 7:54:41 PM11/15/09
to

> > When starting the identification process the system must be at steady
> > state.  The three temperature sensors are at different temperatures.
> > That could be steady state for a combination of heater outputs but it
> > is hard to know.  If all the heaters started at the same ambient
> > temperature then I know the system was at steady state.
>
> > Peter Nachtwey
>
> Peter,
>       Okay, I will post that experiment but it's not as clean.  Since
> I only had shared access to the prototype I couldn't let the machine
> cool down long enough for a real restart, and (of course) the room
> temperature changed.   These thermal systems have really long "tails";
> some sections (plastic) absorb heat and let it out very slowly.
>
> Ray

Well I looked around, while I do have SIMO heater by heater data the
subject heater input is random trying to obtain information the sys-id
routines like.
Incidentally: In case I forget; some of the data was taken has a
problem which I found out after much work and threatening to sue the
programmers; the PWM percentages were rounded down to units not tenths
and such. That's the reason for the second set of PWM data.
Maybe I should have quit when they separated the programming from
engineering (: Endeavour to write and check your own control and
monitoring algorithms; you will have a happier life.

Ray

RRogers

unread,
Nov 15, 2009, 8:16:34 PM11/15/09
to

> > When starting the identification process the system must be at steady
> > state.  The three temperature sensors are at different temperatures.
> > That could be steady state for a combination of heater outputs but it
> > is hard to know.  If all the heaters started at the same ambient
> > temperature then I know the system was at steady state.
>
> > Peter Nachtwey
>
> Peter,
>       Okay, I will post that experiment but it's not as clean.  Since
> I only had shared access to the prototype I couldn't let the machine
> cool down long enough for a real restart, and (of course) the room
> temperature changed.   These thermal systems have really long "tails";
> some sections (plastic) absorb heat and let it out very slowly.
>
> Ray

Well I looked around, while I do have SIMO heater by heater data the

JCH

unread,
Nov 16, 2009, 5:45:06 AM11/16/09
to

"RRogers" <rero...@plaidheron.com> schrieb im Newsbeitrag
news:b489cc63-2964-418f...@z3g2000prd.googlegroups.com...
> simultaneous differential equations...


If you can't find one differential equation (process transfer function) as
part of a set you won't be able to solve anything.

See basics and decoupling of MIMO system:

* http://home.arcor.de/janch/janch/_control/20091117-mimo-system/


--
Regards JCH


RRogers

unread,
Nov 16, 2009, 10:48:44 AM11/16/09
to
clip......

>
> If you can't find one differential equation (process transfer function) as
> part of a set you won't be able to solve anything.
>
> See basics and decoupling of MIMO system:
>
> *http://home.arcor.de/janch/janch/_control/20091117-mimo-system/
>
> --
> Regards JCH

JCH & Peter

This thread is getting long and unfocused. Rather than talk about
doing something perhaps we could start a new thread or blog, and
actually have some contests and results doing system-id on some data
sets? Nothing formal, just trials and analysis to improve our grasp
of the problems. A reference site is NICONET but that is real data
and the "truth" is unknown, although I think they have some results
for comparison. In any case dividing the data up into analysis and
prediction parts allows an objective criteria of tracking accuracy.
In addition we could construct various systems (say circuits or
something like what CLF showed), run simulations, present the data,
and see if the others can reconstruct the source of the data. A
variation is that the subjects of the test can specify the type of
drives and we can see what problems/solutions various experimental
designs present.
Experiment design is a crucial part of system identification.

JCH
*http://home.arcor.de/janch/janch/_control/20091117-mimo-system/
Perhaps my spam filters are blocking but all I see is a complicated
block diagram.
I hope that you don't mind if I disagree with your flat statement.
Among other things the heat equation is quite explicit and succinct;
but the Laplace transform (or any other form of solution) has an
infinite number of poles/states. Having a good differential equation
form doesn't guarantee simplicity. In other cases, i.e. distributed
systems, the situation can also lead to complications.
In any case actually doing some test cases would be more interesting
than abstract talking.

Ray


pnachtwey

unread,
Nov 16, 2009, 10:58:34 AM11/16/09
to
As I said above, the quality of the data is very important. The
initial state must be known and usually that is to ensure the system
is at a steady state where the derivatives are 0. It is too bad the
data got truncated too. I work with motion control systems. It is
easy to make sure the system is stopped before getting into excitation
procedure. If the ambient temperature or C changes during the test
then that must be recorded too. Then it isn't a parameter to be
determined.

Are the time periods seconds? If so then time constants are long.
The long time constants make it hard to find the difference between
small changes in estimated time constants. The optimizing routines
calculate a sum of squared error. I like to think of the SSE as the
elevation in a multiple dimension terrain and the optimizer is seeking
the valley floor or pit. The slope of the SSE is checked in all
dimensions, one for each parameter being optimized. If the slopes
are very flat it is hard for the optimizing routine to get to a final
best position. The truncating data will make this almost impossible
because small difference make a big difference on a flat 'desert'

We have a motor that we put a relatively heavy disk on. The weight
increases the time constant to about 1 minutes which is forever in a
motion control system. This system doesn't ID with a consistent set
of numbers but all seem to work. It appear that there are some steep
hills but once the SSE gets off the steep hills the valley is more
like a flat 'desert'. The long time constants do make finding the
lowest part in the 'desert' difficult. However, any solution on the
desert floor seems to work equally well as the elevations are all the
same but the question is will they work equally well under all motion
conditions. This is odd but true. We have put a lot of effort in to
exciting the motor so that the SSE 'desert floor' has more of a slope.

Motors with short time constants are ID with more consistent
parameters.

Ray, you had the cards stacked against you, but your main problem was
with the data, not optim() or lsqrsolve(). The superimposed method
would work. The MIMO problem was not the biggest problem you had.

Hopefully this explains why auto tuning sometimes doesn't always
work.

Peter Nachtwey

JCH

unread,
Nov 17, 2009, 3:30:19 AM11/17/09
to

"RRogers" <rero...@plaidheron.com> schrieb im Newsbeitrag

[...]


>
> *http://home.arcor.de/janch/janch/_control/20091117-mimo-system/
> Perhaps my spam filters are blocking but all I see is a complicated
> block diagram.


I haven't sent more. Have again a close look to:

* http://home.arcor.de/janch/janch/_control/20091117-mimo-system/

This block diagram shows you eliminating coupling superposition. The MIMO
system 'should act' as if there were just two separate control loops not
interconnecting. Each of them can be optimized separately if the process
transfer functions are known. Therefore you MUST have a good process
identification.

Controller y1 influences process value x2 and vice versa!

Be aware that this is just an EXAMPLE for 2 input and 2 output signals. We
have to find 4 differential equations for this EXAMPLE using process
identification methods.

Physical derivation is difficult. Therefore measured (true) values should be
used for finding the differential equations.


--
Regards JCH


RRogers

unread,
Nov 17, 2009, 10:04:22 AM11/17/09
to

I am sure you don't want to get into this area but....
For a multi-state SISO system with poles on the negative real axis
there is a mathematical theory that is guaranteed to produce a
convergent series of models; I think it was extended to complex poles
as well. It is based up Muntz polynomials being converted to an
orthogonal polynomial series; really posynomials.
Good news and bad news from this theory. Bad news: almost any
sequence of poles will work. Good news: almost any sequence of poles
will work.
That means that you can guide the process and make mistakes; the
mistakes will be washed out later. But the result, while accurate,
will not necessarily be physically meaningful. It also implicitly
means that there are an infinite number of perfectly accurate models;
which might throw optimization procedures off a little. To summarize
you can estimate the two dominant poles, then by subtracting them out
of the data in a precise manner find another pole, use that in a
precise manner, and so on; the "precise manner" is the generation of
the orthogonal polynomial sequence. The most obvious application
would be to use the impulse response; but I think I see how to extend
it to arbitrary inputs.
I have been thinking about applying it to the standard sys-id process
to make successive approximations converge; get better and better as
the order increases. I haven't actually managed to get insight on
how to do this for MIMO approximations.

Enough said
Ray

pnachtwey

unread,
Nov 17, 2009, 7:11:57 PM11/17/09
to

They are all independent but the results at the temperature sensors
will be from the sum of the 3 heaters. This should hold true unless
there is something that I don't know where superimposition doesn't
apply. The states for a system of SOPDT equations would simply have
the temperature and the temperature rate at each of the 3 temperature
sensors. I don't see how the temperature at one sensor will affect
another temperature since the temperature sensors are not heat
sources.

Peter Nachtwey

RRogers

unread,
Nov 18, 2009, 12:01:27 PM11/18/09
to
clip....

> They are all independent but the results at the temperature sensors
> will be from the sum of the 3 heaters.  This should hold true unless
> there is something that I don't know where superimposition doesn't
> apply.  The states for a system of SOPDT equations would simply have
> the temperature and the temperature rate at each of the 3 temperature
> sensors.  I don't see how the temperature at one sensor will affect
> another temperature since the temperature sensors are not heat
> sources.
>
> Peter Nachtwey

Your right; and I see how the dependent/independent has to be handled
for state-space representation. A state space reduction to a non-
singular matrix is required to make A nonsingular in the ABCD
representation.

> > t1'=A1*t1+B11*u1(t-dt11)+B12*u2(t-dt12)+B13*u3(t-dt13)+C
> > t2'=A2*t2+B21*u1(t-dt21)+B22*u2(t-dt22)+B23*u3(t-dt23)+C
> > t3'=A3*t3+B31*u1(t-dt31)+B32*u2(t-dt32)+B33*u3(t-dt33)+C

As a note:
In more complicated systems you need additional terms on the left.
m1*t1''' + n1*t1'' + p1*t1'
m2*t2''' + n2*t2'' + p2*t2'

Ray


Phil Hobbs

unread,
Jan 9, 2020, 10:38:25 AM1/9/20
to
On 2009-11-06 00:55, Tim Wescott wrote:
> On Thu, 05 Nov 2009 18:37:21 +0100, Frank W. wrote:
>
>> Is there a control engineering expert here?
>>
>> I could us a bit of help on how to implement a PID autotune function for
>> a heating application (a small boiler).
>>
>> My current PID autotune function produces no reliable results. I use the
>> relay feedback method (Åstrom and Hägglund) but it seems that the
>> ultimate period Tu - one of the paramters determined with the relay
>> feedback test - is directly correlated to the relay output step u. u is
>> an arbitrary value which makes Tu an arbitary value. Since Tu is
>> required to compute Ti and Td (e.g. Ti = 0.5 Tu = Ziegler-Nichols),
>> autotune is not possible.
>
> Are you sure you haven't just stopped too soon in turning the output step
> down? I see where you're having trouble with this; is there some reason
> you can't turn u down even further? Things may start looking better when
> you get it down to where the temperature curve is roughly symmetrical
> around the average.
>
>> I believe this is because the machine heats very quickly (1500W boiler)
>> and cools down very slowly, so the process value isn't a sinusoid. Here
>> is an example graphics, green line is setpoint, red plot is process
>> value and the grey vertical bars represent heat output u:
>>
>> http://img196.imageshack.us/img196/7409/examplejr.jpg
>>
>> Regardless of u, the temperature always dips the same, small amount
>> below setpoint - because the machine cools slowly, the temperature can
>> not fall far below setpoint before it's reacting to the heat. It then
>> shoots up by an amount that is proportional to u. As a result, the plot
>> resembles mountains above setpoint. If u is large, temperature will
>> shoot high above the setpoint and take very long to cool down. Big
>> mountains, so Tu gets large, Tu ~ u = my problem.
>>
>> Since all PID temperature controllers have Autotune, there must be a
>> solution for this problem. Any ideas?
>
> Why do you want to use autotune? If this is a product that you're
> working on, and if the boiler design is the same for all of the parts
> that you ship, then you should just need one tuning.
>
> Autotune is for when you're selling a shrink-wrapped controller that has
> to work for anything -- and autotune is often considered to be a good way
> to get the tuning in the ballpark so that a human can get involved and
> actually make it work right.
>
> Unless you need to adapt to a variety of different, unexpected boiler
> combinations -- that don't change as the boiler is operating -- then
> there's not much value in autotune, IMHO.
>

For simple temperature control applications such as using small TE
coolers to stabilize diode lasers and optical detectors, we use a very
simple autotuning algorithm that works really well.

We hit the plant with a step function, and fit the T(t) curve to a plant
model consisting of a time delay tau followed by an integrator. That
has a simple transform proportional to exp(-j omega tau)/(j omega), so
we compute P and I to get about a 65-degree phase margin, and wind up
with a nice-looking transient response and decent bandwidth. This is
done at test time using BIST functions in firmware.

The main thing that improves the BW is to put an 0603 thermistor on the
bottom of the cold plate circuit board, right next to the TEC and
connected to the top plate of the TEC by a ground pour. That gets the
response down into 100-ms territory. (It's hard to measure at that
rate, because it's faster than the TEC itself.)

As I say, that's a much simpler plant than a domestic heating system,
but it works great.

Cheers

Phil Hobbs

--
Dr Philip C D Hobbs
Principal Consultant
ElectroOptical Innovations LLC / Hobbs ElectroOptics
Optics, Electro-optics, Photonics, Analog Electronics
Briarcliff Manor NY 10510

http://electrooptical.net
http://hobbs-eo.com

Phil Hobbs

unread,
Jan 9, 2020, 10:44:25 AM1/9/20
to
Doh, just noticed the time stamp. Weirdly this thread came up as unread
in Thunderbird this morning.
0 new messages