gpc_demo1.m

由于自校正控制器在实际中的广泛运用

function varargout = GPC_demo1(varargin)
% GPC_demo1 M-file for GPC_demo1.fig
%      GPC_demo1, by itself, creates alfa new GPC_demo1 or raises the existing
%      singleton*.
%
%      H = GPC_demo1 returns the handle to alfa new GPC_demo1 or the handle to
%      the existing singleton*.
%
%      GPC_demo1('CALLBACK',hObject,eventData,handles,...) calls the local
%      function named CALLBACK in GPC_demo1.M with the given input arguments.
%
%      GPC_demo1('Property','Value',...) creates alfa new GPC_demo1 or raises the
%      existing singleton*.  Starting from the left, property value pairs are
%      applied to the GUI before GPC_demo1_OpeningFunction gets called.  An
%      unrecognized property name or invalid value makes property application
%      stop.  All inputs are passed to GPC_demo1_OpeningFcn via varargin.
%
%      *See GUI Options on GUIDE's Tools menu.  Choose "GUI allows only one
%      instance to run (singleton)".
%
% See also: GUIDE, GUIDATA, GUIHANDLES

% Edit the above text to modify the response to help GPC_demo1

% Last Modified by GUIDE v2.5 17-Dec-2005 18:30:33

% Begin initialization code - DO NOT EDIT
gui_Singleton = 1;
gui_State = struct('gui_Name',       mfilename, ...
                   'gui_Singleton',  gui_Singleton, ...
                   'gui_OpeningFcn', @GPC_demo1_OpeningFcn, ...
                   'gui_OutputFcn',  @GPC_demo1_OutputFcn, ...
                   'gui_LayoutFcn',  [] , ...
                   'gui_Callback',   []);
if nargin & isstr(varargin{1})
    gui_State.gui_Callback = str2func(varargin{1});
end

if nargout
    [varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:});
else
    gui_mainfcn(gui_State, varargin{:});
end
% End initialization code - DO NOT EDIT

% --- Executes just before GPC_demo1 is made visible.
function GPC_demo1_OpeningFcn(hObject, eventdata, handles, varargin)
% This function has no output args, see OutputFcn.
% hObject    handle to figure
% eventdata  reserved - to be defined in alfa future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
% varargin   command line arguments to GPC_demo1 (see VARARGIN)

% Choose default command line output for GPC_demo1
handles.output = hObject;

% Update handles structure
guidata(hObject, handles);

% This sets up the initial plot - only do when we are invisible
% so window can get raised using GPC_demo1.
if strcmp(get(hObject,'Visible'),'off')
    axes(handles.axes1);
    cla;
    text(0.2,0.5,'输入各参数和初始值，按"Start"开始仿真 ');
    
    axes(handles.axes2);
    cla;
    text(0.2,0.5,'仿真结束后，按"End"退出程序');
end

% UIWAIT makes GPC_demo1 wait for user response (see UIRESUME)
% uiwait(handles.figure1);


% --- Outputs from this function are returned to the command line.
function varargout = GPC_demo1_OutputFcn(hObject, eventdata, handles)
% varargout  cell array for returning output args (see VARARGOUT);
% hObject    handle to figure
% eventdata  reserved - to be defined in alfa future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Get default command line output from handles structure
varargout{1} = handles.output;

% --- Executes on button press in pushbutton1.
function pushbutton1_Callback(hObject, eventdata, handles)
% hObject    handle to pushbutton1 (see GCBO)
% eventdata  reserved - to be defined in alfa future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)



% 运行步数 K
data = getappdata(gcbf, 'metricdata');
K=data.K;

% 控制参数 lamda ,alfa
data = getappdata(gcbf, 'metricdata');
r=data.r;

data = getappdata(gcbf, 'metricdata');
alfa=data.alfa;

% 参考轨迹yr的模型
popup_sel_index = get(handles.popupmenu3, 'Value');
switch popup_sel_index
    case 1                  %DC signal
        Yr=zeros((K+2),1);  
        for i=1:(K+2)
            Yr(i,1)=1;
        end
    case 2                  %Sw signal
        Yr=zeros((K+2),1);  
        for i=1:(K+2)
            if K<200
                Yr(i,1)=sawtooth(i/5); 
            else
                Yr(i,1)=sawtooth(i/25);
            end 
        end
    case 3                 %Cos signal
        Yr=zeros((K+2),1);  
        for i=1:(K+2)
            if K<200
               Yr(i,1)=cos(i/5);
            else
               Yr(i,1)=cos(i/25);  
            end   
            
        end
    case 4                 %Square signal
        Yr=zeros((K+2),1);  
        for i=1:(K+2)
            if K<200
               Yr(i,1)=square(i/4);
            else
               Yr(i,1)=square(i/40);  
            end   
        end
    case 5                 %Satrurate signal
        Yr=zeros((K+2),1);  
        for i=1:(K+2)
            if i<140
                Yr(i,1)=0.01*i;
            else
                Yr(i,1)=1.4;
            end             
        end
end

%输入在时间为0和－1时y和u的初始值
data = getappdata(gcbf, 'metricdata');
Y(1,1)=data.yfu1;

data = getappdata(gcbf, 'metricdata');
Y(2,1)=data.y0;

data = getappdata(gcbf, 'metricdata');
U(1,1)=data.ufu1;

data = getappdata(gcbf, 'metricdata');
U(2,1)=data.u0;


syms z

% 输入被控模型阶次
na=2;
nb=1;

% 输入被控模型的系数
a1=-1.5;
a2=0.7;
b0=1;
b1=0.5;

AA=zeros(1,na);
BB=zeros(1,nb);
AA(1,1)=a1*(-1);
AA(1,2)=a2*(-1);
BB(1,1)=b0;
BB(1,2)=b1;

N=4;
Nu=3;

E=zeros(N,N);
F=zeros(N,na+1);
H=zeros(N,1);
G=zeros(N,N);

A=1+a1*z^(-1)+a2*z^(-2);
B=b0+b1*z^(-1);

EE=1;
FF=z*(1-A*(1-z^(-1)))*z^na;

GG=EE*B*z;
GGxishu=sym2poly(GG);
G(1,1)=GGxishu(1,1);
H(1,1)=GGxishu(1,2);

FFxishu=sym2poly(FF);
E(1,1)=1;
for jj=1:(na+1)
    F(1,jj)=FFxishu(1,jj);
end
fj0=F(1,1);
ej=fj0;
EE=(EE+ej*z^(-1))*z;
FF=z*(FF/(z^na)-ej*A*(1-z^(-1)))*z^na ;

GG=(EE*B)*z;
GGxishu=sym2poly(GG);
for jj=1:2
    G(2,jj)=GGxishu(1,jj);
end
H(2,1)=GGxishu(1,3);

for j=2:N
    EExishu=sym2poly(EE);
    for jj=1:j
        E(j,jj)=EExishu(1,jj);
    end
    FFxishu=sym2poly(FF);
    for jj=1:(na+1)
        F(j,jj)=FFxishu(1,jj);
    end
    fj0=F(j,1);
    ej=fj0;
    EE=poly2sym(EExishu);
    EE=subs(EE,'x','z');
    EE=(EE/(z^(j-1))+ej*z^(-j))*(z^j);
    FF=poly2sym(FFxishu);
    FF=subs(FF,'x','z');
    FF=z*(FF/(z^na)-ej*A*(1-z^(-1)))*z^na;
    
    GG=EE/(z^j)*B*z^(j+1);
    GGxishu=sym2poly(GG);
    for jj=1:j
        G(j,jj)=GGxishu(1,jj);
    end
    
    if j<N
        H((j+1),1)=GGxishu(1,(j+2));
    end
end
G1=zeros(N,Nu);  % G 矩阵
for j=1:Nu
    for jj=j:N
        G1(jj,j)=G(N,jj-j+1);
    end
end


Y0=zeros(N,1);
Yd=zeros(N,1);
 
d=[1 zeros(1,Nu-1)];

for i=3:K+2
    aa=i; 
    bb=i;
    for a=1:na
        Y2(a,1)=Y(aa-1,1);
        aa=aa-1;
    end 
    for b=1:nb+1
        U2(b,1)=U(bb-1,1);
        bb=bb-1;
    end
    AAA=AA;
    AAA(:,1)=[]; 
         
    Y(i,1)=AA*Y2+BB*U2;
    Yd(1,1)=Y(i,1);
    for t=1:na+1
        y(t,1)=Y(i+1-t,1);
    end
    
    delta=U(i-1,1)-U(i-2,1);
    for t=1:N 
        Y0(t,1)=F(t,1)*Y(i,1)+F(t,2)*Y((i-1),1)+F(t,3)*Y((i-2),1)+H(t,1)*delta;
    end
    
    for j=2:N       
        Yd(j,1)=alfa.*Yd(j-1,1)+(1-alfa).*Yr(i,1);
    end
         
    deltaUA=G1'*G1+r.*eye(Nu);
    deltaUB=inv(deltaUA);
    deltaUC=deltaUB*G1';
    deltaU=deltaUC*[Yd-Y0];
    
    U(i,1)=U(i-1,1)+d*deltaU;
end
Y(1:1:2)=[]; 




% 描绘出所求得Y的图像
ii=1;

axes(handles.axes1);
cla;

axis([0 K -2.5 2.5]);
hold on;

xlabel('预测步长N'),ylabel('参考输入yr和模型输出y')
title('参考输入yr和模型输出y的轨迹')
    
    
axes(handles.axes2);
cla;
axis([0 K -1.5 1.5]);
hold on;

xlabel('预测步长N'),ylabel('控制输入u')
title('控制输入u的轨迹')

for i=1:K
    
    YR=zeros(i,1);
    for j=1:i
        YR(j,1)=Yr(j,1);
    end
    
    
    YY=zeros(i,1);
    for j=1:i
        YY(j,1)=Y(j,1);
    end
    axes(handles.axes1);
    cla;
    plot(YR,'LineWidth',1.3, 'Color','r'); 
    hold on;
    plot(YY,'LineWidth',1.3, 'Color',[0.439,0.529,0.129]); 
    hold on;

    
    UU=zeros(i,1);
    for j=1:i
        UU(j,1)=U(j,1);
    end
    axes(handles.axes2);
    cla;
    plot(UU,'LineWidth',1.3, 'Color','b'); 
    hold on;
    
     
    mmm(ii)=getframe;
    ii=ii+1;
    
    data = getappdata(gcbf, 'metricdata');
    yr = YR(i,1);
    set(handles.yr, 'String',yr);
    data = getappdata(gcbf, 'metricdata');
    y = YY(i,1);
    set(handles.y, 'String',y);
    data = getappdata(gcbf, 'metricdata');
    u = UU(i,1);
    set(handles.u, 'String',u);
    
end;
movie(mmm,1,1);
hold on;


  

% --------------------------------------------------------------------
function FileMenu_Callback(hObject, eventdata, handles)
% hObject    handle to FileMenu (see GCBO)
% eventdata  reserved - to be defined in alfa future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)


% --------------------------------------------------------------------
function OpenMenuItem_Callback(hObject, eventdata, handles)
% hObject    handle to OpenMenuItem (see GCBO)
% eventdata  reserved - to be defined in alfa future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
file = uigetfile('*.fig');
if ~isequal(file, 0)
    open(file);
end

% --------------------------------------------------------------------
function PrintMenuItem_Callback(hObject, eventdata, handles)
% hObject    handle to PrintMenuItem (see GCBO)
% eventdata  reserved - to be defined in alfa future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
printdlg(handles.figure1)

% --------------------------------------------------------------------
function CloseMenuItem_Callback(hObject, eventdata, handles)
% hObject    handle to CloseMenuItem (see GCBO)
% eventdata  reserved - to be defined in alfa future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
selection = questdlg(['Close ' get(handles.figure1,'Name') '?'],...
                     ['Close ' get(handles.figure1,'Name') '...'],...
                     'Yes','No','Yes');
if strcmp(selection,'No')
    return;
end

delete(handles.figure1)