Menu
▾
▴
[Octave-cvsupdate] SF.net SVN: octave:[10691] trunk/octave-forge/extra/control-devel
|
From: <par...@us...> - 2012-06-25 14:56:59
|
Revision: 10691
http://octave.svn.sourceforge.net/octave/?rev=10691&view=rev
Author: paramaniac
Date: 2012-06-25 14:56:49 +0000 (Mon, 25 Jun 2012)
Log Message:
-----------
control-devel: use kalman predictor for simulation
Modified Paths:
--------------
trunk/octave-forge/extra/control-devel/devel/Destillation.m
trunk/octave-forge/extra/control-devel/devel/Evaporator.m
trunk/octave-forge/extra/control-devel/devel/GlassFurnace.m
trunk/octave-forge/extra/control-devel/devel/pH.m
trunk/octave-forge/extra/control-devel/inst/__slicot_identification__.m
trunk/octave-forge/extra/control-devel/inst/moen4.m
Added Paths:
-----------
trunk/octave-forge/extra/control-devel/devel/HeatingSystemKP.m
trunk/octave-forge/extra/control-devel/devel/PowerPlantKP.m
Modified: trunk/octave-forge/extra/control-devel/devel/Destillation.m
===================================================================
--- trunk/octave-forge/extra/control-devel/devel/Destillation.m 2012-06-25 12:57:13 UTC (rev 10690)
+++ trunk/octave-forge/extra/control-devel/devel/Destillation.m 2012-06-25 14:56:49 UTC (rev 10691)
@@ -66,13 +66,16 @@
Y_dest_n20=Y(:,7:9);
Y_dest_n30=Y(:,10:12);
+Y = {Y_dest; Y_dest_n10; Y_dest_n20; Y_dest_n30};
+U = {U_dest; U_dest_n10; U_dest_n20; U_dest_n30};
-dat = iddata (Y_dest, U_dest)
+dat = iddata (Y, U)
-[sys, x0] = moen4 (dat, 's', 5, 'n', 4) % s=5, n=4
+[sys, x0] = moen4 (dat, 's', 5, 'n', 4, 'noise', 'k') % s=5, n=4
+x0=x0{1};
-[y, t] = lsim (sys, U_dest, [], x0);
+[y, t] = lsim (sys, [U_dest, Y_dest], [], x0);
%[y, t] = lsim (sys, U_dest);
err = norm (Y_dest - y, 1) / norm (Y_dest, 1)
Modified: trunk/octave-forge/extra/control-devel/devel/Evaporator.m
===================================================================
--- trunk/octave-forge/extra/control-devel/devel/Evaporator.m 2012-06-25 12:57:13 UTC (rev 10690)
+++ trunk/octave-forge/extra/control-devel/devel/Evaporator.m 2012-06-25 14:56:49 UTC (rev 10691)
@@ -51,10 +51,10 @@
dat = iddata (Y, U)
-[sys, x0] = moen4 (dat, 's', 10, 'n', 4) % s=10, n=4
+[sys, x0] = moen4 (dat, 's', 10, 'n', 4, 'noise', 'k') % s=10, n=4
-[y, t] = lsim (sys, U, [], x0);
+[y, t] = lsim (sys, [U, Y], [], x0);
err = norm (Y - y, 1) / norm (Y, 1)
Modified: trunk/octave-forge/extra/control-devel/devel/GlassFurnace.m
===================================================================
--- trunk/octave-forge/extra/control-devel/devel/GlassFurnace.m 2012-06-25 12:57:13 UTC (rev 10690)
+++ trunk/octave-forge/extra/control-devel/devel/GlassFurnace.m 2012-06-25 14:56:49 UTC (rev 10691)
@@ -48,10 +48,10 @@
dat = iddata (Y, U)
-[sys, x0, info] = moen4 (dat, 's', 10, 'n', 5) % s=10, n=5
+[sys, x0, info] = moen4 (dat, 's', 10, 'n', 5, 'noise', 'k') % s=10, n=5
-[y, t] = lsim (sys, U, [], x0);
+[y, t] = lsim (sys, [U, Y], [], x0);
err = norm (Y - y, 1) / norm (Y, 1)
Added: trunk/octave-forge/extra/control-devel/devel/HeatingSystemKP.m
===================================================================
--- trunk/octave-forge/extra/control-devel/devel/HeatingSystemKP.m (rev 0)
+++ trunk/octave-forge/extra/control-devel/devel/HeatingSystemKP.m 2012-06-25 14:56:49 UTC (rev 10691)
@@ -0,0 +1,103 @@
+%{
+1. Contributed by:
+
+ Roy Smith
+ Dept. of Electrical & Computer Engineering
+ University of California,
+ Santa Barbara, CA 93106
+ U.S.A.
+ ro...@ec...
+
+2. Process/Description:
+
+ The experiment is a simple SISO heating system.
+ The input drives a 300 Watt Halogen lamp, suspended
+ several inches above a thin steel plate. The output
+ is a thermocouple measurement taken from the back of
+ the plate.
+
+3. Sampling interval:
+
+ 2.0 seconds
+
+4. Number of samples
+
+ 801
+
+5. Inputs:
+
+ u: input drive voltage
+ ...
+6. Outputs:
+
+ y: temperature (deg. C)
+ ...
+7. References:
+
+ The use of this experiment and data for robust
+ control model validation is described in:
+
+ "Sampled Data Model Validation: an Algorithm and
+ Experimental Application," Geir Dullerud & Roy Smith,
+ International Journal of Robust and Nonlinear Control,
+ Vol. 6, No. 9/10, pp. 1065-1078, 1996.
+
+8. Known properties/peculiarities
+
+ The data (and nominal model) is the above paper have the
+ output expressed in 10's deg. C. This has been rescaled
+ to the original units of deg. C. in the DaISy data set.
+ There is also a -1 volt offset in u in the data shown plotted
+ in the original paper. This has been removed in the
+ DaISy dataset.
+
+ The data shows evidence of discrepancies. One of the
+ issues studied in the above paper is the size of these
+ discrepancies - measured in this case in terms of the norm
+ of the smallest perturbation required to account for the
+ difference between the nominal model and the data.
+
+ The steady state input (prior to the start of the experiment)
+ is u = 6.0 Volts.
+
+%}
+
+clear all, close all, clc
+
+load heating_system.dat
+U=heating_system(:,2);
+Y=heating_system(:,3);
+
+
+dat = iddata (Y, U, 2.0, 'inname', 'input drive voltage', \
+ 'inunit', 'Volt', \
+ 'outname', 'temperature', \
+ 'outunit', '°C')
+
+% s=15, n=7
+%[sys1, x0] = moen4 (dat, 's', 15, 'n', 7)
+[sys1, x0] = moen4 (dat, 's', 15, 'n', 7, 'noise', 'k')
+
+%sys2 = arx (dat, 7, 7) % normally na = nb
+[sys2, x02] = arx (dat, 7);
+
+%[y1, t1] = lsim (sys1, U, [], x0);
+[y1, t1] = lsim (sys1, [U, Y], [], x0);
+
+%[y2, t] = lsim (sys2(:, 1), U);
+[y2, t] = lsim (sys2, U, [], x02);
+
+
+err1 = norm (Y - y1, 1) / norm (Y, 1)
+err2 = norm (Y - y2, 1) / norm (Y, 1)
+
+figure (1)
+plot (t, Y, t, y1, t, y2)
+title ('DaISy: Heating System [99-001]')
+xlim ([t(1), t(end)])
+xlabel ('Time [s]')
+ylabel ('Temperature [Degree Celsius]')
+legend ('measured', 'simulated subspace', 'simulated ARX', 'location', 'southeast')
+
+
+
Added: trunk/octave-forge/extra/control-devel/devel/PowerPlantKP.m
===================================================================
--- trunk/octave-forge/extra/control-devel/devel/PowerPlantKP.m (rev 0)
+++ trunk/octave-forge/extra/control-devel/devel/PowerPlantKP.m 2012-06-25 14:56:49 UTC (rev 10691)
@@ -0,0 +1,84 @@
+%{
+This file describes the data in the powerplant.dat file.
+1. Contributed by:
+ Peter Van Overschee
+ K.U.Leuven - ESAT - SISTA
+ K. Mercierlaan 94
+ 3001 Heverlee
+ Pet...@es...
+2. Process/Description:
+ data of a power plant (Pont-sur-Sambre (France)) of 120 MW
+3. Sampling time
+ 1228.8 sec
+4. Number of samples:
+ 200 samples
+5. Inputs:
+ 1. gas flow
+ 2. turbine valves opening
+ 3. super heater spray flow
+ 4. gas dampers
+ 5. air flow
+6. Outputs:
+ 1. steam pressure
+ 2. main stem temperature
+ 3. reheat steam temperature
+7. References:
+ a. R.P. Guidorzi, P. Rossi, Identification of a power plant from normal
+ operating records. Automatic control theory and applications (Canada,
+ Vol 2, pp 63-67, sept 1974.
+ b. Moonen M., De Moor B., Vandenberghe L., Vandewalle J., On- and
+ off-line identification of linear state-space models, International
+ Journal of Control, Vol. 49, Jan. 1989, pp.219-232
+8. Known properties/peculiarities
+
+9. Some MATLAB-code to retrieve the data
+ !gunzip powerplant.dat.Z
+ load powerplant.dat
+ U=powerplant(:,1:5);
+ Y=powerplant(:,6:8);
+ Yr=powerplant(:,9:11);
+
+%}
+
+clear all, close all, clc
+
+% NB: the code from DaISy is wrong:
+% powerplant(:,1) is just the sample number
+% therefore increase indices by one
+% it took me weeks to find that silly mistake ...
+load powerplant.dat
+U=powerplant(:,2:6);
+Y=powerplant(:,7:9);
+Yr=powerplant(:,10:12);
+
+inname = {'gas flow',
+ 'turbine valves opening',
+ 'super heater spray flow',
+ 'gas dampers',
+ 'air flow'};
+
+outname = {'steam pressure',
+ 'main steam temperature',
+ 'reheat steam temperature'};
+
+tsam = 1228.8;
+
+dat = iddata (Y, U, tsam, 'outname', outname, 'inname', inname)
+
+[sys, x0] = moen4 (dat, 's', 10, 'n', 8, 'noise', 'k') % s=10, n=8
+
+[y, t] = lsim (sys, [U, Y], [], x0);
+
+err = norm (Y - y, 1) / norm (Y, 1)
+
+figure (1)
+p = columns (Y);
+for k = 1 : p
+ subplot (3, 1, k)
+ plot (t, Y(:,k), 'b', t, y(:,k), 'r')
+endfor
+%title ('DaISy: Power Plant')
+%legend ('y measured', 'y simulated', 'location', 'southeast')
+
+st = isstable (sys)
+
Modified: trunk/octave-forge/extra/control-devel/devel/pH.m
===================================================================
--- trunk/octave-forge/extra/control-devel/devel/pH.m 2012-06-25 12:57:13 UTC (rev 10690)
+++ trunk/octave-forge/extra/control-devel/devel/pH.m 2012-06-25 14:56:49 UTC (rev 10691)
@@ -51,10 +51,10 @@
dat = iddata (Y, U)
-[sys, x0] = moen4 (dat, 's', 15, 'n', 6) % s=15, n=6
+[sys, x0] = moen4 (dat, 's', 15, 'n', 6, 'noise', 'k') % s=15, n=6
-[y, t] = lsim (sys, U, [], x0);
+[y, t] = lsim (sys, [U, Y], [], x0);
err = norm (Y - y, 1) / norm (Y, 1)
st = isstable (sys)
Modified: trunk/octave-forge/extra/control-devel/inst/__slicot_identification__.m
===================================================================
--- trunk/octave-forge/extra/control-devel/inst/__slicot_identification__.m 2012-06-25 12:57:13 UTC (rev 10690)
+++ trunk/octave-forge/extra/control-devel/inst/__slicot_identification__.m 2012-06-25 14:56:49 UTC (rev 10691)
@@ -154,6 +154,11 @@
in_v = __labels__ (outname, "y");
in_v = cellfun (@(x) ["v@", x], in_v, "uniformoutput", false);
inname = [in_u; in_v];
+ elseif (strncmpi (noise, "k", 1)) # Kalman predictor
+ sys = ss ([a-k*c], [b-k*d, k], c, [d, zeros(p)], tsam);
+ in_u = __labels__ (inname, "u");
+ in_y = __labels__ (outname, "y");
+ inname = [in_u; in_y];
else # no error inputs, default
sys = ss (a, b, c, d, tsam);
endif
Modified: trunk/octave-forge/extra/control-devel/inst/moen4.m
===================================================================
--- trunk/octave-forge/extra/control-devel/inst/moen4.m 2012-06-25 12:57:13 UTC (rev 10690)
+++ trunk/octave-forge/extra/control-devel/inst/moen4.m 2012-06-25 14:56:49 UTC (rev 10691)
@@ -71,6 +71,7 @@
## @table @var
## @item 'n'
## The desired order of the resulting state-space system @var{sys}.
+## @var{s} > @var{n} > 0.
##
## @item 's'
## The number of block rows @var{s} in the input and output
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
