Overview
Comment:Updated to use qp_mu (Adrian Wills/Will Heath)
Downloads: Tarball | ZIP archive | SQL archive
Timelines: family | ancestors | descendants | both | origin/master | trunk
Files: files | file ages | folders
SHA3-256: 6ec2f31f7fa701e875a1fc1e2ef8f1e4a0ae23a01898974cabb332520e80e864
User & Date: gawthrop@users.sourceforge.net on 2002-08-26 10:12:55
Other Links: branch diff | manifest | tags
Context
2002-08-27
10:48:29
Corrected documentation check-in: f60511b828 user: gawthrop@users.sourceforge.net tags: origin/master, trunk
2002-08-26
10:12:55
Updated to use qp_mu (Adrian Wills/Will Heath) check-in: 6ec2f31f7f user: gawthrop@users.sourceforge.net tags: origin/master, trunk
2002-08-23
09:13:58
Sensitivity version check-in: c578ade4ca user: gawthrop@users.sourceforge.net tags: origin/master, trunk
Changes

Modified mttroot/mtt/lib/control/PPP/ppp_qp.m from [e94ddba2d2] to [519bc8a390].

1
2
3
4
5
6
7
8
9

10
11
12
13
14
15
16
17
18
19
20





21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36


37
38
39
40
41
42
43
44
45
46
47
48
49
50
function [u,U,J] = ppp_qp (x,W,J_uu,J_ux,J_uw,Us0,Gamma,gamma)

  ## usage:  [u,U] = ppp_qp (x,W,J_uu,J_ux,J_uw,Gamma,gamma)
  ## INPUTS:
  ##      x: system state    
  ##      W: Setpoint vector
  ##      J_uu,J_ux,J_uw: Cost derivatives (see ppp_lin)
  ##      Us0: value of U* at tau=0 (see ppp_lin)
  ##      Gamma, gamma: U constrained by Gamma*U <= gamma 

  ## Outputs:
  ##      u: control signal
  ##      U: control weight vector
  ##
  ## Predictive pole-placement of linear systems using quadratic programming
  ## Use ppp_input_constraint and ppp_output_constraint to generate Gamma and gamma
  ## Use ppp_lin to generate J_uu,J_ux,J_uw
  ## Use ppp_cost to evaluate resultant cost function

  ## Copyright (C) 1999 by Peter J. Gawthrop
  ## 	$Id$	






  ## Check the sizes
  n_x = length(x);

  [n_U,m_U] = size(J_uu);
  if n_U != m_U
    error("J_uu must be square");
  endif

  [n,m] = size(J_ux);
  if (n != n_U)||(m != n_x)
    error("J_ux should be %ix%i not %ix%i",n_U,n_x,n,m);
  endif


  if length(gamma)>0		# Constraints exist: do the QP algorithm 


    U = qp(J_uu,(J_ux*x - J_uw*W),Gamma,gamma); # QP solution for weights U
    #U = pd_lcp04(J_uu,(J_ux*x - J_uw*W),Gamma,gamma); # QP solution for weights U
    u = Us0*U;			# Control signal
  else			# Do the unconstrained solution
    ## Compute the open-loop gains
    K_w = J_uu\J_uw;
    K_x = J_uu\J_ux;

    ## Closed-loop control
    U = K_w*W - K_x*x;		# Basis functions weights - U(t)
    u = Us0*U;			# Control u(t)
  endif

endfunction
|








>











>
>
>
>
>















|
>
>
|
|












1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
function [u,U,iterations] = ppp_qp (x,W,J_uu,J_ux,J_uw,Us0,Gamma,gamma,mu)

  ## usage:  [u,U] = ppp_qp (x,W,J_uu,J_ux,J_uw,Gamma,gamma)
  ## INPUTS:
  ##      x: system state    
  ##      W: Setpoint vector
  ##      J_uu,J_ux,J_uw: Cost derivatives (see ppp_lin)
  ##      Us0: value of U* at tau=0 (see ppp_lin)
  ##      Gamma, gamma: U constrained by Gamma*U <= gamma 
  ##      mu  Parameter of qp_mu
  ## Outputs:
  ##      u: control signal
  ##      U: control weight vector
  ##
  ## Predictive pole-placement of linear systems using quadratic programming
  ## Use ppp_input_constraint and ppp_output_constraint to generate Gamma and gamma
  ## Use ppp_lin to generate J_uu,J_ux,J_uw
  ## Use ppp_cost to evaluate resultant cost function

  ## Copyright (C) 1999 by Peter J. Gawthrop
  ## 	$Id$	

  if nargin<9
    mu = 0
  endif


  ## Check the sizes
  n_x = length(x);

  [n_U,m_U] = size(J_uu);
  if n_U != m_U
    error("J_uu must be square");
  endif

  [n,m] = size(J_ux);
  if (n != n_U)||(m != n_x)
    error("J_ux should be %ix%i not %ix%i",n_U,n_x,n,m);
  endif


  if length(gamma)>0		# Constraints exist: do the QP algorithm
    [U,iterations] = qp_mu(J_uu,(J_ux*x - J_uw*W),Gamma,gamma,mu); # QP solution for weights U	

    ##U = qp(J_uu,(J_ux*x - J_uw*W),Gamma,gamma); # QP solution for weights U
    ##U = pd_lcp04(J_uu,(J_ux*x - J_uw*W),Gamma,gamma); # QP solution for weights U
    u = Us0*U;			# Control signal
  else			# Do the unconstrained solution
    ## Compute the open-loop gains
    K_w = J_uu\J_uw;
    K_x = J_uu\J_ux;

    ## Closed-loop control
    U = K_w*W - K_x*x;		# Basis functions weights - U(t)
    u = Us0*U;			# Control u(t)
  endif

endfunction

Modified mttroot/mtt/lib/control/PPP/ppp_qp_sim.m from [58fd63cdfa] to [442bf78425].

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16




17
18

19
20
21
22
23
24
25
function [T,y,u,J] = ppp_qp_sim (A,B,C,D,A_u,A_w,t,Q, Tau_u,Min_u,Max_u,Order_u, Tau_y,Min_y,Max_y,Order_y, W,x_0,Delta_ol,movie)

  ## usage: [T,y,u,J] = ppp_qp_sim (A,B,C,D,A_u,A_w,t,Q, Tau_u,Min_u,Max_u,Order_u, Tau_y,Min_y,Max_y,Order_y, W,x_0,movie)
  ## Needs documentation - see ppp_ex11 for example of use.
  ## OUTPUTS
  ## T: Time vector
  ## y,u,J output, input and cost

  ## Copyright (C) 1999 by Peter J. Gawthrop
  ## 	$Id$	
  
  if nargin<19			# No intermittent control
    Delta_ol = 0;
  endif

  if nargin<20			# No movie




    movie = 0;
  endif


  ## Check some sizes
  [n_x,n_u,n_y] = abcddim(A,B,C,D);

  [n_x0,m_x0] = size(x_0);
  if (n_x0 != n_x)||(m_x0 != 1)
    error(sprintf("Initial state x_0 must be %ix1 not %ix%i",n_x,n_x0,m_x0));
|















>
>
>
>


>







1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
function [T,y,u,Iterations] = ppp_qp_sim (A,B,C,D,A_u,A_w,t,Q, Tau_u,Min_u,Max_u,Order_u, Tau_y,Min_y,Max_y,Order_y, W,x_0,Delta_ol,mu,movie)

  ## usage: [T,y,u,J] = ppp_qp_sim (A,B,C,D,A_u,A_w,t,Q, Tau_u,Min_u,Max_u,Order_u, Tau_y,Min_y,Max_y,Order_y, W,x_0,movie)
  ## Needs documentation - see ppp_ex11 for example of use.
  ## OUTPUTS
  ## T: Time vector
  ## y,u,J output, input and cost

  ## Copyright (C) 1999 by Peter J. Gawthrop
  ## 	$Id$	
  
  if nargin<19			# No intermittent control
    Delta_ol = 0;
  endif

  if nargin<20			# No movie
    mu = 0;
  endif

  if nargin<21			# No movie
    movie = 0;
  endif


  ## Check some sizes
  [n_x,n_u,n_y] = abcddim(A,B,C,D);

  [n_x0,m_x0] = size(x_0);
  if (n_x0 != n_x)||(m_x0 != 1)
    error(sprintf("Initial state x_0 must be %ix1 not %ix%i",n_x,n_x0,m_x0));
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
  disp("Designing controller");
  [k_x,k_w,K_x,K_w,Us0,J_uu,J_ux,J_uw,J_xx,J_xw,J_ww] = ppp_lin (A,B,C,D,A_u,A_w,t,Q);

  ## Set up various time vectors
  dt = t(2)-t(1);		# Time increment

  ## Make sure Delta_ol is multiple of dt
  Delta_ol = floor(Delta_ol/dt)*dt

  if Delta_ol>0			# Intermittent control
    T_ol = 0:dt:Delta_ol-dt;	# Create the open-loop time vector
  else
    T_ol = 0;
    Delta_ol = dt;
  endif







|







50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
  disp("Designing controller");
  [k_x,k_w,K_x,K_w,Us0,J_uu,J_ux,J_uw,J_xx,J_xw,J_ww] = ppp_lin (A,B,C,D,A_u,A_w,t,Q);

  ## Set up various time vectors
  dt = t(2)-t(1);		# Time increment

  ## Make sure Delta_ol is multiple of dt
  Delta_ol = floor(Delta_ol/dt)*dt;

  if Delta_ol>0			# Intermittent control
    T_ol = 0:dt:Delta_ol-dt;	# Create the open-loop time vector
  else
    T_ol = 0;
    Delta_ol = dt;
  endif
71
72
73
74
75
76
77

78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111

112
113
114
115
116
117
118
119
120
121
122
123
124
125

126
127
128
129
130
131
132
133
134
135
136
137
  [Ad, Bd] = sys2ss(dsys);

  x = x_0;			# Initialise state

  ## Initialise the saved variable arrays
  X = [];
  u = [];

  du = [];
  J = [];
  tick= time;
  i = 0;
  disp("Simulating ...");
  for t=T_cl			# Outer loop at Delta_ol
    ##disp(sprintf("Time %g", t));
    ## Output constraints
    [Gamma_y,gamma_y] = ppp_output_constraint  (A,B,C,D,x,A_u,Tau_y,Min_y,Max_y,Order_y);
    
    ## Composite constraints 
    Gamma = [Gamma_u; Gamma_y];
    gamma = [gamma_u; gamma_y];
    
    ## Compute U(t)
    [uu U] = ppp_qp (x,W,J_uu,J_ux,J_uw,Us0,Gamma,gamma); # Compute U
 
    ## Compute the cost (not necessary but maybe interesting)
#    [J_t] = ppp_cost (U,x,W,J_uu,J_ux,J_uw,J_xx,J_xw,J_ww); # cost
#    J = [J J_t];

    ## Simulation loop
    i_ol = 0;
    for t_ol=T_ol		# Inner loop at dt

      ## Compute ol control
      i_ol = i_ol+1;
      range = (i_ol-1)*n_U + 1:i_ol*n_U; # Extract current U*
      ut = Ustar_ol(:,range)*U;	# Compute OL control (U* U)

      ## Simulate the system
      i = i+1;
      X = [X x];		# Save state
      u = [u ut];		# Save input

      x = Ad*x + Bd*ut;	# System

#       if movie			# Plot the moving horizon
# 	tau = T(1:n_T-i);	# Tau with moving horizon
# 	tauT = T(i+1:n_T);	# Tau with moving horizon + real time
# 	[ys,us,xs,xu,AA] = ppp_ystar (A,B,C,D,x,A_u,U,tau); # OL response
# 	plot(tauT,ys, tauT(1), ys(1), "*")
#       endif
    endfor
  endfor
  
  ## Save the last values
  X = [X x];		# Save state
  u = [u ut];		# Save input


  tock = time;
  Iterations = length(T_cl)
  Elapsed_Time = tock-tick
  y = C*X + D*u;		# System output

  T = 0:dt:t+Delta_ol;		# Overall time vector

endfunction










>















|


















>














>


<









76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135

136
137
138
139
140
141
142
143
144
  [Ad, Bd] = sys2ss(dsys);

  x = x_0;			# Initialise state

  ## Initialise the saved variable arrays
  X = [];
  u = [];
  Iterations = [];
  du = [];
  J = [];
  tick= time;
  i = 0;
  disp("Simulating ...");
  for t=T_cl			# Outer loop at Delta_ol
    ##disp(sprintf("Time %g", t));
    ## Output constraints
    [Gamma_y,gamma_y] = ppp_output_constraint  (A,B,C,D,x,A_u,Tau_y,Min_y,Max_y,Order_y);
    
    ## Composite constraints 
    Gamma = [Gamma_u; Gamma_y];
    gamma = [gamma_u; gamma_y];
    
    ## Compute U(t)
    [uu, U, iterations] = ppp_qp (x,W,J_uu,J_ux,J_uw,Us0,Gamma,gamma,mu); # Compute U
 
    ## Compute the cost (not necessary but maybe interesting)
#    [J_t] = ppp_cost (U,x,W,J_uu,J_ux,J_uw,J_xx,J_xw,J_ww); # cost
#    J = [J J_t];

    ## Simulation loop
    i_ol = 0;
    for t_ol=T_ol		# Inner loop at dt

      ## Compute ol control
      i_ol = i_ol+1;
      range = (i_ol-1)*n_U + 1:i_ol*n_U; # Extract current U*
      ut = Ustar_ol(:,range)*U;	# Compute OL control (U* U)

      ## Simulate the system
      i = i+1;
      X = [X x];		# Save state
      u = [u ut];		# Save input
      Iterations = [Iterations iterations]; # Save iteration count
      x = Ad*x + Bd*ut;	# System

#       if movie			# Plot the moving horizon
# 	tau = T(1:n_T-i);	# Tau with moving horizon
# 	tauT = T(i+1:n_T);	# Tau with moving horizon + real time
# 	[ys,us,xs,xu,AA] = ppp_ystar (A,B,C,D,x,A_u,U,tau); # OL response
# 	plot(tauT,ys, tauT(1), ys(1), "*")
#       endif
    endfor
  endfor
  
  ## Save the last values
  X = [X x];		# Save state
  u = [u ut];		# Save input
  Iterations = [Iterations iterations]; # Save iteration count

  tock = time;

  Elapsed_Time = tock-tick
  y = C*X + D*u;		# System output

  T = 0:dt:t+Delta_ol;		# Overall time vector

endfunction




MTT: Model Transformation Tools
GitHub | SourceHut | Sourceforge | Fossil RSS ]