### [lecture #1] 2020.9.3 outline of the lecture, review of classical and modern control theory (1/3) †

• review : stabilization of SISO unstable plant by classical and modern control theory
• transfer functions / differential equations
• poles / eigenvalues
• impulse response / initial value response
• ...
s = tf('s')
P = 1/(s-1)
pole(P)
impulse(P)
k = 2
Tyr = feedback(P*k, 1)
step(Tyr)
k = 10
Tyr = feedback(P*k, 1)
step(Tyr)
k = 0.5
Tyr = feedback(P*k, 1)
step(Tyr)

### [lecture #2] 2020.9.10 review of classical and modern control theory (2/3) with introduction of Matlab/Simulink †

1. minute paper
2. introduction of Matlab and Simulink text_fixed.pdf Basic usage of MATLAB and Simulink used for 情報処理演習及び考究II/Consideration and Practice of Information Processing II: Advanced Course of MATLAB
• interactive system (no compilation, no variable definition)
• m file
3. system representation: Transfer Function(TF) / State-Space Representation (SSR)
• example: mass-spring-damper system
• definition of SSR
• from SSR to TF
• from TF to SSR: controllable canonical form
4. open-loop characteristic
5. closed-loop stability
%-- 20/09/10 12:48 --%
a = 1
b = 2
a + b
ex0910_1
P
ctrlpref
ex0910_1
P
P.num
P.num{:}
P.den{:}
ex0910_2
roots(P.den{:})
roots(L.den{:})
ex0910_2
K
ex0910_2
• Q: マトラボの使いかたが分かるオススメの本やサイトが知りたいです
• A: この授業では基本的な使い方ができれば十分のため、まずは今日示したpdfファイルを見てもらえればと思います。

### [lecture #3] 2020.9.17 review of classical and modern control theory (3/3) †

1. LQR problem
%-- 20/09/17 12:42 --%
ex0917_1
A
B
Uc
det(Uc)
R
M
M = [2, 1; 1, 1]
eig(M)
ex0917_1
C
D
F
R
Q
P
A
eig(A)
eig(P)
R\B'
R\B'*P
F
J
J(end)

### [lecture #4] 2020.9.24 relation between LQR and H infinity control problem (1/2) †

• GOAL: to learn difference in concepts between LQR problem and H infinity control problem
1. a simple example relating LQR and H infinity control problems
• For given plant G $G = \left[\begin{array}{c|c:c} a & 1 & b \\ \hline \sqrt{q} & 0 & 0 \\ 0 & 0 & \sqrt{r} \\ \hdashline 1 & 0 & 0 \end{array} \right] = \left\{ \begin{array}{l} \dot x = ax + bu + w\\ z = \left[ \begin{array}{c} \sqrt{q} x \\ \sqrt{r} u \end{array}\right] \\ x = x \end{array}\right.$ with zero initial state value x(0) = 0, find a state-feedback controller $u = -f x$ such that \begin{eqnarray} (i) &&\quad \mbox{closed loop is stable} \\ (ii) &&\quad \mbox{minimize} \left\{\begin{array}{l} \| z \|_2 \mbox{ for } w(t) = \delta(t) \quad \mbox{(LQR)} \\ \| T_{zw} \|_\infty \mbox{($H_\infty$ control problem)}\end{array}\right. \end{eqnarray}
• comparison of norms in (ii) (for a = -1, b = 1, q = 1, r = 1) $\begin{array}{|c||c|c|}\hline & \mbox{LQR}: f=-1+\sqrt{2} & \quad \quad H_\infty: f=1\quad\quad \\ \hline\hline J=\|z\|_2^2 & & \\ \hline \|T_{zw}\|_\infty & & \\ \hline \end{array}$
2. an alternative description to LQR problem
1. J = (L2 norm of z)^2
2. impulse resp. with zero initial value = initial value resp. with zero disturbance
3. definition of H infinity norm (SISO)
s = tf('s');
G1 = 1/(s+1);
bode(G1);
norm(G1, 'inf')
G2 = 1/(s^2 + 0.1*s + 1);
bode(G2);
norm(G2, 'inf')
4. definition of H infinity norm (SIMO)
5. solve the problem by hand
6. solve the problem by tool(hinfsyn) ex0924_1.m
%-- 20/09/24 12:36 --%
s = tf('s');
G1 = 1/(s+1);
bode(G1);
norm(G1, 'inf')
norm(G1)
norm(G1, inf)
bode(G2)
G2 = 1/(s^2 + 0.1*s + 1);
bode(G2);
G2
ctrlpref
bode(G2);
norm(G2, inf)
• Q: 日本語での講義が良いです。
• A: 反対意見が無かったため、次回から日本語で講義を行います。

### [lecture #5] 2020.10.08 relation between LQR and H infinity control problem (2/2) †

1. complete the table in simple example
2. confirm the cost function J for both controllers by simulation mod1008.mdl
• block diagram in the simulink model
• how to approximate impulse disturbance with a step function
• impulse disturbance resp. with zero initial condition = initial condition resp. with zero disturbance
3. confirm the closed-loop H infinity norm for both controllers by simulation
• H infinity norm = L2 induced norm
• review: steady-state response; the worst-case disturbance w(t) which maximizes L2 norm of z(t) ?
• how to make the worst-case disturbance w(t)? w(t) for the simple example ?
4. general state-feedback case: hinf.pdf
%-- 20/10/08 12:36 --%
ex0924_1
mod1008
f
f = -1+sqrt(2)
x0
x0 = 0
h
h = 0.01
f
zz
zz(end)
sqrt(2)-1
h = 0.001
zz(end)
f
f = 1
zz(end)
h
x0
x0 = 1
zz(end)
x0
h
h = 0.00001
zz(end)
x0 = 0
zz(end)
format long e
zz(end)
x0
f
h
h = 10
zz(end)
sqrt(zz(end)/ww(end))
format short
sqrt(zz(end)/ww(end))
h
h = 100
sqrt(zz(end)/ww(end))
f
f = -1+sqrt(2)
sqrt(zz(end)/ww(end))

### [lecture #6] 2020.10.15 Mixed sensitivity problem 1/3 †

1. outline: map_v1.1_mixedsens1.pdf
• sensitivity function S and complementary sensitivity function T
2. H infinity control problem (general case)
• with generalized plant G
• including the state-feedback case
3. reference tracking problem
• how to translate the condition (ii) into one with H infinity norm ?
• corresponding generalized plant G ?
• introduction of weighting function for sensitivity function in (ii)
4. design example ex1015_1.m ex1015_2.m
5. the small gain theorem
• proof: Nyquist stability criterion
%-- 20/10/15 12:54 --%
ex1015_1
P
eig(P)
K
help step
ex1015_2
who
K_hinf
eig(K.a)
eig(K_hinf.a)
help hinfsyn
ex1015_2
P

### [lecture #7] 2020.10.22 Mixed sensitivity problem 2/3 †

1. outline: from point to set map_v1.1_mixedsens2.pdf
2. review: the small gain theorem ... robust stability = H infinity norm condition
3. normalized uncertainty Delta
4. uncertainty model
5. how to determine P0 and WT
6. robust stabilization problem and equivalent problem
%-- 20/10/22 12:55 --%
ex1022_1
ex1022_2
ex1022_3
mod1022
c
c = 0.8
c = 2

### [lecture #8] 2020.10.29 Mixed sensitivity problem 3/3 †

1. mixed sensitivity problem => (1) and (2) : proof
2. generalized plant for mixed senstivity problem
3. design example ex1029_1.m minimize gamma by hand
4. gamma iteration by bisection method ex1029_2.m
5. intro. to RP(problem of NP) ex1029_3.m
%-- 20/10/29 13:00 --%
ex1029_1
ex1029_2
ex1029_3

### [lecture #9] 2020.11.5 robust performance problem 1/3 †

1. review
• mixed sensitivity problem : N.P. but not R.P.
• robust performance problem (R.P.) c.f. the last whiteboard, but can not be solved by tool
• the small gain theorem
2. an equivalent robust stability (R.S.) problem to R.P.
• (i) introduction of a fictitious uncertainty Delta_p (for performance)
• (ii) for 2-by-2 uncertainty block Delta hat which includes Delta and Delta_p
3. definition of H infinity norm for general case (MIMO)
4. proof of ||Delta hat||_inf <= 1
5. design example: ex1105_1.m
• robust performance is achieved but large gap
• non structured uncertainty is considered ... the design problem is too conservative
%-- 20/11/05 13:59 --%
M = [1/sqrt(2); 1i; 1/sqrt(2), -1i]
M'
eig(M'*M)
svd(M)
M = [1/sqrt(2), 1i; 1/sqrt(2), -1i]
M'
M'*M
eig(M'*M)
svd(M)
max(svd(M))
ex1105_1

### [lecture #10] 2020.11.12 Robust performance problem (2/3) †

1. return of mini report #1
2. SVD: singular value decomposition
%-- 20/11/12 13:18 --%
M = [sqrt(2), -1i/sqrt(2); sqrt(2), 1i/sqrt(2)]
help svd
X
svd(M)
[U,Sigma,V] = svd(M);
Sigma
U
U'*U
V'*V
V*V'
ex1112_1
svd(M)
ex1112_1
help rand
ex1112_1
svd(M)

### [lecture #11] 2020.11.19 Robust performance problem (3/3) †

1. review
• H infinity norm (MIMO case)
• R.S. problems for structured and unstructured uncertainty
2. scaled H infinity control problem
3. relation between three problems
4. how to determine structure of scaling matrix
5. design example ex1119_1.m
ex1105_1
gam2 = gam_opt
K2 = K_opt;
ex1119_1
gam_opt
6. mini exam #1 (10 min.)
%-- 20/11/19 12:51 --%
ex1105_1

### [lecture #12] 2020.11.26 Robust performance problem (3/3) (cont.), Control system design for practical system (1/3) †

1. return of mini exam #1
2. review of scaling ex1126_1.m
3. mini report #2 report2.pdf
4. introduction of a practical system: active noise control in duct
%-- 20/11/26 13:09 --%
ex1105_1
gam2 = gam_opt
K2 = K_opt;
ex1119_1
gam_opt
gam2
ex1126_1
gam2
gam3
close all
clear all
clf
ex1105_1; %
gam2 = gam_opt
K2 = K_opt;
ex1119_1
gam3 = gam_opt
K3 = K_opt;
Delta_tilde = [0, 1/sqrt(2); 0, -1/sqrt(2)]; % example of non-structured
uncertainty
fprintf('***1st check for singular values of Delta_tilde:');
svd(Delta_tilde) % less than or equal to 1
fprintf('***2nd check for closed-loop stability of M2(Ghat(gamma2) and K2) and
Delta_tilde:');
M2 = lft(mdiag(1,1/gam2,1)*Ghat, K2, 1, 1);
clp2 = lft(Delta_tilde, M2, 2, 2);
real(eig(clp2.a)) % closed-loop stability regardless the structure of
uncertainty block
clp2
clp2.a
max(eig(clp2.a))
fprintf('***3rd check for closed-loop stability of M3(Ghat(gamma3) and K3) and
Delta_tilde:');
M3 = lft(mdiag(1,1/gam3,1)*Ghat, K3, 1, 1);
clp3 = lft(Delta_tilde, M3, 2, 2);
real(eig(clp3.a)) % closed-loop instability by the non-structured uncertainty
block
fprintf('***4th check for closed-loop H infinity norm of M3:');
norm(M3, 'inf') % larger than 1
fprintf('***5th check for closed-loop H infinity norm of M3 with scaling:');
W = mdiag(d_opt,1);
M3_d = W\M3*W;
norm(M3_d, 'inf') % less than 1
d_opt
ex1126_2
ctrlpref
ex1126_2
346/(4*1.62)

!!! the remaining page is under construction (the contents below are from last year) !!!

### [lecture #13] 2019.12.5 Control system design for practical system (2/3) †

1. return of mini report #2; ... You will have a mini exam #2 related to this report next week
2. review of the experimental system
• closed-loop system of 2-by-2 plant G and controller K
• closed-loop gain is desired to be minimized for constant speed operation
• frequency response data of G can be used; how to handle modeling error of G ?
3. design example (modeling error for Gyu is only considered for simplicity)
• frequency response experiment data
servo1.dat
servo2.dat
• determination of plant model(nominal plant and additive uncertainty weight)
• configuration of generalized plant and controller design by scaled H infinity control problem using one-dimensional search on the scaling d
• comparison of closed-loop gain characteristics with and without control
• result of control experiment
result.dat
4. final report and remote experimental system
1. design your controller(s) so that the system performance is improved compared with the design example
2. Draw the following figures and explain the difference between two control systems (your controller and the design example):
1. bode diagram of controllers
2. gain characteristic of closed-loop system from w to z
3. time response of control experiment
3. Why is the performance of your system improved(or unfortunately deteriorated)?
• due date: 6th(Mon) Jan 17:00
• submit your report(pdf or doc) by e-mail to kobayasi@nagaokaut.ac.jp
• You can use Japanese
• maximum controller order is 20
• the system will be started until next lecture
• You can send up to 10 controllers
• control experimental results will be uploaded here
• freqresp ... frequency response will be measured and uploaded everyday
5. how to improve the performance ?
• accuracy of the nominal(physical) model
• weighting for robust stability
6. specifications of the experimental system
1. program sources for frequency response experiment
• freqresp.h
• freqresp_module.c
• freqresp_app.c
• format of servo1.dat (w is used instead of u for servo2.dat)
1st column ... frequency (Hz)
2nd column ... gain from u(Nm) to y(rad/s)
3rd column ... phase (deg) from u to y
4th column ... gain from u to z
5th column ... phase (deg) from u to z
2. program sources for control experiment
• hinf.h
• hinf_module.c
• hinf_app.c
• format of result.dat
1st column: time (s)
4th column: u (Nm)
5th column: w (Nm)
3. configuration of control experiment
• disturbance signal w is specified as described in hinf.h and hinf_module.c:
w = 0; // disturbance torque for driven motor
if((t > 2)&&(t < 3)){
w = RATED_TORQ * -0.15;
}
if((t > 4)&&(t < 5)){
w = RATED_TORQ * -0.1 * sin(2*M_PI*5.0 * (t-4.0));
}
da_conv(torq_volt_conv_1(w), 1);
• control signal u is limited as specified in hinf.h and hinf_module.c:
#define U_MAX (RATED_TORQ / 3.0)

if(u > U_MAX) u = U_MAX;
if(u < -U_MAX) u = -U_MAX;
u is generated by PI control for t < 1(s). Your designed controller is started at t = 1(s).
4. calculation of rotational speed
• The rotational speed is approximately calculated by using difference for one sampling period in hinf_module.c and freqresp_module.c like:
theta_rad = (double)read_theta(0) / (double)Pn212 * 2.0 * M_PI;
theta_rad_before = theta_rad;
where the sampling period is given as 0.25 ms.
%-- 2019/12/05 13:03 --%
pwd
nominal
weight
cont
compare
ctrlpref
compare
perf
• Q: 実機パラメータのイナーシャはどうやって試験をして導出しているか。
• A: 周波数応答実験の結果に物理モデルが近くなるように調整しました。
• Q: C_S, C_L はダンパなのか。C_S はどこの部分のダンパというかたちになるのか。
• A: C_S は軸の捩りに伴う減衰係数です（これが0の場合、二つの慣性の捩り振動は減衰することなく永久に持続する）。C_L は従動側慣性の回転に伴う減衰係数です（これが0の場合、二慣性系全体は与えられた初速に応じて永久に回転し続ける）。

### [lecture #14] 2019.12.12 Control system design for practical system (3/3) †

• review & supplemental explanations
• final report
• control objective is to suppress speed fluctuation of the driven-side motor not the drive-side motor (a large vibration at the drive-side motor is allowed)
• no strict control objective is given ( there is a freedom to define what is good performance; a frequency dependent weighting function can be introduced to evaluate the performance )
• c2d() is used to discretize the resultant continuous-time controller
• web based remote experiment system
• now you can login after registration
• room temperature is displayed and stored in temp.txt (Bosch Sensortec BME280 is used)
• preparation of your own controller(s) by using the remote experiment system
• mini exam #2 (14:00〜)

### [lecture #15] 2019.12.19 Control system design for practical system (cont.) †

• return of mini exam #2
• preparation of your own controller(s)
• use design example #3 for comparison with your design
• discuss relationship between the required figures (a), (b), and (c)
• driving and driven torque are set to zero when driven-motor speed z or driving-motor speed y exceeds a limit as in hinf.h:
#define SPEED_MAX 60.0 // driving torque becomes zero when rotational speed is out of range from -60.0 rad/s to 60.0 rad/s in order to avoid shutting down the driving motor by the maximum speed excess alarm (2018.12.19)
and hinf_module.c:
       int flag_speed_excess = 0; // 2018.12.19

if(flag_speed_excess == 1) w = 0; // 2018.12.20
da_conv(torq_volt_conv_1(w), 1);

if((fabs(z) > SPEED_MAX) || (fabs(y) > SPEED_MAX)) flag_speed_excess = 1; // 2018.12.20 further modified after today's lecture
		if(flag_speed_excess == 1) u = 0; // 2018.12.19

#ifndef NO_CONTROL
da_conv(torq_volt_conv_0(u), 0);
#endif
• questionnaires
• to university
• for web-based experimental environment
%-- 2019/12/19 13:31 --%
nominal
weight
cont
compare
perf
nominal
weight
cont
compare
perf
weight