制御系CAD

制御を意識したモノつくりを目指して

線形系の応答計算

投稿日時: 2022年8月17日投稿者: cacsd

[1]自由系の時間応答
a) 次の１次自由系の時間応答を重ねて描け。
$\dot{x}=-\frac{1}{T}x,\quad x(0)=1$
$\quad(T=0.5,1,2; K=1)$
b) 次の２次自由系の時間応答を重ねて描け。
$\dot{x}= \left[\begin{array}{cc} 0 & 1\\ -\omega_n^2 & -2\zeta\omega_n \end{array}\right]x ,\quad x(0)= \left[\begin{array}{cc} 1 \\ 0 \end{array}\right]$
$y=\left[\begin{array}{cc} 1 & 0 \end{array}\right]x$
$\quad(\zeta=0.01,\frac{1}{\sqrt{2}},1; \omega_n=1)$

//resp1.sce
//—–
function dx=f(t,x),dx=a*x,endfunction
x0=1; t0=0;t=0:0.1:10;
a=-1/0.5;
x1=ode(x0,t0,t,f);
a=-1/1;
x2=ode(x0,t0,t,f);
a=-1/2;
x3=ode(x0,t0,t,f);
scf(1);
plot(t,x1,t,x2,t,x3)
mtlb_grid
//—–
function dx=f(t,x),dx=A*x,endfunction
x0=[1;0];t0=0;t=0:0.1:10;
A=[0 1;-1 -0.02]; C=[1 0];
x1=ode(x0,t0,t,f); y1=C*x1;
A=[0 1;-1 -sqrt(2)]; C=[1 0];
x2=ode(x0,t0,t,f); y2=C*x2;
A=[0 1;-1 -2]; C=[1 0];
x3=ode(x0,t0,t,f); y3=C*x3;
scf(2);
plot(t,y1,t,y2,t,y3)
mtlb_grid
//—–
//eof

[2]インパルス応答
a) 次の１次系のインパルス応答を重ねて描け。
$\dot{x}=-\frac{1}{T}x+\frac{K}{T}u$
$\quad(T=0.5,1,2; K=1)$
b) 次の２次系のインパルス応答を重ねて描け。
$\dot{x}= \left[\begin{array}{cc} 0 & 1\\ -\omega_n^2 & -2\zeta\omega_n \end{array}\right]x + \left[\begin{array}{cc} 0 \\ \omega_n^2 \end{array}\right]u$
$y=\left[\begin{array}{cc} 1 & 0 \end{array}\right]x$
$\quad(\zeta=0.01,\frac{1}{\sqrt{2}},1; \omega_n=1)$

//resp2.sce
//—–
function dx=f(t,x),dx=a*x,endfunction
t0=0; t=0:0.01:5;
a=-1/0.5; b=1/0.5;
x1=ode(b,t0,t,f);
a=-1; b=1;
x2=ode(b,t0,t,f);
a=-1/2; b=1/2;
x3=ode(b,t0,t,f);
scf(1);
plot(t,x1,t,x2,t,x3)
mtlb_grid
//—–
function dx=f(t,x),dx=A*x,endfunction
t0=0;t=0:0.1:10;
A=[0 1;-1 -0.02]; B=[0;1]; C=[1 0];
x1=ode(B,t0,t,f); y1=C*x1;
A=[0 1;-1 -sqrt(2)]; B=[0;1]; C=[1 0];
x2=ode(B,t0,t,f); y2=C*x2;
A=[0 1;-1 -2]; B=[0;1]; C=[1 0];
x3=ode(B,t0,t,f); y3=C*x3;
scf(2);
plot(t,y1,t,y2,t,y3,t,1)
mtlb_grid
//—–
//eof

[3]正弦波入力に対する時間応答
a) 次の１次自由系の $u=\sin\,10t$ に対する時間応答を，零入力応答と零状態応答と重ねて描け。
$\dot{x}=-\frac{1}{T}x+\frac{K}{T}u,\quad x(0)=1$
$\quad(T=1; K=1)$
b) 次の２次自由系の $u=\sin\,10t$ に対する時間応答を，零入力応答と零状態応答と重ねて描け。
$\dot{x}= \left[\begin{array}{cc} 0 & 1\\ -\omega_n^2 & -2\zeta\omega_n \end{array}\right]x + \left[\begin{array}{cc} 0 \\ \omega_n^2 \end{array}\right]u ,\quad x(0)= \left[\begin{array}{cc} 1 \\ 0 \end{array}\right]$
$y=\left[\begin{array}{cc} 1 & 0 \end{array}\right]x$
$\quad(\zeta=\frac{1}{\sqrt{2}}; \omega_n=1)$

//resp3.sce
//—–
function dx=f1(t,x),dx=a*x,endfunction
function dx=f2(t,x),dx=a*x+b*u(t),endfunction
function ut=u(t),ut=sin(10*t),endfunction
a=-1; b=1;
x0=1; t0=0; t=0:0.01:5;
x1=ode(x0,t0,t,f1);
x0=0; t0=0; t=0:0.01:5;
x2=ode(x0,t0,t,f2);
x0=1; t0=0; t=0:0.01:5;
x3=ode(x0,t0,t,f2);
scf(1);
plot(t,x1,t,x2,t,x3)
mtlb_grid
//—–
function dx=f1(t,x),dx=A*x,endfunction
function dx=f2(t,x),dx=A*x+B*u(t),endfunction
function ut=u(t),ut=sin(10*t),endfunction
A=[0 1;-1 -sqrt(2)]; B=[0;1]; C=[1 0];
x0=[1;0]; t0=0; t=0:0.1:10;
x1=ode(x0,t0,t,f1); y1=C*x1;
x0=[0;0]; t0=0; t=0:0.1:10;
x2=ode(x0,t0,t,f2); y2=C*x2;
x0=[1;0]; t0=0; t=0:0.1:10;
x3=ode(x0,t0,t,f2); y3=C*x3;
scf(2);
plot(t,y1,t,y2,t,y3,t,1)
mtlb_grid
//—–
//eof

[4]ステップ応答
a) 次の１次系のステップ応答を重ねて描け。
$\dot{x}=-\frac{1}{T}x+\frac{K}{T}u$
$\quad(T=-0.5,-1,2; K=1)$
b) 次の２次系のステップ応答を重ねて描け。
$\dot{x}= \left[\begin{array}{cc} 0 & 1\\ -\omega_n^2 & -2\zeta\omega_n \end{array}\right]x + \left[\begin{array}{cc} 0 \\ \omega_n^2 \end{array}\right]u$
$y=\left[\begin{array}{cc} 1 & 0 \end{array}\right]x$
$\quad(\zeta=0.01,\frac{1}{\sqrt{2}},1; \omega_n=1)$

//resp4.sce
//—–
function dx=f(t,x),dx=a*x+b*u(t),endfunction
function ut=u(t),ut=1,endfunction
x0=0; t0=0; t=0:0.01:5;
a=-1/0.5; b=1/0.5;
x1=ode(x0,t0,t,f);
a=-1; b=1;
x2=ode(x0,t0,t,f);
a=-1/2; b=1/2;
x3=ode(x0,t0,t,f);
scf(1);
plot(t,x1,t,x2,t,x3,t,1-1/%e)
mtlb_grid
//—–
function dx=f(t,x),dx=A*x+B*u(t),endfunction
function ut=u(t),ut=1,endfunction
x0=[0;0];t0=0;t=0:0.1:10;
A=[0 1;-1 -0.02]; B=[0;1]; C=[1 0];
x1=ode(x0,t0,t,f); y1=C*x1;
A=[0 1;-1 -sqrt(2)]; B=[0;1]; C=[1 0];
x2=ode(x0,t0,t,f); y2=C*x2;
A=[0 1;-1 -2]; B=[0;1]; C=[1 0];
x3=ode(x0,t0,t,f); y3=C*x3;
scf(2);
plot(t,y1,t,y2,t,y3,t,1)
mtlb_grid
//—–
//eof

[5]周波数応答
a) 次の１次系の周波数応答を重ねて描け。
$\dot{x}=-\frac{1}{T}x+\frac{K}{T}u$
$\quad(T=-0.5,-1,2; K=1)$
b) 次の２次系の周波数応答を重ねて描け。
$\dot{x}= \left[\begin{array}{cc} 0 & 1\\ -\omega_n^2 & -2\zeta\omega_n \end{array}\right]x + \left[\begin{array}{cc} 0 \\ \omega_n^2 \end{array}\right]u$
$y=\left[\begin{array}{cc} 1 & 0 \end{array}\right]x$
$\quad(\zeta=0.01,\frac{1}{\sqrt{2}},1; \omega_n=1)$

//resp5.sce
//—–
function g=f(a,b,w), g=b./(%i*w-a), endfunction
w=logspace(-1,1);
a=-1/2; b=1/2; g=f(a,b,w); ga1=20*log10(abs(g));
ph1=180/%pi*atan(imag(g),real(g));
a=-1/1; b=1/1; g=f(a,b,w); ga2=20*log10(abs(g));
ph2=180/%pi*atan(imag(g),real(g));
a=-1/0.5; b=1/0.5; g=f(a,b,w); ga3=20*log10(abs(g));
ph3=180/%pi*atan(imag(g),real(g));
scf(1);
subplot(211),plot2d(w,ga1,logflag=’ln’),
plot2d(w,ga2,logflag=’ln’),
plot2d(w,ga3,logflag=’ln’),mtlb_grid
subplot(212),plot2d(w,ph1,logflag=’ln’),
plot2d(w,ph2,logflag=’ln’),
plot2d(w,ph3,logflag=’ln’),mtlb_grid
//—–
function g=f(A,B,C,w),
g=[]; N=length(w)
for i=1:N, g=[g C*((%i*w(i)*eye(A)-A)\B)]; end
endfunction
B=[0;1]; C=[1 0]; w=logspace(-1,1,200);
A=[0 1;-1 -2*0.01]; g=f(A,B,C,w); ga1=20*log10(abs(g));
ph1=180/%pi*atan(imag(g),real(g));
A=[0 1;-1 -2/sqrt(2)]; g=f(A,B,C,w); ga2=20*log10(abs(g));
ph2=180/%pi*atan(imag(g),real(g));
A=[0 1;-1 -2*1]; g=f(A,B,C,w); ga3=20*log10(abs(g));
ph3=180/%pi*atan(imag(g),real(g));
scf(2);
subplot(211),plot2d(w,ga1,logflag=’ln’),
plot2d(w,ga2,logflag=’ln’),
plot2d(w,ga3,logflag=’ln’),mtlb_grid
subplot(212),plot2d(w,ph1,logflag=’ln’),
plot2d(w,ph2,logflag=’ln’)
plot2d(w,ph3,logflag=’ln’),mtlb_grid
//—–
//eof

演習

投稿日時: 2022年8月17日投稿者: cacsd

[1] 次の行列を定義せよ。

$1^\circ\ A=\left[\begin{array}{cc} 0 & 1 \\ -1 & -0.02 \end{array}\right],\ B=\left[\begin{array}{cc} 0 \\ 1 \end{array}\right],\ C=\left[\begin{array}{cc} 1 & 0 \end{array}\right]$

$2^\circ\ A=\left[\begin{array}{cccc} 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 1 \\ -1 & 0 & 0 & 0 \\ 0 & -1 & 0 & 0 \end{array}\right],\ B=\left[\begin{array}{cc} 0 & 0 \\ 0 & 0 \\ 1 & 1 \\ 1 & -1 \end{array}\right]$

$3^\circ\ A=\left[\begin{array}{cccc} 0 & 1 & 0 & 0 \\ 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 1 \\ -1 & -4 & -6 & -4 \end{array}\right],\ B=\left[\begin{array}{cc} 0 \\ 0 \\ 0 \\ 1 \end{array}\right]$

[2] 行列

$E=\left[\begin{array}{ccccc} 1 & 2 & 1 & 0 & 0\\ 3 & 4 & 0 & 1 & 0\\ 0 & 0 & 1 & 0 & 0\\ 0 & 0 & 0 & 2 & 0\\ 0 & 0 & 0 & 0 & 3 \end{array}\right]$

に対して

$1^\circ$ 行列Eの列ベクトルを逆順に並べよ。

$2^\circ$ 行列Eの行ベクトルを逆順に並べよ。

$3^\circ$ 行列Eの奇数番目の列ベクトルと偶数番目の行ベクトルからなる小行列を求めよ。

[3]

$1^\circ$ 行列

$A=\left[\begin{array}{ccccc} 0 & 1 \\ -2 & -3 \end{array}\right]$

の固有値分解 $VRV^{-1}$ を求め，つぎを確かめよ。

$A=VRV^{-1}$

$2^\circ$ 行列

$C=\left[\begin{array}{ccccc} 1 & 0 & 0 & 0 \\ 1 & 1 & 0 & 0 \end{array}\right]$

の特異値分解 $USV^T$ を求め，つぎを確かめよ。

$C=USV^T$
$U^TU=I$
$V^TV=I$

$3^\circ$ x=[3 1 5 4 2] の要素を大きい順に並べよ。

[4] 平均 $\mu$ ，分散 $\sigma^2$ の正規分布密度関数

$p(x)=\frac{1}{\sqrt{2\pi}\sigma}e^{-\frac{(x-\mu)^2}{2\sigma^2}} \quad(-3\le x\le3)$

に対して

(1) $\mu=0$ ， $\sigma=1$ ,
(2) $\mu=1$ ， $\sigma=0.5$ ,
(3) $\mu=-1$ ， $\sigma=2$

の場合のグラフを重ねて描け。

[5] 与えられた中心の座標 $(x,y)$ と半径 $r$ に対して，円を描く関数を作成し

(1) $(x,y)=(0,0)$ ， $r=1$
(2) $(x,y)=(1,1)$ ， $r=2$
(3) $(x,y)=(-1,-1)$ ， $r=3$

の場合の円を重ねて描け。

SCILABの基礎事項

投稿日時: 2022年8月17日投稿者: cacsd

[0] はじめに
SCILABは対話型数値計算プログラムで、次からフリーで入手できます。
http://www.scilab.org/products/scilab/download

[1] 行列データの入力
１）サイズ２×２の行列
$A=\left[\begin{array}{cc}1&2\\3&4\end{array}\right]$
を定義するために、次を入力してみます。
–> A=[1 2;3 4]
ここで、行列データは[ ]で囲み、要素の区切りはスペースまたはカンマ、行替えはセミコロン；を用いることに注意してください。次のように最後に；を置くと、出力が抑制されます。
–> A=[1 2;3 4];
すでに入力した行列データへのアクセスは、次のように行います。
–> A
また、特定の要素へのアクセスや修正は、次のように行います。
–> A(2,1)
–> A(2,1)=-5

２）２次の単位行列
$B=\left[\begin{array}{cc}1&0\\0&1\end{array}\right]$
は、次により入力できます。
–> B=eye(2,2)
便利なのは、次のように異なる行数と列数が許されることです。
–> B=eye(2,3)
これは、次のサイズ２×３の行列
$B=\left[\begin{array}{ccc}1&0&0\\0&1&0\end{array}\right]$
の定義を意味します。

３）サイズ３×２の零行列
$C=\left[\begin{array}{cc}0&0\\0&0\\0&0\end{array}\right]$
は、次により定義できます。
–> C=zeros(3,2)

４）３次の対角行列
$D=\left[\begin{array}{ccc}1&0&0\\0&2&0\\0&0&3\end{array}\right]$
は、次により定義できます。
–> D=diag([1 2 3])
すなわち、ベクトルにdiagを作用させているわけですが、行列にdiagを作用させ、対角成分を抽出することもできます。
–> d=diag(D)

５）行列A,B,C,Dを小行列としてもつブロック行列
$E=\left[\begin{array}{cc}A&B\\C&D\end{array}\right]$
を定義するためには、次のように入力します。
–> E=[A B;C D]

６）サイズ２×２の複素行列
$F=\left[\begin{array}{cc}1+i&2+3i\\3+2i&4+4i\end{array}\right]$
も、次により定義できます。
–> F=[1+%i, 2+3*%i;3+2*%i, 4+4*%i]
ここで、%iは虚数単位を意味します。

７）SCILABでは、空行列を定義できます。
–> G=[]

８）すべての要素が1である行列
$H=\left[\begin{array}{cc}1&1\\1&1\end{array}\right]$
は、次により定義できます。
–> H=ones(2,2)

[2] 行列データの管理
１）現時点でどのような行列を定義したかを調べるためには、次のように入力します。
–> who_user
また、たとえば入力した行列データdを消去するためには、次のように入力します。
–> clear d
すべてのデータを消去するときは、次のように入力します。
–> clear

２）行列データEのサイズを調べるには、次のように入力します。
–> [m,n]=size(E)

３）浮動小数点表示を行うためには、次のように入力します。
–> format(’e’)
これを固定小数点表示に戻すためには、次のように入力します。
–> format(’v’)

４）カレントディレクトリのチェックを行うためには、次のコマンドが使えます
–> pwd
–> ls
カレントディレクトリの変更のためには、次のコマンドを与えます。
–> cd(’dname’)

５）メモリ上のデータのセーブとロードのためには、次のコマンドを与えます。
–> save(’fname’)
–> load(’fname’)

[3] コロン記法
１）まず、次を入力してみてください。
–> h=1:3
これは、「form 1 to 3」を意味します。すなわちサイズ1×5の行列h=[1 2 3]が生成されます。コロン「:」はtoの意味で使われ、colon notationと呼ばれています。

２）それでは、次はどうでしょう。
–> h=0:0.5:3
これは、「form 1 to 3 by 0.5」を意味します。すなわち、２つのコロンで挟まれた0.5は増分を意味し、サイズ1×5の行列h=[1 1.5 2 2.5 3]が生成されます。ちなみに、h=1:3はh=1:1:3と同じと考えておきましょう。

３）それでは、次はどうでしょう。
–> h=0:-2
これは、h=1:1:-2の意味になりますが、1から-2まで、増分1では到達できないので、空行列となります。

[4] 小行列の抽出
１）コロン記法の利点は，行列の添え字部分に適用できることにあります。たとえば，行列Eの(2,3)要素には
–> E(2,3)
でアクセスできるのですが，この2と3は，スカラーである必要はなく，ベクトルが使えます。たとえば，次を確かめてください。
–> E([1 3],: )
これは第１行と弟２行からなるサイズ2×5の小行列が取り出されることに注意してください。列を指定する「：」はすべての列を意味しています。

２）次はどうでしょう。
–> E(:,2:4 )
ここで、2:4は[2 3 4]というベクトルでしたよね。第2～4列からなるサイズ5×3の小行列が取り出されることに注意してください。行を指定する「：」はすべての行を意味しています。

３）それでは，次はどうでしょう。
–> E([1 3],2:4 )

[5] 行列の分割
１）いろいろな応用では、ある行列を分割する必要が出てきます。これは小行列の抽出の観点から、たとえば次のようにできます。
–> [m,n]=size(E); i=2; j=2;
–> A=E(1:i,1:j);
–> B=E(1:i,j+1:n);
–> C=E(i+1:m,1:j);
–> D=E(i+1:m,j+1:n);
ここで、i=0のときは1:i=[]、i=mのときはi+1:m=[]などに注意してください。

[6] 行列演算
いま、行列Aが次のように定義されていることを確認してください。
$A=\left[\begin{array}{cc}1&2\\-5&4\end{array}\right]$
また、行列Bを新しく定義します。
–> B=eye(A)

１）行列A,Bの和、差、乗算は次のように計算できます。
–> X=A+B
–> X=A-B
–> X=A*B

２）行列Aが正則のとき、その逆行列は
–> X=inv(A)
で計算されます。しかし、これは逆行列の要素が欲しいとき以外は使わないようにしましょう。
多くの場合、逆行列は次の連立一次方程式
$AX=B$
または
$XA=B$
を解くときに用いますが、SCILABでは、次のように計算します。
–> X=A\B
–> X=B/A
それぞれ左からの割り算、右からの割り算と呼ばれます。
これは数値的に安定な数値計算アルゴリズムを得るためは直交行列だけの演算ですませることが必要だからです。

３）行列Aのべき乗と転置共役は次のように計算でされます。
–> X=A^2
–> X=A’

４）それでは、次の演算結果と比較してみましょう。
–> X=A.*B
–> X=A.\B
–> X=B./A
–> X=A.^2
–> X=A.’
すなわち、演算子の前に「ドット．」をおくと、要素ごとの演算やただの転置になります。

５）次に、組み込み関数としては、次が使えます。
sin,cos,tan
asin,acos,atan
sinh,cosh,tanh
asinh,acosh,atanh
sqrt
exp,log,log10
abs
real,imag

６）行列関数としては次が準備されています。
–> X=expm(A)
–> X=norm(A)
–> X=rank(A)
–> X=det(A)
ぞれぞれ、行列指数関数、ノルム、階数、行列式を計算します。

７）また、固有値問題を解くには次の関数を用います。
–> [V,R]=spec(A)
単に固有値だけを求めるのであれば
–> r=spec(A)
でＯＫです。

８）さらに、特異値問題を解くには次の関数を用います。
–> [U,S,V]=svd(A)
単に特異値だけを求めるのであれば
–> s=svd(A)
でＯＫです。

９）多項式 $x+1$ と $x+2$ の掛け算は、次のように行います。
–> r=convol([1 1],[1 2])

１０）いま、次のベクトルを定義します。
–> x=[3 1 5 4 2]
このとき、次はxの要素を大きい順番に並べ替えます。
–> [y,i]=gsort(x)
yがソーティングされた結果、iが元のどの位置にあったかを示すインデックスです。たとえば
–> x（i)
としてみるとyと同じ結果が出ます。

１１）次は、ベクトルの要素和、最大値、最小値を求める関数です。
–> y=sum(x)
–> y=max(x)
–> y=min(x)

[7] 関数のグラフ
１）次の関数のグラフを描いてみます。
$p(x)=\frac{1}{2\pi}e^{-\frac{x^2}{2}}\quad(-3\le x\le3)$
次のコマンドを与えてください。
–> x=-3:0.1:3;
–> y=1/sqrt(2*%pi)*exp(-x.*x/2)；
–> plot(x,y)
ここで，なぜx*xでなく，x.*xとしているか考えてみてください。

２）描画範囲（-3<=x<=3, 0<=ｙ<=0.4）を指定し，グリッドを描画するには，次のようなコマンドを追加します。
–> mtlb_grid
–> mtlb_axis([-3 3 0 0.5])
–> xtitle(’N(0,1)’)
ここで，「mtlb_」はMATLABのコマンド相当という意味です。

３）次の１次系の周波数応答を片対数グラフを描いてみましょう。
$|\hat{g}(j\omega)|=20\log_{10}\frac{1}{\sqrt{1+\omega^2}}\quad(0.1\le x\le10)$
次のコマンドを与えてください。
–> w=logspace(-1,1);
–> ga=20*log10(sqrt(1+w.^2).^(-1));
–> plot2d(w,ga,logflag=’ln’)
–> mtlb_grid

[8] プログラミング
これまでは対話型入力を行ってきましたが，一連のコマンドをエディタを用いてあらかじめ編集しておいて一括処理する，すなわちscilabでのプログラミングについて学びます。まず，白紙を示すアイコンをクリックして，エディタscipadを起動してください。

１） 1から9までの奇数の和を求める次のプログラムを入力し，ファイル名ex1.sceでセーブします。
//ex1.sce
//—–
s=0;
for i=1:2:9, s=s+i; end
disp(s,’sum=’)
//—–
//eof
scipadのexecuteコマンドをクリックして実行します。答えは25のはずです。

２）マシンイプシロン(machine epsilon)と呼ばれる特別な定数を求めてみましょう。次のプログラムをファイル名ex2.sceでセーブします。
//ex2.sce
//—–
s=1;
while 1+s>1, s=s/2; end
s=s*2;
disp(s,’eps=’),
//—–
//eof
scipadのexecuteコマンドをクリックして実行します。epsは次の%epsと等しくなります。
–-> %eps
次を実行してみてください。
–> 1+%eps/2==1
すなわち，
マシンイプシロン%eps=2^{-k}は，1+2^{-k-1}>1が成り立たなくなる最小のkを用いて定義されます。浮動小数点による数値計算では，1+%eps=1となる正数が存在することに注意してください。

３）次は，if文を使った例です。次のプログラムをファイル名ex3.sceでセーブします。
//ex3.sce
//—–
x=input(’x= ‘);
if imag(x)~=0, disp(’complex number’)
elseif x>0, disp(’positive real number’)
elseif x<0, disp(’negative real number’)
else disp(’zero’)
end
//—–
//eof
scipadのexecuteコマンドをクリックして実行します。

４）関数の例として，行列の階数を計算するものを示しましょう。次のプログラムをファイル名rank2.sceでセーブします。
//rank2.sce
//—–
function r=rank2(A)
s=svd(A); tol=max(size(A))*s(1)*%eps; r=sum(s>tol);
endfunction
//—–
//eof
このとき，次を入力してください。
–> exec(‘rank2.sce’)
–> r=rank2([1 2;3 4])

５）次の微分方程式の解を求めることを考えます。
$\dot{x}(t)=-x(t)+\sin(10t)$
次のプログラムをファイル名ex4.sceでセーブします。
//ex4.sce
//—–
function dx=f(t,x), dx=-x+sin(10*t); endfunction
x0=1; t0=0; t=0:0.01:5;
x=ode(x0,t0,t,f);
plot(t,x), mtlb_grid
//—–
//eof
scipadのexecuteコマンドをクリックして実行します。

MATLABの基礎事項

投稿日時: 2022年8月17日投稿者: cacsd

[1] 行列データの入力 　
１）サイズ２×２の行列
$A=\left[\begin{array}{cc}1&2\\3&4\end{array}\right]$
を定義するために、次を入力してみます。
>> A=[1 2;3 4]
ここで、行列データは[ ]で囲み、要素の区切りはスペースまたはカンマ、行替えはセミコロン；を用いることに注意してください。次のように最後に；を置くと、出力が抑制されます。
>> A=[1 2;3 4];
すでに入力した行列データへのアクセスは、次のように行います。
>> A
また、特定の要素へのアクセスや修正は、次のように行います。
>> A(2,1)
>> A(2,1)=-5

２）２次の単位行列
$B=\left[\begin{array}{cc}1&0\\0&1\end{array}\right]$
は、次により入力できます。
>> B=eye(2,2)
便利なのは、次のように異なる行数と列数が許されることです。
>> B=eye(2,3)
これは、次のサイズ２×３の行列
$B=\left[\begin{array}{ccc}1&0&0\\0&1&0\end{array}\right]$
の定義を意味します。

３）サイズ３×２の零行列
$C=\left[\begin{array}{cc}0&0\\0&0\\0&0\end{array}\right]$
は、次により定義できます。
>> C=zeros(3,2)

４）３次の対角行列
$D=\left[\begin{array}{ccc}1&0&0\\0&2&0\\0&0&3\end{array}\right]$
は、次により定義できます。
>> D=diag([1 2 3])
すなわち、ベクトルにdiagを作用させているわけですが、行列にdiagを作用させ、対角成分を抽出することもできます。
>> d=diag(D)

５）行列A,B,C,Dを小行列としてもつブロック行列
$E=\left[\begin{array}{cc}A&B\\C&D\end{array}\right]$
を定義するためには、次のように入力します。
>> E=[A B;C D]

６）サイズ２×２の複素行列
$F=\left[\begin{array}{cc}1+i&2+3i\\3+2i&4+4i\end{array}\right]$
も、次により定義できます。
>> F=[1+i, 2+3*i;3+2*i, 4+4*i]
ここで、iは虚数単位として定義されています。もちろん再定義可能です。

７）MATLABでは、空行列を定義できます。
>> G=[]

８）すべての要素が1である行列
$H=\left[\begin{array}{cc}1&1\\1&1\end{array}\right]$
は、次により定義できます。
>> H=ones(2,2)

[2] 行列データの管理
１）現時点でどのような行列を定義したかを調べるためには、次のように入力します。
>> who
行列のザイズを合わせて表示するには、次のように入力します。
>> whos
また、たとえば入力した行列データdを消去するためには、次のように入力します。
>>clear d
すべてのデータを消去するときは、次のように入力します。
>>clear

２）行列データEのサイズを調べるには、次のように入力します。
>> [m,n]=size(E)

３）浮動小数点表示を行うためには、次のように入力します。
>> format short e
これを固定小数点表示に戻すためには、次のように入力します。
>> format

４）カレントディレクトリのチェックを行うためには、次のコマンドが使えます
>> pwd
>> ls
カレントディレクトリの変更のためには、次のコマンドを与えます。
>> cd dname

５）メモリ上のデータのセーブとロードのためには、次のコマンドを与えます。
>> save filename
>> load filename

[3] コロン記法
１）まず、次を入力してみてください。
>> h=1:3
これは、「form 1 to 3」を意味します。すなわちサイズ1×5の行列h=[1 2 3]が生成されます。コロン「:」はtoの意味で使われ、colon notationと呼ばれています。

２）それでは、次はどうでしょう。
>> h=0:0.5:3
これは、「form 1 to 3 by 0.5」を意味します。すなわち、２つのコロンで挟まれた0.5は増分を意味し、サイズ1×5の行列h=[1 1.5 2 2.5 3]が生成されます。ちなみに、h=1:3はh=1:1:3と同じと考えておきましょう。

３）それでは、次はどうでしょう。
>> h=0:-2
これは、h=1:1:-2の意味になりますが、1から-2まで、増分1では到達できないので、空行列となります。

[4] 小行列の抽出
１）コロン記法の利点は，行列の添え字部分に適用できることにあります。たとえば，行列Eの(2,3)要素には
>> E(2,3)
でアクセスできるのですが，この2と3は，スカラーである必要はなく，ベクトルが使えます。たとえば，次を確かめてください。
>> E([1 3],: )
これは第１行と弟２行からなるサイズ2×5の小行列が取り出されることに注意してください。列を指定する「：」はすべての列を意味しています。

２）次はどうでしょう。
>> E(:,2:4 )
ここで、2:4は[2 3 4]というベクトルでしたよね。第2～4列からなるサイズ5×3の小行列が取り出されることに注意してください。行を指定する「：」はすべての行を意味しています。

３)それでは，次はどうでしょう。
>> E([1 3],2:4 )

[5] 行列の分割
１）いろいろな応用では、ある行列を分割する必要が出てきます。これは小行列の抽出の観点から、たとえば次のようにできます。
>> [m,n]=size(E); i=2; j=2;
>> A=E(1:i,1:j);
>> B=E(1:i,j+1:n);
>> C=E(i+1:m,1:j);
>> D=E(i+1:m,j+1:n);
ここで、i=0のときは1:i=[]、i=mのときはi+1:m=[]などに注意してください。

[6] 行列演算
いま、行列Aが次のように定義されていることを確認してください。
$A=\left[\begin{array}{cc}1&2\\-5&4\end{array}\right]$
また、行列Bを新しく定義します。
>> B=eye(A)

１）行列A,Bの和、差、乗算は次のように計算できます。
>> X=A+B
>> X=A-B
>> X=A*B

２）行列Aが正則のとき、その逆行列は
>> X=inv(A)
で計算されます。しかし、これは逆行列の要素が欲しいとき以外は使わないようにしましょう。
多くの場合、逆行列は次の連立一次方程式
$AX=B$
または
$XA=B$
を解くときに用いますが、MATLABでは、次のように計算します。
>> X=A\B
>> X=B/A
それぞれ左からの割り算、右からの割り算と呼ばれます。
これは数値的に安定な数値計算アルゴリズムを得るためは直交行列だけの演算ですませることが必要だからです。

３）行列Aのべき乗と転置共役は次のように計算でされます。
>> X=A^2
>> X=A’

４）それでは、次の演算結果と比較してみましょう。
>> X=A.+B
>> X=A.-B
>> X=A.*B
>> X=A.\B
>> X=B./A
>> X=A.^2
>> X=A.’
すなわち、演算子の前に「ドット．」をおくと、要素ごとの演算やただの転置になります。

５）次に、組み込み関数としては、次が使えます。
sin,cos,tan,asin,acos,atan,sinh,cosh,tanh,asinh,acosh,atanh,sqrt,exp,log,log10,abs,real,imag

６）行列関数としては次が準備されています。
>> X=expm(A)
>> X=norm(A)
>> X=rank(A)
>> X=det(A)
ぞれぞれ、行列指数関数、ノルム、階数、行列式を計算します。

７）また、固有値問題を解くには次の関数を用います。
–>[V,R]=eig(A)
単に固有値だけを求めるのであれば
–>r=eig(A)
でＯＫです。

８）さらに、特異値問題を解くには次の関数を用います。
–>[U,S,V]=svd(A)
単に特異値だけを求めるのであれば
–>s=svd(A)
でＯＫです。

９）多項式 $x+1$ と $x+2$ の掛け算は、次のように行います。
–>r=convol([1 1],[1 2])

１０）いま、次のベクトルを定義します。
–>x=[3 1 5 4 2]
このとき、次はxの要素を小さい順番に並べ替えます。
–>[y,i]=sort(x)
yがソーティングされた結果、iが元のどの位置にあったかを示すインデックスです。たとえば
–>x（i)
としてみるとyと同じ結果が出ます。

１１）次は、ベクトルの要素和、最大値、最小値を求める関数です。
–>y=sum(x)
–>y=max(x)
–>y=min(x)

[7] 関数のグラフ
１）次の関数のグラフを描いてみます。
$p(x)=\frac{1}{2\pi}e^{-\frac{x^2}{2}}\quad(-3\le x\le3)$
次のコマンドを与えてください。
>> x=-3:0.1:3;
>> y=1/sqrt(2*pi)*exp(-x.*x/2)；
>> plot(x,y)
ここで，なぜx*xでなく，x.*xとしているか考えてみてください。またpiは $\pi$ として定義されています。

２）描画範囲（-3<=x<=3, 0<=ｙ<=0.4）を指定し，グリッドを描画するには，次のようなコマンドを追加します。
>> grid
>> axis([-3 3 0 0.5])
>> xtitle(’N(0,1)’)

３）次の１次系の周波数応答を片対数グラフび描いてみましょう。
$|\hat{g}(j\omega)|=20\log_{10}\frac{1}{\sqrt{1+\omega^2}}\quad(0.1\le x\le10)$
次のコマンドを与えてください。
>> w=logspace(-1,1);
>> ga=20*log10(sqrt(1+w.^2).^(-1));
>> semilogx(w,ga)
>> grid

[8] プログラミング
これまでは対話型入力を行ってきましたが，一連のコマンドをエディタを用いてあらかじめ編集しておいて一括処理する，すなわちMATLABでのプログラミングについて学びます。まず，新規スクリプトを示すアイコンをクリックして，エディタを起動してください。

１） 1から9までの奇数の和を求める次のプログラムを入力し，ファイル名ex1.mでセーブします。
//ex1.m
//—–
s=0;
for i=1:2:9,
s=s+i;
end
disp(s,’sum=’)
//—–
//eof
エディタの実行ボタンをクリックします。答えは25のはずです。

２）マシンイプシロン(machine epsilon)と呼ばれる特別な定数を求めてみましょう。次のプログラムをファイル名ex2.mでセーブします。
//ex2.m
//—–
s=1;
while 1+s>1,
s=s/2;
end
s=s*2;
disp(s,’eps=’),
//—–
//eof
エディタの実行ボタンをクリックします。このsは次のepsと等しくなります。
>> eps
次を実行してみてください。
>> 1+eps/2==1
すなわち，
マシンイプシロンeps=a^{-k}は，1+a^{-k-1}>1が成り立たなくなる最小のkを用いて定義されます。浮動小数点による数値計算では，1+eps=1となる正数が存在することに注意してください。

３）次は，if文を使った例です。次のプログラムをファイル名ex3.mでセーブします。
//ex3.m
//—–
x=input(’x= ‘);
if imag(x)~=0,
disp(’complex number’)
elseif x>0,
disp(’positive real number’)
elseif x<0,
disp(’negative real number’)
else
disp(’zero’)
end
//—–
//eof
エディタの実行ボタンをクリックします。。

４）関数の例として，行列の階数を計算するものを示しましょう。次のプログラムをファイル名rank2.mでセーブします。
//rank2.m
//—–
function r=rank2(A)
s=svd(A);
tol=max(size(A))*s(1)*eps; r=sum(s>tol);
endfunction
//—–
//eof
このとき，次を入力してください。
>> r=rank([1 2;3 4])

５）次の微分方程式の解を求めることを考えます。
$\dot{x}(t)=-x(t)+\sin(10t)$
次のプログラムをファイル名ex4.mでセーブします。
//ex4.m
//—–
function dx=f(t,x),dx=-x+sin(10*t), endfunction
x0=1; t0=0;
t=0:0.01:5;
x=ode(x0,t0,t,f);
plot(t,x),
mtlb_grid
//—–
//eof
エディタの実行ボタンをクリックします。

SVD

投稿日時: 2022年8月14日投稿者: cacsd

制御系CADに関わる数値計算では実行列の特異値分解(Singular-Value Decomposition, SVD)を多用しますが、データ駆動分野の計算でも基盤技術となっています。

●いまデータは次の複素行列にまとめられているとします。

$\displaystyle{(1)\quad X=\left[\begin{array}{ccccc} x_1 & x_2 & \cdots & x_m \end{array}\right]\in{\bf C}^{n\times m} \quad(x_1,x_2,\cdots,x_n\in{\bf C}^n) }$

ここで列ベクトル $x_1,x_2,\cdots,x_n$ はn次元数ベクトルですが、データ駆動分野では、nの値は数個ではなく、数千、数万のオーダーになります。

● $X$ をサイズ $n\times m$ の複素行列（ ${\rm rank}\,X=k$ ）とします。このとき、ユニタリ行列 $U$ と $V$ が存在して

$\displaystyle{(2.1)\quad X=U\Sigma V^H\quad(U\in{\bf C}^{n\times n},\Sigma\in{\bf R}^{n\times m},V\in{\bf C}^{m\times m}) }$

が成り立ちます（ $V^H=V^{*T}$ ）。ここで、実行列 $\Sigma$ は次のように表されます。

$\displaystyle{(2.2)\quad \Sigma= \left[\begin{array}{cc} \widehat{\Sigma} & 0_{k\times m-k}\\ 0_{n-k\times k} & 0_{n-k\times m-k} \end{array}\right]\in{\bf R}^{n\times m} }$

$\displaystyle{(2.3)\quad \widehat{\Sigma} ={\rm diag}\{\sigma_1,\cdots,\sigma_k\}\quad(\sigma_1\ge\cdots\ge\sigma_k>0) }$

MATLABのコマンドでは、次のように求めます。

>> [U,S,V]=svd(X)

●いま $\boxed{n\ge m}$ とします。このとき(2)は次式となります。

$\displaystyle{(3.1)\quad X= \underbrace{\left[\begin{array}{cc} U_1 & U_2 \end{array}\right]}_{U} \underbrace{\left[\begin{array}{cc} \Sigma_1 \\ 0_{n-m\times m} \end{array}\right]}_{\Sigma}V^H =U_1\Sigma_1V^H }$

$\displaystyle{(3.2)\quad \Sigma_1={\rm diag}\{\sigma_1,\cdots,\sigma_k, 0,\cdots,0\}\in{\bf R}^{m\times m} }$

MATLABのコマンドでは、次のように求めます。

>> [U1,S1,V]=svd(X,’econ’)

(2)の場合は $U$ は大規模な行列となりますが、(3)の場合はスリムな行列となりますので、Economic SVDと呼ぶことがあります。

●いま実効階数を $r$ とします。すなわち微小な $\epsilon>0$ に対して

$\displaystyle{(4)\quad \sigma_1\ge\cdots\ge\sigma_r> \epsilon >\sigma_{r+1}\ge\cdots\ge\sigma_k>0 }$

とします。このとき(3)から次式のように近似できます。

$\displaystyle{(5.1)\quad X\approx \widetilde{X}= \underbrace{\left[\begin{array}{cc} \widetilde{U} & \cdot \end{array}\right]}_{U} \underbrace{\left[\begin{array}{cc} \widetilde{\Sigma} & 0_{r\times m-r} \\ 0_{n-r\times r} & 0_{n-r\times m-r} \end{array}\right]}_{\Sigma'} \underbrace{\left[\begin{array}{c} \widetilde{V}^H\\ \cdot \end{array}\right]}_{V^H} =\widetilde{U}\widetilde{\Sigma}\widetilde{V}^H }$

$\displaystyle{(5.2)\quad \widetilde{\Sigma} ={\rm diag}\{\sigma_1,\cdots,\sigma_r\}\quad(\sigma_1\ge\cdots\ge\sigma_r>\epsilon) }$

●(2.1)、(3.1)、(5.1)の関係はダイアディック形式を用いて次のように表されます。

$\displaystyle{(6)\quad \underbrace{\sum_{i=1}^{m}\sigma_iu_iv_i^H}_{U\Sigma V^H} =\underbrace{\sum_{i=1}^{k}\sigma_iu_iv_i^H}_{U_1\Sigma_1 V^H} \approx \underbrace{\sum_{i=1}^{r}\sigma_iu_iv_i^H}_{\widetilde{U}\widetilde{\Sigma}\widetilde{V}^H} }$

●相関行列 $XX^H$ または $X^HX$ は、(3.1)より、次式のように表されます。

$\displaystyle{(7.1)\quad \begin{array}{l} XX^H=U\left[\begin{array}{cc} \Sigma_1 \\ 0_{n-m\times m} \end{array}\right] \underbrace{V^HV}_{I_m}\left[\begin{array}{cc} \Sigma_1 & 0_{m\times n-m} \end{array}\right] U^H\\=U\left[\begin{array}{cc} \Sigma_1^2 & 0_{m\times n-m} \\ 0_{n-m\times m} & 0_{n-m\times n-m} \end{array}\right] U^H \end{array} }$

$\displaystyle{(7.2)\quad X^HX=V\left[\begin{array}{cc} \Sigma_1 & 0_{m\times n-m} \end{array}\right] \underbrace{U^HU}_{I_n}\left[\begin{array}{cc} \Sigma_1 \\ 0_{n-m\times m} \end{array}\right] V^H=V\Sigma_1^2 V^H }$

これらより

$\displaystyle{(8.1)\quad \begin{array}{l} XX^HU=U\left[\begin{array}{cc} \Sigma_1^2 & 0_{m\times n-m} \\ 0_{n-m\times m} & 0_{n-m\times n-m} \end{array}\right] \end{array} }$

$\displaystyle{(8.2)\quad X^HXV=V\Sigma_1^2 }$

すなわち、非零の特異値は相関行列 $XX^H$ または $X^HX$ の固有値の正の平方根です。もし $X=X^H$ の場合は $X$ の特異値は $X$ の固有値の絶対値となります。これはスカラの場合、次式を意味しています。

$\displaystyle{(9)\quad \sigma=\sqrt{xx^*}=\sqrt{(x_R+jx_I)(x_R-jx_I)}=\sqrt{x_R^2+x_I^2}=|x| }$

データ駆動分野では、nの値は数個ではなく、数千、数万のオーダーになります。一方、ｍの値はそのように大きな値とはなりません。したがって、(8.2)の小規模な固有値問題を解いて $\Sigma_1$ と $V$ を求め、(3.1)に基づいて、(8.1)の大規模な固有値問題の固有ベクトルの一部 $U_1$ を、次式により求めることが行われます。

$\displaystyle{(10)\quad U_1=XV\Sigma_1^{-1} }$

(5.1)に基づく場合は次式となります。

$\displaystyle{(11)\quad \widetilde{U}=X\widetilde{V}\widetilde{\Sigma}^{-1} }$

このような手法はスナップショットと呼ばれています。

●観測データを行ベクトルとしてもつ $X$ を考えます。

$\displaystyle{(12)\quad X=\left[\begin{array}{c} X_{11},\cdots,X_{1m}\\ \vdots\\ %X_{i1},\cdots,X_{im}\\ %\vdots\\ X_{n1},\cdots,X_{nm}\\ \end{array}\right] }$

ここで、 $n$ 回の観測が行われ、各観測では $m$ 個のデータを得るものとしています。これから

$\displaystyle{(13)\quad \bar{X}=\left[\begin{array}{c} 1\\ \vdots\\ 1 \end{array}\right] \left[\begin{array}{c} \displaystyle{\frac{1}{n}\sum_{i=1}^nX_{i1},\cdots,\frac{1}{n}\sum_{i=1}^nX_{im}} \end{array}\right] }$

を引いた行列 $B$ の特異値分解を

$\displaystyle{(14)\quad B=X-\bar{X}=U\Sigma V^H }$

とします。このとき、共分散行列は

$\displaystyle{(15)\quad C=\frac{1}{n-1}BB^H=\frac{1}{n-1}V\Sigma^2 V^H\Rightarrow D=\frac{1}{n-1}\Sigma^2 }$

のように表されます。すなわち共分散行列の主成分ベクトルと、その方向に沿う不変分散値が、行列 $B$ の特異値分解から得られることに留意します。

補遺：出力FB型追従SMI制御

投稿日時: 2021年11月2日投稿者: cacsd

補遺：出力FB型追従SMI制御…Homework

[1] 定値外乱を受ける制御対象

$\displaystyle{(1)\quad \begin{array}{l} \dot{x}(t)=Ax(t)+Bu(t)+w\\ y(t)=Cx(t)\\ (x(t)\in{\rm\bf R}^n, u(t)\in{\rm\bf R}^m, y(t)\in{\rm\bf R}^p, w\in{\rm\bf R}^n, m=p) \end{array} }$

の出力を、次のコマンド（定値目標）

$\displaystyle{(2)\quad r\in{\rm\bf R}^m %\begin{array}{l} %\dot{r}(t)=\Gamma(r(t)-r_c)\\ %(r(t),r_c\in{\rm\bf R}^m) %\end{array} %\Gammaは安定行列 }$

に追従させることを考えます。そのために、積分動作

$\displaystyle{(3)\quad \begin{array}{l} \dot{x}_I(t)=y(t)-r\\ %\dot{x}_I(t)=r(t)-y(t)\\ (x_I(t)\in{\rm\bf R}^m) \end{array} }$

を考え、次の拡大系を構成します。

$\displaystyle{(4)\quad %\underbrace{ \left[\begin{array}{c} \dot{x}(t) \\ \dot{x}_I(t) \end{array}\right] %}_{\dot{x}_E(t)} = %\underbrace{ \left[\begin{array}{cc} A & 0 \\ C & 0 \end{array}\right] %}_{A_E} %\underbrace{ \left[\begin{array}{c} x(t) \\ x_I(t) \end{array}\right] %}_{x_E(t)} + %\underbrace{ \left[\begin{array}{c} B \\ 0 \end{array}\right] %}_{B_E} u(t) + %\underbrace{ \left[\begin{array}{c} w \\ -r \end{array}\right] %}_{w_E} }$

定常状態では

$\displaystyle{(5)\quad \left[\begin{array}{c} 0 \\ 0 \end{array}\right] = \left[\begin{array}{cc} A & 0 \\ C & 0 \end{array}\right] \left[\begin{array}{c} x_\infty \\ x_{I\infty} \end{array}\right] + \left[\begin{array}{c} B \\ 0 \end{array}\right] u_\infty + \left[\begin{array}{c} w \\ -r \end{array}\right] }$

を得ます( $x_{\infty}$ , $x_{I\infty}$ , $u_{\infty}$ は定数ベクトル)。

●さて、関係式

$\displaystyle{(10)\quad \underbrace{ \left[\begin{array}{c} {\dot x}(t) \\ y(t)-r \end{array}\right] }_{x_{E2}(t)} = \underbrace{ \left[\begin{array}{cc} A & B \\ C & 0 \end{array}\right] }_{S_E} \underbrace{ \left[\begin{array}{c} x(t)-x_\infty \\ u(t)-u_\infty \end{array}\right] }_{x_{E3}(t)} }$

を用いて、偏差系

$\displaystyle{(15)\quad \begin{array}{l} \boxed{ \underbrace{ \left[\begin{array}{c} \dot x_1(t)\\ \dot x_2(t) \end{array}\right] }_{\dot{x}_{E3}(t)} = \underbrace{ \left[\begin{array}{cc} A & B \\ 0 & 0 \end{array}\right] }_{A_{E3}} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x_{E3}(t)} + \underbrace{ \left[\begin{array}{c} 0\\ I_m \end{array}\right] }_{B_{E3}} {\dot u}(t)}\\ (x_1(t)=x(t)-x_\infty, x_2(t)=u(t)-u_\infty) \end{array} }$

を得ます。これはSM標準形となっていることに注意して、SM制御系を設計します。

●スイッチング関数として、次式を考えます。

$\displaystyle{(16)\quad s(t)= \underbrace{ \left[\begin{array}{cc} S_1C & S_2 \\ \end{array}\right] }_{S} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x_{E3}(t)} = \underbrace{S_2 \left[\begin{array}{cc} M & I \\ \end{array}\right] }_{S} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x_{E3}(t)} \ (M=S_2^{-1}S_1C) }$

(15)に対して、座標変換

$\displaystyle{(17)\quad \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{x'_{E3}(t)} = \underbrace{ \left[\begin{array}{cc} I & 0 \\ S_1C & S_2 \\ \end{array}\right] }_{T_s} \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x_{E3}(t)}\\ \Leftrightarrow \underbrace{ \left[\begin{array}{c} x_1(t)\\ x_2(t) \end{array}\right] }_{x_{E3}(t)} = \underbrace{ \left[\begin{array}{cc} I & 0 \\ -S_2^{-1}S_1C & S_2^{-1} \\ \end{array}\right] }_{T_s^{-1}} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{x'_{E3}(t)} }$

を行って、次式を得ます。

$\displaystyle{(18)\quad \underbrace{ \left[\begin{array}{c} \dot x_1(t)\\ \dot s(t) \end{array}\right] }_{\dot{x}'_{E3}(t)} = \underbrace{ \left[\begin{array}{cc} \bar{A}_{11} & \bar{A}_{12} \\ \bar{A}_{21} & \bar{A}_{22} \\ \end{array}\right] }_{T_sA_{E3}T_s^{-1}} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{x'_{E3}(t)} + \underbrace{ \left[\begin{array}{cc} 0\\ S_2 \end{array}\right] }_{T_sB_{E3}} {\dot u}(t) }$

$\displaystyle{(19)\quad \begin{array}{l} \bar{A}_{11}=A-BM\quad(M=S_2^{-1}S_1C)\\ \bar{A}_{12}=BS_2^{-1}\\ \bar{A}_{21}=S_1(A-BM)\\ \bar{A}_{22}=S_1BS_2^{-1} \end{array} }$

以下では、 $\bar{A}_{11}$ が安定行列となるようにスイッチング関数が選ばれていると仮定します。

●このとき、スライディングモード制御則（SM制御則、SMC則）

$\displaystyle{(20)\quad {\dot u}(t)=\underbrace{{\dot u}_\ell(t)}_{linear\ control}+\underbrace{{\dot u}_n(t)}_{switching\ component}} }$

を、２次安定性

$\displaystyle{(21)\quad \begin{array}{lll} V(\bar{x})= \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right]^T }_{\bar{x}^T(t)} \underbrace{ \left[\begin{array}{cc} P_1 & 0\\ 0 & P_2 \end{array}\right] }_{P} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{\bar{x}(t)}\\ \Rightarrow \dot{V}(\bar{x})\le - \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right]^T }_{\bar{x}^T(t)} \underbrace{ \left[\begin{array}{cc} Q_1 & 0\\ 0 & I \end{array}\right] }_{Q} \underbrace{ \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] }_{\bar{x}(t)} \end{array}} }$

すなわち

$\displaystyle{(21')\quad \begin{array}{lll} V(\bar{x})= \underbrace{x_1^T(t)P_1x_1(t)}_{V(x_1)}+\underbrace{s^T(t)P_2s(t)}_{V(s)}\\ \Rightarrow \dot{V}(\bar{x})\le -x_1^T(t)Q_1x_1(t)-s^T(t)s(t) \end{array}} }$

が成り立つように決定します（ $P_1>0$ , $P_2>0$ , $Q_1>0$ ）。

可到達性の検討

等価制御は

$\displaystyle{(22)\quad \begin{array}{l} s(t)=0\Rightarrow\dot{s}(t)=0 \Rightarrow 0=\bar{A}_{21}x_1(t)+\bar{A}_{22}s(t)+S_2{\dot u}(t)\\ \Rightarrow {\dot u}_{eq}(t)=-\underbrace{S_2^{-1}}_{(SB_{E3})^{-1}} \underbrace{\left[\begin{array}{cc} \bar{A}_{21} & \bar{A}_{22} \\ \end{array}\right]x'_{E3}(t)}_{SA_{E3}x_{E3}(t)}} \end{array} }$

のように得られます。(20)の第１項 ${\dot u}_\ell$ は、この等価制御をベースして

$\displaystyle{(23)\quad \begin{array}{l} {\dot u}_\ell(t)=-\underbrace{S_2^{-1}}_{(SB_{E3})^{-1}} \underbrace{(\left[\begin{array}{cc} \bar{A}_{21} & \bar{A}_{22} \\ \end{array}\right]x'_{E3}(t)-\Phi s(t))}_{(SA_{E3}-\Phi S)x_{E3}(t)}} \end{array}} }$

のように構成します（ $\Phi$ は安定行列）。このとき閉ループ系は次式で与えられます。

$\displaystyle{(24)\quad \begin{array}{l} \left[\begin{array}{c} \dot{x}_1(t)\\ \dot{s}(t) \end{array}\right] = \left[\begin{array}{cc} \bar{A}_{11} & \bar{A}_{12} \\ 0 & \Phi \\ \end{array}\right] \left[\begin{array}{c} x_1(t)\\ s(t) \end{array}\right] + \left[\begin{array}{c} 0\\ S_2{\dot u}_n(t) \end{array}\right] \end{array }$

すなわち

$\displaystyle{(25)\quad \begin{array}{l} \dot{x}_1(t)=\bar{A}_{11}x_1(t)+\bar{A}_{12}s(t)\\ \dot{s}(t)=\Phi s(t)+S_2{\dot u}_n(t) \end{array} }$

ここで、 $\Phi$ は安定行列なので

$\displaystyle{(26)\quad \begin{array}{l} P_2\Phi+\Phi^TP_2=-I \end{array} }$

を満たす $P_2>0$ を選ぶことができます。これを用いて

$\displaystyle{(27)\quad {\dot u}_n(t)=-\underbrace{S_2^{-1}}_{(SB_{E3})^{-1}}\rho\frac{P_2s(t)}{||P_2s(t)||}} }$

と選びます（ $\rho>0$ は定数）。このとき次式が成り立ちます。

$\displaystyle{(28)\quad \begin{array}{l} \dot{V}(s)=2s^T(t)P_2\dot{s}(t)\\ =2s^T(t)P_2(\Phi s(t)-\rho\frac{P_2s(t)}{||P_2s(t)||})\\ =s^T(t)(P_2\Phi+\Phi^TP_2)s(t)+2s^T(t)P_2(-\rho\frac{P_2s(t)}{||P_2s(t)||})\\ = -||s(t)||^2-2\rho||P_2s(t)||\\ \le -||s(t)||^2 \end{array} }$

スライディングモードの検討

$\bar{A}_{11}$ は安定行列なので

$\displaystyle{(29)\quad \begin{array}{l} P_1\bar{A}_{11}+\bar{A}_{11}^TP_1=-Q_1<0 \end{array} }$

を満たす $P_1>0$ を選ぶことができます。

$\displaystyle{(30)\quad \begin{array}{l} \dot{V}(x_1)=2x_1^T(t)P_1\dot{x}_1(t)\\ =2x_1^T(t)P_1( \bar{A}_{11}x_1(t)+\bar{A}_{12}\underbrace{s(t)}_{0})\\ =x_1^T(t)(P_1\bar{A}_{11}+\bar{A}_{11}^TP_1)x_1(t)\\ =-x_1^T(t)Q_1x_1(t) \end{array} }$

積分動作をもつSMC

上で求めた偏差系E3に対するSMCは次式で与えられました。

$\displaystyle{(31)\quad {\dot u}(t)={\dot u}_\ell(t)+{\dot u}_n(t) }$

$\displaystyle{(32)\quad {\dot u}_\ell(t)=-(SB_{E3})^{-1}(SA_{E3}-\Phi S)x_{E3}(t) }$

$\displaystyle{(33)\quad {\dot u}_n(t)=-(SB_{E3})^{-1}\rho\frac{P_2s(t)}{||P_2s(t)||}}=-(SB_{E3})^{-1}\rho\,{\rm sgn}(P_2Sx_{E3}(t)) }$

これらを積分して、制御対象(1)に対する積分動作をもつSMCを導出します。

$\displaystyle{(34)\quad {u}(t)={u}_\ell(t)+{u}_n(t) }$

まず(32)は(10)を用いて次式のように書けます。

$\displaystyle{(35)\quad {\dot u}_\ell(t)=- \underbrace{ (SB_{E3})^{-1}(SA_{E3}-\Phi S)S_E^{-1} }_{\left[\begin{array}{cc} F & F_I \end{array}\right]} \underbrace{ \left[\begin{array}{c} {\dot x}(t) \\ y(t)-r \end{array}\right] }_{x_{E2}(t)} }$

これを積分して

$\displaystyle{(36)\quad u_\ell(t)=-Fx(t)+F_I\int_0^t(r-y(\tau))d\tau }$

次に(33)は(10)を用いて次式のように書けます。

$\displaystyle{(37)\quad {\dot u}_n(t) =-S_2^{-1}\rho\, {\rm sgn}( \underbrace{ P_2SS_E^{-1} }_{\left[\begin{array}{cc} G & G_I \end{array}\right]} \underbrace{ \left[\begin{array}{c} {\dot x}(t) \\ y(t)-r \end{array}\right] }_{x_{E2}(t)}) }$

これを積分すれば

$\displaystyle{(38)\quad u_n(t)=S_2^{-1}\rho\,\int_0^t\underbrace{{\rm sgn}(-G{\dot x}(\tau)+G_I(r-y(\tau)))}_{0,\pm 1}d\tau }$

水中線状構造物：閉ループ系

投稿日時: 2018年8月15日投稿者: cacsd

● $\displaystyle{\bf Derivation \#6}$
$\displaystyle{ \left[\begin{array}{cc|cc} M_{qq}+M_{qq}^a & M_{qp} & 0 & 0\\ M_{pq} & M_{pp}+M_{pp}^b & 0 & 0\\\hline -c_0 I_N & 0 & I_N & 0\\ 0 & -d_0 I_N & 0 & I_N \end{array}\right] \left[\begin{array}{c} \ddot{q} \\ \ddot{p} \\\hline \ddot{q}_a \\ \ddot{p}_b \end{array}\right]}$
$\displaystyle{+\left[\begin{array}{cc|cc} D_{qq}+D_{qq}^a & D_{qp} & 0 & 0\\ D_{pq} & D_{pp}+D_{qq}^b & 0 & 0\\\hline 0 & 0 & {\cal D}_{qq}^a & 0 \\ 0 & 0 & 0 & {\cal D}_{pp}^b \end{array}\right] \left[\begin{array}{c} \dot{q} \\ \dot{p} \\\hline \dot{q}_a \\ \dot{p}_b \end{array}\right]}$
$\displaystyle{ +\left[\begin{array}{cc|cc} K_{qq} & 0 & -K_{q_aq_a} & -K_{q_ap_b} \\ 0 & K_{pp} & -K_{p_bq_a} & -K_{p_bp_b} \\\hline 0 & 0 & c_2 I_N & 0\\ 0 & 0 & 0 & d_2 I_N \end{array}\right] \left[ \begin{array}{c} q \\ p \\\hline q_a \\ p_b \end{array} \right] = \left[\begin{array}{c} f_1 \\ 0 \\\hline 0\\ 0 \end{array}\right] }$

$\displaystyle{M_{qq}+M_{qq}^a=[m^y_{ij} + \displaystyle\sum_{l,k=1}^NF^y_{ilkj}q_lq_k]_{i,j=1,\cdots,N}+\chi_2I_N}$
$\displaystyle{M_{qp}=[\displaystyle\sum_{l,k=1}^N N^y_{ilkj} q_l p_k]_{i,j=1,\cdots,N}}$
$\displaystyle{M_{pq}=[\displaystyle\sum_{l,k=1}^N N^z_{ilkj} p_l q_k]_{i,j=1,\cdots,N}}$
$\displaystyle{M_{pp}+M_{qq}^b=[m^z_{ij} + \displaystyle\sum_{l,k=1}^NF^z_{ilkj}p_lp_k]_{i,j=1,\cdots,N}+\chi_2I_N}$

for i=1:N, for j=1:N
Mqq(i,j)=m_y(i,j);
Mqp(i,j)=0;
Mpq(i,j)=0;
Mpp(i,j)=m_z(i,j);
for l=1:N, for k=1:N
Mqq(i,j)=Mqq(i,j)+F_y(i,l,k,j)*q(l)*q(k);
Mqp(i,j)=Mqp(i,j)+N_y(i,l,k,j)*q(l)*p(k);
Mpq(i,j)=Mpq(i,j)+N_z(i,l,k,j)*p(l)*q(k);
Mpp(i,j)=Mpp(i,j)+F_z(i,l,k,j)*p(l)*p(k);
end, end, end, end
MM=eye(4*N,4*N);
MM(1:2*N,1:2*N)=[Mqq+chi2*eye(N,N) Mqp;Mpq Mpp+chi2*eye(N,N)];
MM(2*N+1:3*N,1:N)=-c_0*eye(N,N);
MM(3*N+1:4*N,N+1:2*N)=-d_0*eye(N,N);

$\displaystyle{D_{qq}+D_{qq}^a=[c^y_{ij} +\chi_3\sum_{j=1}^N\int_0^1|v-\dot{\eta}|\phi_i\phi_jd\xi + \displaystyle\sum_{l,k=1}^ND^y_{ilkj}q_lq_k+ \sum_{l,k=1}^NE^y_{ilkj}q_l\dot{q}_k]_{i,j=1,\cdots,N}}$
$\displaystyle{D_{qp}=[\displaystyle\sum_{l,k=1}^NL^y_{ilkj}q_lp_k+\displaystyle\sum_{l,k=1}^NM^y_{ilkj}q_l\dot{p}_k]_{i,j=1,\cdots,N}}$
$\displaystyle{D_{pq}=[\displaystyle\sum_{l,k=1}^NL^z_{ilkj}p_lq_k+\displaystyle\sum_{l,k=1}^NM^z_{ilkj}p_l\dot{q}_k]_{i,j=1,\cdots,N} }$
$\displaystyle{D_{pp}+D_{qq}^b=[c^z_{ij} +\chi_3\sum_{j=1}^N\int_0^1|\dot{\eta}|\psi_i\psi_jd\xi + \displaystyle\sum_{l,k=1}^ND^z_{ilkj}p_lp_k+ \displaystyle\sum_{l,k=1}^NE^z_{ilkj}p_l\dot{p}_k]_{i,j=1,\cdots,N} }$
$\displaystyle{{\cal D}_{qq}^a=[c_1\displaystyle\sum_{l,k=1}^NS^y_{ijkl}q_{ak}q_{al}]_{i,j=1,\cdots,N}-c_1I_N}$
$\displaystyle{{\cal D}_{pp}^b=[d_1\displaystyle\sum_{l,k=1}^NS^z_{ijkl}p_{bk}p_{bl}]_{i,j=1,\cdots,N}-d_1I_N}$

for i=1:N, for j=1:N
Dqq(i,j)=c_y(i,j)+chi3*c_simp*(abs(V-dy).*f_y(:,i).*f_y(:,j));
Dqp(i,j)=0;
Dpq(i,j)=0;
Dpp(i,j)=c_z(i,j)+chi3*c_simp*(abs(dy).*f_z(:,i).*f_z(:,j));
Dqqa(i,j)=0;
Dppb(i,j)=0;
for l=1:N, for k=1:N
Dqq(i,j)=Dqq(i,j)+D_y(i,l,k,j)*q(l)*q(k)+E_y(i,l,k,j)*q(l)*dq(k);
Dqp(i,j)=Dqp(i,j)+L_y(i,l,k,j)*q(l)*p(k)+M_y(i,l,k,j)*q(l)*dp(k);
Dpq(i,j)=Dpq(i,j)+L_z(i,l,k,j)*p(l)*q(k)+M_z(i,l,k,j)*p(l)*dq(k);
Dpp(i,j)=Dpp(i,j)+D_z(i,l,k,j)*p(l)*p(k)+E_z(i,l,k,j)*p(l)*dp(k);
Dqqa(i,j)=Dqqa(i,j)+c_1*S_y(i,j,k,l)*qa(k)*qa(l);
Dppb(i,j)=Dppb(i,j)+d_1*S_z(i,j,k,l)*pb(k)*pb(l);
end, end, end, end
DD=eye(4*N,4*N);
DD(1:2*N,1:2*N)=[Dqq Dqp;Dpq Dpp];
DD(2*N+1:3*N,2*N+1:3*N)=Dqqa-c_1*eye(N,N);
DD(3*N+1:4*N,3*N+1:4*N)=Dppb-d_1*eye(N,N);

$\displaystyle{K_{qq}=[k^y_{ij} + \displaystyle\sum_{l,k=1}^NB^y_{ilkj}q_kq_l+ \displaystyle\sum_{l,k=1}^NH^y_{ilkj}p_kp_l]_{i,j=1,\cdots,N}}$
$\displaystyle{K_{pp}=[k^z_{ij} + \displaystyle\sum_{l,k=1}^NB^z_{ilkj}p_kp_l+ \displaystyle\sum_{l,k=1}^NH^z_{ilkj}q_kq_l]_{i,j=1,\cdots,N}}$
$\displaystyle{K_{q_aq_a}=[\chi_4 \displaystyle\sum_{l,k=1}^NS^y_{ijkl}q_{vk}q_{vl}]_{i,j=1,\cdots,N}}$
$\displaystyle{K_{q_ap_b}=[-\chi_5 \displaystyle\sum_{l,k=1}^NS^y_{ijkl}q_{vk}\dot{q}_{l}]_{i,j=1,\cdots,N}}$
$\displaystyle{K_{p_bq_a}=[\chi_4 \displaystyle\sum_{l,k=1}^NS^z_{ijkl}p_{vk}p_{vl}]_{i,j=1,\cdots,N}}$
$\displaystyle{K_{p_bp_b}=[\chi_5 \displaystyle\sum_{l,k=1}^NS^z_{ijkl}p_{vk}\dot{p}_{l}]_{i,j=1,\cdots,N}}$

for i=1:N, for j=1:N
Kqq(i,j)=k_y(i,j);
Kpp(i,j)=k_z(i,j);
Kqqa(i,j)=0;
Kqpa(i,j)=0;
Kpqb(i,j)=0;
Kppb(i,j)=0;
for l=1:N, for k=1:N
Kqq(i,j)=Kqq(i,j)+B_y(i,l,k,j)*q(l)*q(k)+H_y(i,l,k,j)*p(l)*p(k);
Kpp(i,j)=Kpp(i,j)+B_z(i,l,k,j)*p(l)*p(k)+H_z(i,l,k,j)*q(l)*q(k);
Kqqa(i,j)=Kqqa(i,j)+chi4*S_y(i,j,k,l)*q(k)*q(l);
Kqpa(i,j)=Kqpa(i,j)-chi5*S_y(i,j,k,l)*q(k)*dp(l);
Kpqb(i,j)=Kpqb(i,j)+chi4*S_z(i,j,k,l)*p(k)*p(l);
Kppb(i,j)=Kppb(i,j)+chi5*S_z(i,j,k,l)*p(k)*dp(l);
end, end, end, end
KK=eye(4*N,4*N);
KK(1:N,1:N)=Kqq;
KK(N+1:2*N,N+1:2*N)=Kpp;
KK(1:2*N,2*N+1:4*N)=-[Kqqa Kqpa;Kpqb Kppb];
KK(2*N+1:3*N,2*N+1:3*N)=c_2*eye(N,N);
KK(3*N+1:4*N,3*N+1:4*N)=d_2*eye(N,N);

$\displaystyle{f_1=[\chi_3\int_0^1|v-\dot{\eta}|\phi_iv_id\xi]_{i=1,\cdots,N}}$

f1=zeros(N,1);
for i=1:N
f1(i)=chi3*c_simp*(abs(V-dy).*f_y(:,i).*V);
end

●for $i,j,k,l=1,\cdots,N$

$\displaystyle{m^y_{ij}=\int_0^1\phi_i\phi_jd\xi+ \Gamma_m\phi_i(1)\phi_j(1) }$
$\displaystyle{m^z_{ij}=\int_0^1\psi_i\psi_jd\xi+ \Gamma_m\psi_i(1)\psi_j(1) }$

for i=1:N, for j=1:N
m_y(i,j)=c_simp*(f_y(:,i).*f_y(:,j))+G_m*f_y(n_d+1,i)*df_y(n_d+1,j);
end, end
m_z=m_y;

$\displaystyle{E^y_{ijk\ell}=\int_0^1 \phi_i\phi'_j\int_0^\xi \phi'_k \phi'_\ell d\xi d\xi +\Gamma_m \phi_i(1)\phi'_j(1) \int_0^1 \phi'_k \phi'_\ell d\xi }$
$\displaystyle{-\int_0^1 \phi_i\phi''_j\int_\xi^1\int_0^\xi \phi'_k \phi'_\ell d\xi d\xi -\Gamma_m \int_0^1\phi_i\phi''_j\int_0^1 \phi'_k \phi'_\ell d\xi d\xi}$
$\displaystyle{E^z_{ijk\ell}=\int_0^1 \psi_i\psi'_j\int_0^\xi \psi'_k \psi'_\ell d\xi d\xi +\Gamma_m \psi_i(1)\psi'_j(1) \int_0^1 \psi'_k \psi'_\ell d\xi }$
$\displaystyle{-\int_0^1 \psi_i\psi''_j\int_\xi^1\int_0^\xi \psi'_k \psi'_\ell d\xi d\xi -\Gamma_m \int_0^1\psi_i\psi''_j\int_0^1 \psi'_k \psi'_\ell d\xi d\xi}$

for i=1:N, for j=1:N, for k=1:N, for l=1:N
E_y(i,j,k,l)=c_simp*(f_y(:,i).*df_y(:,j).*(simp_u*(df_y(:,k).*df_y(:,l))))…
+G_m*f_y(n_d+1,i)*df_y(n_d+1,j)*c_simp*(df_y(:,k).*df_y(:,l))…
-c_simp*(f_y(:,i).*ddf_y(:,j).*(simp_b*(simp_u*(df_y(:,k).*df_y(:,l)))))…
-G_m*c_simp*(f_y(:,i).*ddf_y(:,j).*(simp_u*(df_y(:,k).*df_y(:,l))));
end, end, end, end
E_z=E_y;

$\displaystyle{F^y_{ijk\ell}=\int_0^1 \phi_i\phi'_j\int_0^\xi \phi'_k \phi'_\ell d\xi d\xi +\Gamma_m \phi_i(1)\phi'_j(1) \int_0^1 \phi'_k \phi'_\ell d\xi }$
$\displaystyle{-\int_0^1 \phi_i\phi''_j\int_\xi^1\int_0^\xi \phi'_k \phi'_\ell d\xi d\xi -\Gamma_m \int_0^1\phi_i\phi''_j\int_0^1 \phi'_k \phi'_\ell d\xi d\xi}$
$\displaystyle{F^z_{ijk\ell}=\int_0^1 \psi_i\psi'_j\int_0^\xi \psi'_k \psi'_\ell d\xi d\xi +\Gamma_m \psi_i(1)\psi'_j(1) \int_0^1 \psi'_k \psi'_\ell d\xi }$
$\displaystyle{-\int_0^1 \psi_i\psi''_j\int_\xi^1\int_0^\xi \psi'_k \psi'_\ell d\xi d\xi -\Gamma_m \int_0^1\psi_i\psi''_j\int_0^1 \psi'_k \psi'_\ell d\xi d\xi}$

for i=1:N, for j=1:N, for k=1:N, for l=1:N
F_y(i,j,k,l)=E_y(i,j,k,l);
end, end, end, end
F_z=F_y;

$\displaystyle{M^y_{ijk\ell}= \int_0^1 \phi_i\phi'_j\int_0^\xi \psi'_k \psi'_\ell d\xi d\xi +\Gamma_m \phi_i(1)\phi'_j(1) \int_0^1 \psi'_k \psi'_\ell d\xi }$
$\displaystyle{-\int_0^1 \phi_i\phi''_j\int_\xi^1\int_0^\xi \psi'_k \psi'_\ell d\xi d\xi -\Gamma_m \int_0^1\phi_i\phi''_j\int_0^1 \psi'_k \psi'_\ell d\xi d\xi }$
$\displaystyle{M^z_{ijk\ell}= \int_0^1 \psi_i\psi'_j\int_0^\xi \phi'_k \phi'_\ell d\xi d\xi +\Gamma_m \psi_i(1)\psi'_j(1) \int_0^1 \phi'_k \phi'_\ell d\xi }$
$\displaystyle{-\int_0^1 \psi_i\psi''_j\int_\xi^1\int_0^\xi \phi'_k \phi'_\ell d\xi d\xi -\Gamma_m \int_0^1\psi_i\psi''_j\int_0^1 \phi'_k \phi'_\ell d\xi d\xi }$

for i=1:N, for j=1:N, for k=1:N, for l=1:N
M_y(i,j,k,l)=c_simp*(f_y(:,i).*df_y(:,j).*(simp_u*(df_z(:,k).*df_z(:,l))))…
+G_m*f_y(n_d+1,i)*df_y(n_d+1,j)*c_simp*(df_z(:,k).*df_z(:,l))…
-c_simp*(f_y(:,i).*ddf_y(:,j).*(simp_b*(simp_u*(df_z(:,k).*df_z(:,l)))))…
-G_m*c_simp*(f_y(:,i).*ddf_y(:,j).*(simp_u*(df_z(:,k).*df_z(:,l))));
end, end, end, end
M_z=M_y;

$\displaystyle{N^y_{ijk\ell}= \int_0^1 \phi_i\phi'_j\int_0^\xi \psi'_k \psi'_\ell d\xi d\xi +\Gamma_m \phi_i(1)\phi'_j(1) \int_0^1 \psi'_k \psi'_\ell d\xi }$
$\displaystyle{-\int_0^1 \phi_i\phi''_j\int_\xi^1\int_0^\xi \psi'_k \psi'_\ell d\xi d\xi -\Gamma_m \int_0^1\phi_i\phi''_j\int_0^1 \psi'_k \psi'_\ell d\xi d\xi }$
$\displaystyle{+\int_0^1 (\phi_i \phi_j''\psi_k''\psi_\ell''+\phi_i \phi_j''\psi_k'\psi_\ell'''+3\phi_i \phi_j'\psi_k''\psi_\ell'''+\phi_i \phi_j'^T\psi_k'\psi_\ell'''' )d\xi}$
$\displaystyle{N^z_{ijk\ell}= \int_0^1 \psi_i\psi'_j\int_0^\xi \phi'_k \phi'_\ell d\xi d\xi +\Gamma_m \psi_i(1)\psi'_j(1) \int_0^1 \phi'_k \phi'_\ell d\xi }$
$\displaystyle{-\int_0^1 \psi_i\psi''_j\int_\xi^1\int_0^\xi \phi'_k \phi'_\ell d\xi d\xi -\Gamma_m \int_0^1\psi_i\psi''_j\int_0^1 \phi'_k \phi'_\ell d\xi d\xi }$
$\displaystyle{+\int_0^1 (\psi_i \psi_j''\phi_k''\phi_\ell''+\psi_i \psi_j''\phi_k'\phi_\ell'''+3\psi_i \psi_j'\phi_k''\phi_\ell'''+\psi_i \psi_j'\phi_k'\phi_\ell'''' )d\xi}$

for i=1:N, for j=1:N, for k=1:N, for l=1:N
N_y(i,j,k,l)=c_simp*(f_y(:,i).*df_y(:,j).*(simp_u*(df_z(:,k).*df_z(:,l))))…
+G_m*f_y(n_d+1,i)*df_y(n_d+1,j)*c_simp*(df_z(:,k).*df_z(:,l))…
-c_simp*(f_y(:,i).*ddf_y(:,j).*(simp_b*(simp_u*(df_z(:,k).*df_z(:,l)))))…
-G_m*c_simp*(f_y(:,i).*ddf_y(:,j).*(simp_u*(df_z(:,k).*df_z(:,l))))…
+c_simp*(f_y(:,i).*ddf_y(:,j).*ddf_z(:,k).*ddf_z(:,l)…
+f_y(:,i).*ddf_y(:,j).*df_z(:,k).*dddf_z(:,l)…
+3*f_y(:,i).*df_y(:,j) .*ddf_z(:,k).*dddf_z(:,l)…
+f_y(:,i).*df_y(:,j) .*df_y(:,k).*ddddf_y(:,l));
end, end, end, end
N_z=N_y;

$\displaystyle{S^y_{ijkl}=\int_0^1 \phi_i\phi_j\phi_k\phi_l d\xi}$
$\displaystyle{S^z_{ijkl}=\int_0^1 \psi_i\psi_j\psi_k\psi_l d\xi}$

for i=1:N, for j=1:N, for k=1:N, for l=1:N
S_y(i,j,k,l)=c_simp*(f_y(:,i).*f_y(:,j).*f_y(:,k).*f_y(:,l));
end, end, end, end
S_z=S_y;

$\displaystyle{c^y_{ij}= \int_0^1 2 \sqrt{\beta}{u} \phi_i \phi'_j d\xi }$
$\displaystyle{c^z_{ij}= \int_0^1 2 \sqrt{\beta}{u} \psi_i \psi'_j d\xi }$

for i=1:N, for j=1:N
c_y(i,j)=2*sqrt(beta)*c_simp*(U.*f_y(:,i).*df_y(:,j));
end, end
c_z=c_y;

$\displaystyle{k^y_{ij}= \int_0^1\gamma\phi_i\phi'_jd\xi+\gamma\Gamma_p\phi_i(1)\phi'_j(1) }$
$\displaystyle{- \int_0^1 (\gamma(1-\xi+\Gamma_p)}$ $\displaystyle{-\sqrt{\beta}{\dot{u}}(1-\xi) }$ $\displaystyle{) \phi_i \phi_j''d\xi }$ $\displaystyle{+\int_0^1 {u^2} \phi_i\phi''_j d\xi }$ $\displaystyle{+ \int_0^1 \phi_i \phi''''_j d\xi}$
$\displaystyle{k^z_{ij}= \int_0^1\gamma\psi_i\psi'_jd\xi+\gamma\Gamma_p\psi_i(1)\psi'_j(1) }$
$\displaystyle{- \int_0^1 (\gamma(1-\xi+\Gamma_p)}$ $\displaystyle{-\sqrt{\beta}{\dot{u}}(1-\xi) }$ $\displaystyle{) \psi_i \psi_j''d\xi }$ $\displaystyle{+\int_0^1 {u^2} \psi_i\psi''_j d\xi }$ $\displaystyle{+ \int_0^1 \psi_i \psi''''_j d\xi}$

for i=1:N, for j=1:N
k_y_0(i,j)=c_simp*(gamma*f_y(:,i).*df_y(:,j))+gamma*G_p*f_y(n_d+1,i).*df_y(n_d+1,j)…
-c_simp*(gamma*(1-xi+G_p))+c_simp*(f_y(:,i).*ddddf_y(:,j));
end, end

for i=1:N, for j=1:N
k_y(i,j)= k_y_0(i,j)+sqrt(beta)*c_simp*(dU.*(1-xi).*f_y(:,i).*ddf_y(:,j))…
+c_simp*(U.^2.*f_y(:,i).*ddf_y(:,j));
end, end
k_z=k_y;

$\displaystyle{B^y_{ijk\ell}={1\over 2}\gamma \int_0^1 \phi_i\phi'_j\phi'_k \phi'_\ell d\xi +{1\over 2}\gamma\Gamma_p \phi_i(1)\phi'_j(1) \phi'_k(1) \phi'_\ell(1)}$
$\displaystyle{-{3\over 2} \int_0^1 (\gamma(1-\xi+\Gamma_p)}$ $\displaystyle{-\sqrt{\beta}{\dot{u}}(1-\xi) }$ $\displaystyle{) \phi_i \phi_j'' \phi'_k \phi'_\ell d\xi }$
$\displaystyle{+\int_0^1 \sqrt{\beta}{\dot{u}}\phi_i \phi_j'' \int_\xi^1 \phi'_k \phi'_\ell d\xi d\xi +\int_0^1 {u^2} ( \phi_i\phi_j' \phi'_k \phi''_\ell- \phi_i\phi''_j\int_\xi^1 \phi'_k \phi''_\ell d\xi) d\xi }$
$\displaystyle{+\int_0^1 (\phi_i \phi_j''\phi_k''\phi_\ell''+4\phi_i \phi_j'\phi_k''\phi_\ell'''+\phi_i \phi_j''\phi_k''\phi_\ell'''' )d\xi}$
$\displaystyle{B^z_{ijk\ell}={1\over 2}\gamma \int_0^1 \psi_i\psi'_j\psi'_k \psi'_\ell d\xi +{1\over 2}\gamma\Gamma_p \psi_i(1)\psi'_j(1) \psi'_k(1) \psi'_\ell(1)}$
$\displaystyle{-{3\over 2} \int_0^1 (\gamma(1-\xi+\Gamma_p)}$ $\displaystyle{-\sqrt{\beta}{\dot{u}}(1-\xi) }$ $\displaystyle{) \psi_i \psi_j'' \psi'_k \psi'_\ell d\xi }$
$\displaystyle{+\int_0^1 \sqrt{\beta}{\dot{u}}\psi_i \psi_j'' \int_\xi^1 \psi'_k \psi'_\ell d\xi d\xi +\int_0^1 {u^2} ( \psi_i\psi_j' \psi'_k \psi''_\ell- \psi_i\psi''_j\int_\xi^1 \psi'_k \psi''_\ell d\xi) d\xi }$
$\displaystyle{+\int_0^1 (\psi_i \psi_j''\psi_k''\psi_\ell''+4\psi_i \psi_j'\psi_k''\psi_\ell'''+\psi_i \psi_j''\psi_k''\psi_\ell'''' )d\xi}$

for i=1:N, for j=1:N, for k=1:N, for l=1:N
B_y_0(i,j,k,l)=1/2*gamma*c_simp*(f_y(:,i).*df_y(:,j).*df_y(:,k).*df_y(:,l))…
+1/2*gamma*G_p*f_y(n_d+1,i).*df_y(n_d+1,j).*df_y(n_d+1,k).*df_y(n_d+1,l)…
-3/2*c_simp*(gamma*(1-xi+G_p))…
+c_simp*(f_y(:,i).*ddf_y(:,j).*ddf_y(:,k).*ddf_y(:,l)…
+4*f_y(:,i).*df_y(:,j) .*ddf_y(:,k).*dddf_y(:,l)…
+f_y(:,i).*ddf_y(:,j).*ddf_y(:,k).*ddddf_y(:,l));
end, end, end, end

for i=1:N, for j=1:N, for k=1:N, for l=1:N
B_y(i,j,k,l)= B_y_0(i,j,k,l)…
+3/2*sqrt(beta)*c_simp*(dU.*(1-xi).*f_y(:,i).*ddf_y(:,j).*df_y(:,k).*df_y(:,l))…
+sqrt(beta)*c_simp*(dU.*f_y(:,i).*ddf_y(:,j).*(simp_b*(df_y(:,k).*df_y(:,l))))…
+c_simp*(U.^2.*(f_y(:,i).*df_y(:,j).*df_y(:,k).*ddf_y(:,l)…
-f_y(:,i).*ddf_y(:,j).*(simp_b*(df_y(:,k).*ddf_y(:,l)))));
end, end, end, end
B_z=B_y;

$\displaystyle{H^y_{ijk\ell}={1\over 2}\gamma \int_0^1 \phi_i\phi'_j\psi'_k \psi'_\ell d\xi +{1\over 2}\gamma\Gamma_p \phi_i(1)\phi'_j(1) \psi'_k(1) \psi'_\ell(1)}$
$\displaystyle{-{1\over 2} \int_0^1(\gamma(1-\xi+\Gamma_p)}$ $\displaystyle{-\sqrt{\beta}{\dot{u}}(1-\xi) }$ $\displaystyle{) \phi_i \phi_j'' \psi'_k \psi'_\ell d\xi }$
$\displaystyle{- \int_0^1 (\gamma(1-\xi+\Gamma_p)}$ $\displaystyle{-\sqrt{\beta}{\dot{u}}(1-\xi) }$ $\displaystyle{) \phi_i \phi_j' \psi'_k \psi''_\ell d\xi }$
$\displaystyle{-\int_0^1 \sqrt{\beta}{\dot{u}}\phi_i \phi_j'' \int_\xi^1 \psi'_k \psi'_\ell d\xi d\xi +\int_0^1 {u^2} ( \phi_i\phi_j' \psi'_k \psi''_\ell - \phi_i\phi''_j\int_\xi^1 \psi'_k \psi''_\ell d\xi) d\xi }$
$\displaystyle{H^z_{ijk\ell}={1\over 2}\gamma \int_0^1 \psi_i\psi'_j\phi'_k \phi'_\ell d\xi +{1\over 2}\gamma\Gamma_p \psi_i(1)\psi'_j(1) \phi'_k(1) \phi'_\ell(1)}$
$\displaystyle{-{1\over 2} \int_0^1(\gamma(1-\xi+\Gamma_p)}$ $\displaystyle{-\sqrt{\beta}{\dot{u}}(1-\xi) }$ $\displaystyle{) \psi_i \psi_j'' \phi'_k \phi'_\ell d\xi }$
$\displaystyle{- \int_0^1 (\gamma(1-\xi+\Gamma_p)}$ $\displaystyle{-\sqrt{\beta}{\dot{u}}(1-\xi) }$ $\displaystyle{) \psi_i \psi_j' \phi'_k \phi''_\ell d\xi }$
$\displaystyle{-\int_0^1 \sqrt{\beta}{\dot{u}}\psi_i \psi_j'' \int_\xi^1 \phi'_k \phi'_\ell d\xi d\xi +\int_0^1 {u^2} ( \psi_i\psi_j' \phi'_k \phi''_\ell - \psi_i\psi''_j\int_\xi^1 \phi'_k \phi''_\ell d\xi) d\xi }$

for i=1:N, for j=1:N, for k=1:N, for l=1:N
H_y_0(i,j,k,l)=1/2*gamma*c_simp*(f_y(:,i).*df_y(:,j).*df_y(:,k).*df_y(:,l))…
+1/2*gamma*G_p*f_y(n_d+1,i).*df_y(n_d+1,j).*df_y(n_d+1,k).*df_y(n_d+1,l)…
-1/2*c_simp*(gamma*(1-xi+G_p).*f_y(:,i).*ddf_y(:,j).*df_y(:,k).*df_y(:,l))…
-c_simp*(gamma*(1-xi+G_p).*f_y(:,i).*df_y(:,j).*df_z(:,k).*ddf_z(:,l));
end, end, end, end

for i=1:N, for j=1:N, for k=1:N, for l=1:N
H_y(i,j,k,l)= H_y_0(i,j,k,l)…
+1/2*sqrt(beta)*c_simp*(dU.*(1-xi).*f_y(:,i).*ddf_y(:,j).*df_y(:,k).*df_y(:,l)…
+sqrt(beta)*(1-xi).*f_y(:,i).*df_y(:,j).*df_z(:,k).*ddf_z(:,l)…
+sqrt(beta)*f_y(:,i).*ddf_y(:,j).*(simp_b*(df_z(:,k).*df_z(:,l))))…
+c_simp*(U.^2.*(f_y(:,i).*df_y(:,j).*df_z(:,k).*ddf_z(:,l)…
-f_y(:,i).*ddf_y(:,j).*(simp_b*(df_y(:,k).*ddf_y(:,l)))));
end, end, end, end
H_z=H_y;

$\displaystyle{D^y_{ijk\ell}=\int_0^1 2 \sqrt{\beta}{u} ( \phi_i\phi_j' \phi'_k \phi'_\ell- \phi_i\phi''_j\int_\xi^1 \phi'_k \phi'_\ell d\xi) d\xi}$
$\displaystyle{D^z_{ijk\ell}=\int_0^1 2 \sqrt{\beta}{u} ( \psi_i\psi_j' \psi'_k \psi'_\ell- \psi_i\psi''_j\int_\xi^1 \psi'_k \psi'_\ell d\xi) d\xi}$

for i=1:N, for j=1:N, for k=1:N, for l=1:N
D_y(i,j,k,l)=2*sqrt(beta)*c_simp*(U.*(f_y(:,i).*df_y(:,j).*df_y(:,k).*df_y(:,l)…
-f_y(:,i).*ddf_y(:,j).*(simp_b*(df_y(:,k).*df_y(:,l)))));
end, end, end, end
D_z=D_y;

$\displaystyle{L^y_{ijk\ell}=\int_0^1 2 \sqrt{\beta}{u} ( \phi_i\phi_j' \psi'_k \psi'_\ell- \phi_i\phi''_j\int_\xi^1 \psi'_k \psi'_\ell d\xi) d\xi}$
$\displaystyle{L^z_{ijk\ell}=\int_0^1 2 \sqrt{\beta}{u} ( \psi_i\psi_j' \phi'_k \phi'_\ell- \psi_i\psi''_j\int_\xi^1 \phi'_k \phi'_\ell d\xi) d\xi}$

for i=1:N, for j=1:N, for k=1:N, for l=1:N
L_y(i,j,k,l)=2*sqrt(beta)*c_simp*(U.*(f_y(:,i).*df_y(:,j).*df_z(:,k).*df_z(:,l)…
-f_y(:,i).*ddf_y(:,j).*(simp_b*(df_z(:,k).*df_z(:,l)))));
end, end, end, end
L_z=L_y;

水中線状構造物：振動子

投稿日時: 2018年8月15日投稿者: cacsd

●
$\displaystyle{\left[\begin{array}{cc} M_{qq} & M_{qp} \\ M_{pq} & M_{pp} \\ \end{array}\right] \left[\begin{array}{c} \ddot{q} \\ \ddot{p} \end{array}\right] +\left[\begin{array}{cc} D_{qq} & D_{qp} \\ D_{pq} & D_{pp} \\ \end{array}\right] \left[\begin{array}{c} \dot{q} \\ \dot{p} \end{array}\right]}$ $\displaystyle{ +\left[\begin{array}{cc} K_{qq} & K_{qp} \\ K_{pq} & K_{pp} \\ \end{array}\right] \left[ \begin{array}{c} q \\ p \end{array} \right]}$
$\displaystyle{{= \left[\begin{array}{c} \int_0^1\phi(\xi)\frac{L^3}{EI} f_yd\xi \\ \int_0^1\psi(\xi)\frac{L^3}{EI} f_zd\xi \end{array}\right]　}}$

●
$\displaystyle{ V_R(s,t)= \left[\begin{array}{c} V(s,t)-\dot{y}(s,t)\\ -\dot{z}(s,t) \end{array}\right] }$
$\displaystyle{||V_R(s,t)||^2=(V(s,t)-\dot{y}(s,t))^2+\dot{z}^2(s,t)}$

$\displaystyle{ \left[\begin{array}{c} f_y(s,t) \\ f_z(s,t) \end{array}\right] =-C_a\rho_eA_e \left[\begin{array}{c} \ddot{y}(s,t) \\ \ddot{z}(s,t) \end{array}\right] }$
$\displaystyle{ +C_{Ds}\cdot\frac{1}{2}\rho_eD_e||V_R(s,t)|| \left[\begin{array}{cc} V(s,t)-\dot{y}(s,t) \\ -\dot{z}(s,t) \end{array}\right] }$
$\displaystyle{+ \frac{1}{||V_R(s,t)||} \left[\begin{array}{cc} V(s,t)-\dot{y}(s,t) & \dot{z}(s,t)\\ -\dot{z}(s,t)& V(s,t)-\dot{y}(s,t) \end{array}\right] \left[\begin{array}{c} f_{Dv}(s,t) \\ f_{L}(s,t) \end{array}\right] }$
$\displaystyle{ \left[\begin{array}{c} f_{Dv}(s,t) \\ f_{L}(s,t) \end{array}\right] = \frac{1}{2}\rho_eD_e||V_R(s,t)||^2 \left[\begin{array}{c} \frac{1}{2}C_{D0}a(s,t) \\ \frac{1}{2}C_{L0}b(s,t) \end{array}\right] }$

$V(s,t)>>\sqrt{\dot{y}^2(s,t)+\dot{z}^2(s,t)}$
$\displaystyle{||V_R(s,t)||^2=(V(s,t)-\dot{y}(s,t))^2+\dot{z}^2(s,t)\simeq V^2(s,t)}$
$\displaystyle{ \frac{1}{||V_R(s,t)||} \left[\begin{array}{cc} V(s,t)-\dot{y}(s,t) & \dot{z}(s,t)\\ -\dot{z}(s,t)& V(s,t)-\dot{y}(s,t) \end{array}\right] \simeq \left[\begin{array}{cc} 1 & \frac{\dot{z}(s,t)}{V(s,t)}\\ -\frac{\dot{z}(s,t)}{V(s,t)} & 1 \end{array}\right] }$

$\displaystyle{ \left[\begin{array}{c} f_y(s,t) \\ f_z(s,t) \end{array}\right] =-C_a\rho_eA_e \left[\begin{array}{c} \ddot{y}(s,t) \\ \ddot{z}(s,t) \end{array}\right] }$
$\displaystyle{+C_{Ds}\cdot\frac{1}{2}\rho_eD_e \left[\begin{array}{c} |V(s,t)-\dot{y}(s,t)|(V(s,t)-\dot{y}(s,t)) \\ -|\dot{z}(s,t)|\dot{z}(s,t) \end{array}\right] }$
$\displaystyle{+ \left[\begin{array}{cc} 1 & \frac{\dot{z}(s,t)}{V(s,t)}\\ -\frac{\dot{z}(s,t)}{V(s,t)} & 1 \end{array}\right] \left[\begin{array}{c} f_{Dv}(s,t) \\ f_{L}(s,t) \end{array}\right] }$
$\displaystyle{ \left[\begin{array}{c} f_{Dv}(s,t) \\ f_{L}(s,t) \end{array}\right] = \frac{1}{2}\rho_eD_eV^2(s,t) \left[\begin{array}{c} \frac{1}{2}C_{D0}a(s,t) \\ \frac{1}{2}C_{L0}b(s,t) \end{array}\right] }$

●
$\displaystyle{\ddot{a}(s,t)+2\epsilon_y\Omega_f(a(s,t)^2-1)\dot{a}(s,t)+4\Omega_f^2a(s,t)=\frac{A_y}{D}\ddot{y}(s,t)}$
$\displaystyle{\ddot{b}(s,t)+ \epsilon_z\Omega_f(b(s,t)^2-1)\dot{b}(s,t)+ \Omega_f^2b(s,t)=\frac{A_z}{D}\ddot{z}(s,t)}$

$\displaystyle{\alpha(\xi,\tau) =\displaystyle\sum_{i=1}^\infty q_{ai}(\tau)\phi_i(\xi)\simeq q_a^T(\tau)\phi(\xi)}$
$\displaystyle{ q_a(\tau)= \left[\begin{array}{c} q_{a1}(\tau) \\ \vdots\\ q_{aN}(\tau) \\ \end{array}\right],\ \phi(\xi)= \left[\begin{array}{c} \phi_1(\xi) \\ \vdots\\ \phi_N(\xi) \\ \end{array}\right] }$

$\displaystyle{\beta(\xi,\tau)=\displaystyle\sum_{i=1}^\infty p_{bi}(\tau)\psi_i(\xi)\simeq p_b^T(\tau)\psi(\xi)}$
$\displaystyle{ p_b(\tau)= \left[\begin{array}{c} p_{b1}(\tau) \\ \vdots\\ p_{bN}(\tau) \\ \end{array}\right],\ \psi(\xi)= \left[\begin{array}{c} \psi_1(\xi) \\ \vdots\\ \psi_N(\xi) \\ \end{array}\right] }$

●
$\displaystyle{{(9)_y:f_y(s,t)=-C_a\rho_eA_e(s,t)\ddot{y}(s,t) }}$
$\displaystyle{{+\frac{1}{2}C_{Ds}\rho_eD_e|V(s,t)-\dot{y}(s,t)|(V(s,t)-\dot{y}(s,t)) }}$
$\displaystyle{{+\frac{1}{4}C_{D0}\rho_eD_eV^2(s,t)a(s,t)+\frac{1}{4}C_{L0}\rho_eD_eV(s,t)\dot{z}(s,t)b(s,t) }}$

$\displaystyle{{(9)_\eta:\frac{L^3}{EI}f_y(s,t) =-\underbrace{C_a\frac{\rho_eA_e}{m+m_i}}_{\chi_2}\underbrace{\frac{d^2\frac{y}{L}}{d\frac{t}{L^4{\frac{m+m_i}{EI}}}}}_{\ddot{\eta}} }}$
$\displaystyle{{+\underbrace{\frac{1}{2}C_{Ds}\frac{\rho_eD_e}{m+m_i}}_{\chi_3} |\underbrace{\frac{V}{\frac{1}{L}\sqrt{\frac{EI}{m+m_i}}}}_{v}- \underbrace{\frac{d\frac{y}{L}}{d\frac{t}{L^2\sqrt{\frac{m+m_i}{EI}}}}}_{\dot{\eta}}| (\underbrace{\frac{V}{\frac{1}{L}\sqrt{\frac{EI}{m+m_i}}}}_{v}- \underbrace{\frac{d\frac{y}{L}}{d\frac{t}{L^2\sqrt{\frac{m+m_i}{EI}}}}}_{\dot{\eta}}) }}$
$\displaystyle{{+\underbrace{\frac{1}{4}C_{D0}\frac{\rho_eD_e}{m+m_i}}_{\chi_4}\underbrace{\frac{V^2}{\frac{1}{L^2}\frac{EI}{m+m_i}}}_{v^2}{a} +\underbrace{\frac{1}{4}C_{L0}\frac{\rho_eD_e}{m+m_i}}_{\chi_5}\underbrace{\frac{V}{\frac{1}{L}\sqrt{\frac{EI}{m+m_i}}}}_{v} \underbrace{\frac{d\frac{z}{L}}{d\frac{t}{L^2\sqrt{\frac{m+m_i}{EI}}}}}_{\dot{\zeta}}{b} }}$
$\displaystyle{{=-\chi_2\ddot{\eta}+\chi_3|v-\dot{\eta}|(v-\dot{\eta})+\chi_4v^2\alpha+\chi_5v\dot{\eta}\beta }}$
$\displaystyle{{=-\chi_2\ddot{\eta}-\chi_3|v-\dot{\eta}|\dot{\eta}+\chi_3|v-\dot{\eta}|v+\chi_4v^2\alpha+\chi_5v\dot{\eta}\beta }}$

$\displaystyle{{\int_0^1\phi(\xi)\times(9)_\eta d\xi=\int_0^1\phi(-\chi_2\ddot{\eta}-\chi_3|v-\dot{\eta}|\dot{\eta}+\chi_3|v-\dot{\eta}|v+\chi_4v^2\alpha+\chi_5v\dot{\eta}\beta)d\xi }}$
$\displaystyle{{=\int_0^1\phi(-\chi_2\underbrace{\phi^T\ddot{q}}_{\ddot{\eta}}-\chi_3|v-\dot{\eta}|\underbrace{\phi^T\dot{q}}_{\dot{\eta}} +\chi_3|v-\dot{\eta}|v+\chi_4\underbrace{(\phi^Tq_v)^2}_{v^2}\underbrace{\phi^T{q}_a}_{\alpha} +\chi_5\underbrace{\phi^Tq_v}_{v}\underbrace{\phi^T\dot{q}}_{\dot{\eta}}\underbrace{\psi^T{p}_b}_{\beta})d\xi }}$
$\displaystyle{{=-\chi_2\int_0^1\phi\phi^T\ddot{q}d\xi-\chi_3\int_0^1|v-\dot{\eta}|\phi\phi^T\dot{q}d\xi+\chi_3\int_0^1\phi|v-\dot{\eta}|vd\xi }}$
$\displaystyle{{+\chi_4\int_0^1 q_v^T\phi\phi^Tq_v\phi\phi^T{q}_ad\xi+\chi_5\int_0^1 q_v^T\phi\phi^T\dot{q}\phi\psi^T{p}_bd\xi }}$

for $i=1,\cdots,N$
$\displaystyle{{\int_0^1\phi_i(\xi)\times(9)_\eta d\xi }}$
$\displaystyle{{= -\chi_2\underbrace{\sum_{j=1}^N \int_0^1 \phi_i \phi_j d\xi\times \ddot{q}_j}_{\ddot{q}_i} -\chi_3\sum_{j=1}^N\int_0^1|v-\dot{\eta}|\phi_i\phi_jd\xi\times \dot{q}_j+\chi_3\int_0^1\phi_i|v-\dot{\eta}|vd\xi }}$
$\displaystyle{{+\chi_4\sum_{j,k,l=1}^N \underbrace{\int_0^1 \phi_k\phi_l\phi_i\phi_j d\xi}_{S^y_{ijkl}}\times {q}_{aj}q_{vk}q_{vl} +\chi_5\sum_{j,k,l=1}^N \underbrace{\int_0^1 \phi_k\phi_l\phi_i \psi_j d\xi}_{S^y_{ijkl}}\times {p}_{bj}q_{vk}\dot{q}_{l} }}$
$\displaystyle{{=-\chi_2\ddot{q}_j-\chi_3\sum_{j=1}^N\int_0^1|v-\dot{\eta}|\phi_i\phi_jd\xi\times \dot{q}_j +\chi_3\int_0^1\phi_i|v-\dot{\eta}|vd\xi }}$
$\displaystyle{{+\chi_4\sum_{j,k,l=1}^N {S^y_{ijkl}}\times {q}_{aj}q_{vk}q_{vl}+\chi_5\sum_{j,k,l=1}^N {S^y_{ijkl}}\times {p}_{bj}q_{vk}\dot{q}_{l} }}$

●
$\displaystyle{{(9)_z:f_z(s,t)=-C_a\rho_eA_e(s,t)\ddot{z}(s,t) }}$
$\displaystyle{{-\frac{1}{2}C_{Ds}\rho_eD_e|\dot{z}(s,t)|\dot{z}(s,t) }}$
$\displaystyle{{-\frac{1}{4}C_{D0}\rho_eD_eV(s,t)\dot{z}(s,t)a(s,t)+\frac{1}{4}C_{L0}\rho_eD_eV^2(s,t)b(s,t) }}$

$\displaystyle{{(9)_\zeta:\frac{L^3}{EI}f_y(s,t) =-\underbrace{C_a\frac{\rho_eA_e}{m+m_i}}_{\chi_2}\underbrace{\frac{d^2\frac{z}{L}}{d\frac{t}{L^4{\frac{m+m_i}{EI}}}}}_{\ddot{\zeta}} }}$
$\displaystyle{{-\underbrace{\frac{1}{2}C_{Ds}\frac{\rho_eD_e}{m+m_i}}_{\chi_3}|\underbrace{\frac{d\frac{z}{L}}{d\frac{t}{L^2\sqrt{\frac{m+m_i}{EI}}}}}_{\dot{\zeta}}|\underbrace{\frac{d\frac{z}{L}}{d\frac{t}{L^2\sqrt{\frac{m+m_i}{EI}}}}}_{\dot{\zeta}} }}$
$\displaystyle{{-\underbrace{\frac{1}{4}C_{D0}\frac{\rho_eD_e}{m+m_i}}_{\chi_4}\underbrace{\frac{V}{\frac{1}{L}\sqrt{\frac{EI}{m+m_i}}}}_{v} \underbrace{\frac{d\frac{z}{L}}{d\frac{t}{L^2\sqrt{\frac{m+m_i}{EI}}}}}_{\dot{\zeta}}{a}+\underbrace{\frac{1}{4}C_{L0}\frac{\rho_eD_e}{m+m_i}}_{\chi_5}\underbrace{\frac{V^2}{\frac{1}{L^2}\frac{EI}{m+m_i}}}_{v^2}{b} }}$
$\displaystyle{{=-\chi_2\ddot{\zeta}-\chi_3|\dot{\zeta}|\dot{\zeta}-\chi_4v\dot{\zeta}\alpha+\chi_5v^2\beta }}$

$\displaystyle{{\int_0^1\psi(\xi)\times(9)_\zeta d\xi=\int_0^1\psi(-\chi_2\ddot{\zeta}-\chi_3|\dot{\zeta}|\dot{\zeta}-\chi_4v\dot{\zeta}\alpha+\chi_5v^2\beta)d\xi }}$
$\displaystyle{{=\int_0^1\psi(-\chi_2\underbrace{\psi^T\ddot{p}}_{\ddot{\zeta}}-\chi_3|\dot{\zeta}|\underbrace{\psi^T\dot{p}}_{\dot{\zeta}} -\chi_4\underbrace{\psi^Tp_v}_{v}\underbrace{\psi^T\dot{p}}_{\dot{\zeta}}\underbrace{\phi^T{q}_a}_{\alpha}+\chi_5\underbrace{(\psi^Tp_v)^2}_{v^2}\underbrace{\psi^T{p}_b}_{\beta})d\xi }}$
$\displaystyle{{=\int_0^1(-\chi_2\psi\psi^T\ddot{p}-\chi_3|\dot{\zeta}|\psi^T\dot{p} -\chi_4p_v^T\psi\psi^T\dot{p}\psi\phi^T{q}_a+\chi_5p_v^T\psi\psi^Tp_v\psi\psi^T{p}_b)d\xi }}$
$\displaystyle{{=-\chi_2\int_0^1\psi\psi^T\ddot{p}d\xi-\chi_3\int_0^1|\dot{\zeta}|\psi^T\dot{p}d\xi}}$
$\displaystyle{{-\chi_4\int_0^1p_v^T\psi\psi^T\dot{p}\psi\phi^T{q}_ad\xi+\chi_5\int_0^1p_v^T\psi\psi^Tp_v\psi\psi^T{p}_bd\xi }}$

for $i=1,\cdots,N$
$\displaystyle{{\int_0^1\psi_i(\xi)\times(9)_\zeta d\xi }}$
$\displaystyle{{= -\chi_2 \sum_{j=1}^N\int_0^1 \psi_i\psi_jd\xi\times \ddot{p}_j -\chi_3\sum_{j=1}^N\int_0^1|\dot{\zeta}|\psi_i\psi_jd\xi\times \dot{p}_j }}$
$\displaystyle{{-\chi_4\sum_{j,k,l=1}^N \int_0^1 \psi_k\psi_l\psi_i\phi_j d\xi\times {q}_{aj}p_{vk}\dot{p}_{l} +\chi_5\sum_{j,k,l=1}^N \int_0^1 \psi_k\psi_l\psi_i \psi_j d\xi\times {p}_{bj}p_{vk}p_{vl} }}$
$\displaystyle{{= -\chi_2 \underbrace{\sum_{j=1}^N\int_0^1 \psi_i\psi_jd\xi\times \ddot{p}_j}_{\ddot{p}_i} -\chi_3\sum_{j=1}^N\int_0^1|\dot{\zeta}|\psi_i\psi_jd\xi\times \dot{p}_j }}$
$\displaystyle{{-\chi_4\sum_{j,k,l=1}^N \underbrace{\int_0^1 \psi_k\psi_l\psi_i\phi_j d\xi}_{S^z_{ijkl}}\times {q}_{aj}p_{vk}\dot{p}_{l} +\chi_5\sum_{j,k,l=1}^N \underbrace{\int_0^1 \psi_k\psi_l\psi_i \psi_j d\xi}_{S^z_{ijkl}}\times {p}_{bj}p_{vk}p_{vl} }}$
$\displaystyle{{= -\chi_2 {\ddot{p}_i}-\chi_3\sum_{j=1}^N\int_0^1|\dot{\zeta}|\psi_i\psi_jd\xi\times \dot{p}_j }}$
$\displaystyle{{-\chi_4\sum_{j,k,l=1}^N {S^z_{ijkl}}\times {q}_{aj}p_{vk}\dot{p}_{l}+\chi_5\sum_{j,k,l=1}^N {S^z_{ijkl}}\times {p}_{bj}p_{vk}p_{vl} }}$

●
$\displaystyle{(10)_a:\ddot{a}(s,t)+2\epsilon_y\Omega_f(a(s,t)^2-1)\dot{a}(s,t)+4\Omega_f^2a(s,t)=\frac{A_y}{D}\ddot{y}(s,t)}$
$\displaystyle{(10)_\alpha:\underbrace{\frac{d^2a}{d\frac{t}{L^4{\frac{m+m_i}{EI}}}}}_{\ddot{\alpha}} +2\epsilon_ym^*\underbrace{\frac{\Omega_f}{\frac{m^*}{L^2\sqrt{\frac{m+m_i}{EI}}}}}_{\omega_f} (\underbrace{a^2}_{\alpha^2}-1)\underbrace{\frac{da}{d\frac{t}{L^2\sqrt{\frac{m+m_i}{EI}}}}}_{\dot{\alpha}}}$
$\displaystyle{+4m^{*2}\underbrace{\frac{\Omega_f^2}{\frac{m^{*2}}{L^4{\frac{m+m_i}{EI}}} }}_{\omega_f^2} \underbrace{a}_\alpha={A_y}\frac{L}{D}\underbrace{\frac{d^2\frac{y}{L}}{d\frac{t}{L^4{\frac{m+m_i}{EI}}}}}_{\ddot{\eta}}}$
$\displaystyle{\ddot{\alpha}+\underbrace{2\epsilon_y m^*\omega_f}_{c_1}(\alpha^2-1)\dot{\alpha}+\underbrace{4m^{*2}\omega_f^2}_{c_2}\alpha=\underbrace{A_y\frac{L}{D}}_{c_0}\ddot{\eta}}$

$\displaystyle{\int_0^1\phi(\xi)\times(10)_\alpha d\xi:\int_0^1\phi(\ddot{\alpha}+c_1(\alpha^2-1)\dot{\alpha}+c_2\alpha- c_0\ddot{\eta}) d\xi=0}$
$\displaystyle{\int_0^1\phi(\underbrace{\phi^T\ddot{q}_a}_{\ddot{\alpha}}+c_1(\underbrace{(\phi^Tq_a)^2}_{\alpha^2}-1)\underbrace{\phi^T\dot{q}_a}_{\dot{\alpha}}+c_2\underbrace{\phi^Tq_a}_{\alpha}- c_0\underbrace{\phi^T\ddot{q}}_{\ddot{\eta}}) d\xi=0}$
$\displaystyle{\int_0^1(\phi\phi^T\ddot{q}_a+c_1(q_a^T\phi\phi^Tq_a-1)\phi\phi^T\dot{q}_a+c_2\phi\phi^Tq_a- c_0\phi\phi^T\ddot{q}) d\xi=0}$

for $i=1,\cdots,N$
$\displaystyle{\underbrace{\sum_{j=1}^N\int_0^1 \phi_i\phi_jd\xi\times \ddot{q}_{aj}}_{\ddot{q}_{ai}} +c_1\sum_{j,k,l=1}^N\underbrace{\int_0^1 \phi_i\phi_j\phi_k\phi_l d\xi}_{S^y_{ijkl}}\times \dot{q}_{aj}q_{ak}q_{al} -c_1\underbrace{\sum_{j=1}^N\int_0^1 \phi_i\phi_j d\xi\times \dot{q}_{aj}}_{\dot{q}_{ai}} }$
$\displaystyle{+c_2\underbrace{\sum_{j=1}^N\int_0^1 \phi_i\phi_jd\xi\times {q}_{aj}}_{{q}_{ai}} -c_0\underbrace{\sum_{j=1}^N\int_0^1 \phi_i\phi_jd\xi\times \ddot{q}_{j}}_{\ddot{q}_{i}}=0 }$
$\displaystyle{{\ddot{q}_{ai}} +c_1\sum_{j,k,l=1}^N S^y_{ijkl}\times \dot{q}_{aj}q_{ak}q_{al}-c_1{\dot{q}_{ai}} }$ $\displaystyle{+c_2{{q}_{ai}} -c_0{\ddot{q}_{i}}=0 }$

●
$\displaystyle{(10)_b:\ddot{b}(s,t)+ \epsilon_z\Omega_f(b(s,t)^2-1)\dot{b}(s,t)+ \Omega_f^2b(s,t)=\frac{A_z}{D}\ddot{z}(s,t)}$
$\displaystyle{(10)_\beta:\underbrace{\frac{d^2{b}}{d\frac{t}{L^4{\frac{m+m_i}{EI}}}}}_{\ddot{\beta}} +\epsilon_z m^*\underbrace{\frac{\Omega_f}{\frac{m^*}{L^2\sqrt{\frac{m+m_i}{EI}}}}}_{\omega_f} (\underbrace{{b^2}}_{\beta^2}-{1})\underbrace{\frac{d{b}}{d\frac{m^*}{L^2\sqrt{\frac{m+m_i}{EI}}}}}_{\dot{\beta}}}$
$\displaystyle{+m^{*2}\underbrace{\frac{\Omega_f^2}{\frac{m^{*2}}{L^4{\frac{m+m_i}{EI}}}}}_{\omega_f^2} \underbrace{{b}}_\beta ={A_z}\frac{L}{D}\underbrace{\frac{d^2\frac{z}{L}}{d\frac{t}{L^4{\frac{m+m_i}{EI}}}}}_{\ddot{\zeta}}}$
$\displaystyle{\ddot{\beta}+\underbrace{\epsilon_z m^*\omega_f}_{d_1}(\beta^2-1)\dot{\beta}+\underbrace{m^{*2}\omega_f^2}_{d_2}\beta=\underbrace{A_z\frac{L}{D}}_{d_0}\ddot{\zeta}}$

$\displaystyle{\int_0^1\psi(\xi)\times(10)_\beta d\xi:\int_0^1\psi(\ddot{\beta}+d_1(\beta^2-1)\dot{\beta}+d_2\beta- d_0\ddot{\zeta}) d\xi=0}$
$\displaystyle{\int_0^1\psi(\underbrace{\psi^T\ddot{p}_b}_{\ddot{\beta}}+d_1(\underbrace{(\psi^Tp_b)^2}_{\beta^2}-1)\underbrace{\psi^T\dot{p}_b}_{\dot{\beta}} +d_2\underbrace{\psi^Tq_a}_{\beta}- d_0\underbrace{\psi^T\ddot{q}}_{\ddot{\zeta}}) d\xi=0}$
$\displaystyle{\int_0^1(\psi\psi^T\ddot{p}_b+d_1(p_b^T\psi\psi^Tp_b-1)\psi\psi^T\dot{p}_b+d_2\psi\psi^Tp_b- d_0\psi\psi^T\ddot{p}) d\xi=0}$

for $i=1,\cdots,N$
$\displaystyle{\underbrace{\sum_{j=1}^N\int_0^1 \psi_i\psi_jd\xi\times \ddot{p}_{bj}}_{\ddot{p}_{bi}} +d_1\sum_{j,k,l=1}^N\underbrace{\int_0^1 \psi_i\psi_j\psi_k\psi_l d\xi}_{S^z_{ijkl}}\times \dot{p}_{bj}p_{bk}p_{bl} -d_1\underbrace{\sum_{j=1}^N\int_0^1 \psi_i\psi_j d\xi\times \dot{p}_{bj}}_{\dot{p}_{bi}} }$
$\displaystyle{+d_2\underbrace{\sum_{j=1}^N\int_0^1 \psi_i\psi_jd\xi\times {p}_{bj}}_{{p}_{bi}} -d_0\underbrace{\sum_{j=1}^N\int_0^1 \psi_i\psi_jd\xi\times \ddot{p}_{j}}_{\ddot{p}_{i}}=0 }$
$\displaystyle{{\ddot{p}_{bi}} +d_1\sum_{j,k,l=1}^N S^z_{ijkl}\times \dot{p}_{bj}p_{bk}p_{bl} -d_1{\dot{p}_{bi}} }$ $\displaystyle{+d_2{{p}_{bi}} -d_0{\ddot{p}_{i}}=0 }$

●
$\displaystyle{ \left[\begin{array}{cc} M_{qq}+M_{qq}^a & M_{qp} \\ M_{pq} & M_{pp}+M_{pp}^b \\ \end{array}\right] \left[\begin{array}{c} \ddot{q} \\ \ddot{p} \end{array}\right] +\left[\begin{array}{cc} D_{qq}+D_{qq}^a & D_{qp} \\ D_{pq} & D_{pp}+D_{qq}^b \\ \end{array}\right] \left[\begin{array}{c} \dot{q} \\ \dot{p} \end{array}\right]}$
$\displaystyle{ +\left[\begin{array}{cc} K_{qq} & 0 \\ 0 & K_{pp} \\ \end{array}\right] \left[ \begin{array}{c} q \\ p \end{array} \right] = \left[\begin{array}{c} f_1 \\ 0 \end{array}\right] +\left[\begin{array}{cc} K_{q_aq_a} & 0 \\ 0 & K_{p_bp_b} \\ \end{array}\right] \left[ \begin{array}{c} q_a \\ p_b \end{array} \right]}$
$\displaystyle{ \left[\begin{array}{c} \ddot{q}_a \\ \ddot{p}_b \end{array}\right] +\left[\begin{array}{cc} {\cal D}_{qq}^a & 0 \\ 0 & {\cal D}_{pp}^b \\ \end{array}\right] \left[\begin{array}{c} \dot{q}_a \\ \dot{p}_b \end{array}\right] + \left[\begin{array}{cc} c_2I_N & 0 \\ 0 & d_2I_N \\ \end{array}\right] \left[ \begin{array}{c} q_a \\ p_b \end{array} \right]}$
$\displaystyle{= \left[\begin{array}{cc} c_oI_N & 0 \\ 0 & d_0I_N \\ \end{array}\right] \left[\begin{array}{c} \ddot{q} \\ \ddot{p} \end{array}\right] }$

$\displaystyle{M_{qq}^a=\chi_2I_N}$
$\displaystyle{M_{qq}^b=\chi_2I_N}$
$\displaystyle{D_{qq}^a=[\chi_3\sum_{j=1}^N\int_0^1|v-\dot{\eta}|\phi_i\phi_jd\xi]_{i,j=1,\cdots,N}}$
$\displaystyle{D_{qq}^b=[\chi_3\sum_{j=1}^N\int_0^1|\dot{\eta}|\psi_i\psi_jd\xi]_{i,j=1,\cdots,N} }$
$\displaystyle{f_1=[\chi_3\int_0^1\phi_i|v-\dot{\eta}|vd\xi]_{i=1,\cdots,N}}$
$\displaystyle{K_{q_aq_a}=[\chi_4 \displaystyle\sum_{l,k=1}^NS^y_{ijkl}q_{vk}q_{vl}]_{i,j=1,\cdots,N}}$
$\displaystyle{K_{q_ap_b}=[-\chi_5 \displaystyle\sum_{l,k=1}^NS^y_{ijkl}q_{vk}\dot{q}_{l}]_{i,j=1,\cdots,N}}$
$\displaystyle{K_{p_bq_a}=[\chi_4 \displaystyle\sum_{l,k=1}^NS^z_{ijkl}p_{vk}p_{vl}]_{i,j=1,\cdots,N}}$
$\displaystyle{K_{p_bp_b}=[\chi_5 \displaystyle\sum_{l,k=1}^NS^z_{ijkl}p_{vk}\dot{p}_{l}]_{i,j=1,\cdots,N}}$
$\displaystyle{{\cal D}_{qq}^a=[c_1\displaystyle\sum_{l,k=1}^NS^y_{ijkl}q_{ak}q_{al}]_{i,j=1,\cdots,N}-c_1I_N}$
$\displaystyle{{\cal D}_{pp}^b=[d_1\displaystyle\sum_{l,k=1}^NS^z_{ijkl}p_{bk}p_{bl}]_{i,j=1,\cdots,N}-d_1I_N}$

おわりに

投稿日時: 2017年8月15日投稿者: cacsd

5. おわりに
　当分科会、大学国際戦略分科会（第１分科会）は６回開催された。
　第１回では、当分科会の目的を、次のように定めた。
　「８大学は，『国際戦略本部強化事業』[33]の中で，それぞれの国際戦略本部を設置している。どのような大学国際戦略のもとで，どのような工学系人材育成の方策を計画・実施しているかを調査し，工学系部局における優れた国際展開モデルを模索する。」
また、東北大学准教授・米澤彰純先生をお招きして、「知の大競争時代とは？―理念と現実―」についてご講演をお願いし、次のコメントをいただいた。
　「頭脳循環の議論、とても大切だと思います。ぜひ、八大学それぞれが、それぞれの個性を活かす形で、大学発の戦略を立てていただくことが、学生にとっても喜ばしいことだと思います。」
　第２回では、国際共同教育の意義や実施体制について、文献[34]～[40]に基づいて議論した。
　第３回では、各大学の国際戦略について議論した。
　第４回では、中間報告[41]のアンケート結果（付録参照）について議論した。
　第５回では、オープンイノベーション[17]について議論した。
　第６回では、当分科会からの提言事項を取りまとめた。
　以上の議論に基づいて、工学系部局における優れた国際展開モデルを模索した結果が、本節の内容である。
　なお、第４回開催当日、当分科会が企画して、シンポジウム「理工系学部における国際戦略」を開催した。講師は
　文部科学省顧問・木村孟氏
　政策研究大学院大学・角南篤准教授
　東北大学高等教育開発推進センター・米澤彰純准教授
にお願いした。角南先生と米澤先生には、付録に示すように、中間報告[41]のアンケートにもご回答いただいた。両先生からは次のコメントをいただいた。
　『日本と国際社会をむすび、相互に行き来できる人材をそだてるということは、その分、我が国がいままで経験しなかった大規模な頭脳流出のリスクに直面することを意味します。したがって、国際的に魅力的な教育と研究が今まで以上に目に見える形で行われる必要が生じますね。』（米澤）
　『日本での産学連携は「シリコンバレーモデル」以外に多様な形態があるので、それらを生かした国際戦略を幅広くとらえること。例えば、産業革新機構など日本版SWFの大学による活用など、新たな取り組みを利用する。』（角南）
　社会の仕組みなど変わるはずがないと誰もが思っていたが、2011年になって、ネットの力は社会体制を揺らし始めた。その力は、アナリー・サクセニアン女史の助言
　『日本にとっての課題は、企業や金融機関、教育機関の垣根を取り払うこと』
に沿うような影響を与え得るのだろうか。
（文責・第１分科会主査）

参考文献
[1]「日本の大学『ガラパゴス化』」―中韓との連携通じ脱却を―、日経2010/06/14
[2] AnnaLee Saxenian: The New Argonauts — Regional Advantage in a Global Economy, Harvard University Press, 2006
[3] 酒井泰介：「最新・経済地理学」―グローバル経済と地域の優位性―（[1]の訳書）、日経ＢＰ社、2008
[4]* AnnaLee Saxenian: The New Argonauts, Word into Action, pp.99-110, 2008
[5]* AnnaLee Saxenian: Brain circulation network, 2006
[6]* AnnaLee Saxenian: Global mobility of skill and the competition for
talent,2007
[7] OECD: http://www.oecd.org/home/0,2987,en_2649_201185_1_1_1_1_1,00.html
[8]* OECD:Education at a Glance 2006
[9]* OECD:Education at a Glance 2007
[10]* OECD:Education at a Glance 2008
[11]* OECD:Education at a Glance 2009
[12]* OECD:Education at a Glance 2010
[13] 中国教育部：http://www.spc.jst.go.jp/organization/org_07.html
[14]* 国の人口順リスト：http://ja.wikipedia.org/wiki/
[15]* 国の国内総生産順リスト：http://ja.wikipedia.org/wiki/
[16] 国際化拠点整備事業（グローバル３０）：http://www.jsps.go.jp/j-kokusaika/index.html
[17]* 内閣府：「オープンイノベーション」を再定義する、2010
[18] 文部科学省：平成21年度版「科学技術白書」―世界の大転換期を乗り越える日本発の革新的科学技術を目指して―、日経印刷㈱、2009
[19]「日本人学生、留学生との競争萎縮」―大学国際化「内向き」に拍車―、教育、日経2011/02/28
[20]* 新成長戦略～「元気な日本」復活のシナリオ～、2010、http://www.meti.go.jp/topic/data/growth_strategy/index.html
[21] 世界トップレベル研究拠点プログラム（WPI）：http://www.jsps.go.jp/j-toplevel/index.html
[22]「科学立国のつまづき」―既得権残り政策生かせず―、日経2011/02/27
[23]「中国経済の展望」―起業家精神が原動力に―、経済教室、日経2011/02/07
[24] 安田聡子：「外国人高度人材のグローバル移動とイノベーション」―brain circulation(頭脳循環)の世界的潮流にわが中小企業はどう向き合うか―，中小企業総合研究第６号、pp.21-42 、2007
[25] 産業革新機構：http:// http://www.incj.co.jp/
[26] UM-SJTU: http://umji.sjtu.edu.cn/en/?basic_type=93
[27]「加速する大学の国際化」―留学生争奪、新興国熱く―、日経2010/08/15
[28]「留学生アジアで奪い合う」―ライバルは欧米＋中韓―、日経産業新聞2010/09/21
[29]「留学生誘致は大学改革で」、インタビュー領空侵犯、日経2008/09/22
[30]* 日本学生支援機構：平成２２年度外国人留学生在籍状況調査結果、2010
[31] JUNBA – サンフランシスコ・ベイエリア大学間連携ネットワーク：http://www.junba.org/index_j.html
[32] 九州大学ロバート・ファン／アントレプレナーシップ・センター：http://www.qrec.kyushu-u.ac.jp/
[33] 日本学術振興会：「グローバル社会における大学の国際展開について」～日本の大学も国際化を推進するための提言～、研究環境国際化の手法開発（大学国際戦略強化事業）最終報告書、2010
[34]* 米澤彰純:「世界大競争時代の大学人材養成～国際共同教育の最前線～」、
第1回　序に代えて― 世界の中の日本の大学教育―
[35]* 米澤彰純: 「世界大競争時代の大学人材養成～国際共同教育の最前線～」、
第２回　ダイナミックに変化する中国と日本のトップ大学
[36]* 米澤彰純:「世界大競争時代の大学人材養成～国際共同教育の最前線～」、
第３回　AUN/SEED-Netと東海大学
[37]* 米澤彰純:「世界大競争時代の大学人材養成～国際共同教育の最前線～」、
第４回　慶應義塾大学と延世大学
[38]* 米澤彰純:「世界大競争時代の大学人材養成～国際共同教育の最前線～」、
第５回　豊田工業大学シカゴ校
[39]* 米澤彰純:「世界大競争時代の大学人材養成～国際共同教育の最前線～」、
第６回　早稲田大学―ナンヤン理工大学ダブルMBAプログラム
[40] 渡邊公一郎：九州大学のグローバル30への取り組み、工学・工業教育研究講演会講演論文集、pp.224-225、2010
[41] 梶原宏之：工学系部局の国際戦略の枠組みに関する一考察、工学・工業教育研究講演会講演論文集、pp.226-227、2010

※ []* はネット経由で入手できる電子データまたは当該URLを示す

国際戦略

投稿日時: 2017年8月15日投稿者: cacsd

4 ８大学工学系部局の国際戦略
4.1 環境の変化
　大学国際戦略分科会が設定すべき戦略を定めるためには、まず８大学工学系部局をとりまく環境の変化を的確に認識しなければならない。1節において、大学の国際化の流れの中で、８大学工学系部局は次のような問題に直面していると述べた。
・世界的な人材獲得競争の下で、優秀な工学系留学生をどう確保するか。
・グローバル経済社会における問題の迅速な解決のためにイノベーションを生み出すことのできる博士課程学生の教育をどう行うか。
　ここでは、2節と3節の調査結果に基づいて、これらの問題の背景を分析し、８大学工学系部局が対応すべき環境の変化について述べる。
　第１番目に、世界的な人材獲得競争の背景にあるのは、新興国の経済発展に伴う留学生数の急激な増加である。2節でみたように、留学生数は、現時点（2007年）での約300万人が、約15年で２倍の約600万人に増加すると予測できる。また、現時点（2008年）で、最も多く留学生を受け入れている国は米国でその数は約60万人（構成比18.7%）である。これに対して日本が受け入れた留学生数約12万人（構成比は3.8%）である。3節でみたように、現在のグローバル経済社会の成立に、米国が優れた留学生を数多く引き寄せたことが貢献していることを考えるとき、今後、世界的な人材獲得競争の激しさが増すとの認識を新たにせざるを得ない。ちなみに、中国の留学生受け入れ数は22万人、また中国は留学生を最も多く、しかも偏りなく各国に出しており、その数は約50万人である。人材獲得競争の相手として意識すべきは米国と将来の中国である。当然のことながら、大学の一部局である（８大学）工学系部局は、このような留学生の急激な増加への対応を免れない。
　ところで、「国際化拠点整備事業―グローバル30（G30）」[16]では、600万人の５%、すなわち30万人を確保したいとの意図から、数値目標30万人が立てられたと聞く。その目的は、次のように掲げられている。
　『世界的な人材獲得競争が厳しくなっている状況の下、我国の高等教育の国際競争力の強化及び留学生等に魅力的な水準の教育等を提供するとともに、留学生と切磋琢磨する環境の中で国際的に活躍できる高度な人材の養成を図ることを目的としています。』
　第２番目に、イノベーションの新たな潮流として、内閣府が作成した資料[17]から、次の記述に注目する。
　『グローバル市場での競争の激化、消費者ニーズの早い変化に対応するために従来以上の速いスピードでイノベーションを実現することが求められ、従来型の自前主義の閉鎖的方法ではなく、必要となる研究開発能力、技術的知見、人的資源及び資金を広くオープンな市場から調達し、効率的なイノベーションを目指す、「オープン・イノベーション」が世界の潮流となってきている。』
　これが意味するものは、3節の頭脳循環に関する調査結果からも、理解できるところである。平成21年度の科学技術白書[18]でも、イノベーションのオープン化・グローバル化に触れている。当然のことながら、技術革新を起こすのに最も近い立場にある８大学工学系部局は、このようなオープン・イノベーションへの対応を免れない。
　第３番目に、日本人学生の理工系離れや海外留学を避ける内向き志向が顕著になってきたことが挙げられる。これは、低成長期に入った我国の少子化傾向にあるにも拘わらず、高度成長期の定員規模が今なお維持されているため「大学全入」の状況となり、学生間の競争や切磋琢磨が以前ほどには行われず、学力の低下が目立ってきていることにもよる。このままでは、各産業界への優秀な人材供給が危うくなる。また、東京大学・園田茂人教授は、大学の国際化が日本人学生の内向き志向に拍車をかけているという皮肉な現象を「国際交流のウィンブルドン化」と名付けている[19]。内向き志向は留学試験への合格率が落ちていることも一因と聞く。大学の一部局である（８大学）工学系部局は、このような日本人学生特有の事情への対応を免れない。

4.2 新成長戦略など
【新成長戦略】
　現政権は「新成長戦略」[20]を打ち出しており、上に述べた８大学工学系部局を取り巻く環境の変化への対応策として参考になる部分がある。
　まず「新成長戦略」の15番目に、「『リーディング大学院』構想等による国際競争力強化と人材育成」があり、次の記述がある（図9参照）。
　『拠点形成と集中投資により、我が国の研究開発・人材育成における国際競争力を強化する。すなわち我が国が強みを持つ学問分野を結集したリーディング大学院を構築し、成長分野などで世界をけん引するリーダーとなる博士人材を国際ネットワークの中で養成する。最先端研究施設・設備や支援体制等の環境整備により国内外から優秀な研究者を引き付けて国際頭脳循環の核となる研究拠点や、つくばナノテクアリーナ等世界的な産学官集中連携拠点を形成する。また、「国立研究開発機関（仮称）」制度の検討を進める。特定分野で世界トップ５０に入る研究・教育拠点を１００以上構築し、イノベーション創出環境を整備するとともに、博士課程修了者の完全雇用と社会での活用を実現する。』
　図9に、「ＷＰＩ世界トップレベル研究拠点プログラム」[21]が参照されているが、そのプログラム委員長井村裕夫氏から次のメッセージがある。
　『近年、優れた頭脳の獲得競争が激化してきている中で、我が国が科学技術の力で政界をリードしていくためには、優秀な人材のグロバールな流動の「環」の中に位置づけられ、世界中から人材が集まる開かれた研究拠点を作っていく必要があります。・・・　今後も、これらのＷＰＩ拠点がさらなる発展を遂げ、真に世界の頭脳循環のハブとなるような世界トップレベルの拠点として卓越した研究成果を生み出すことを期待し、・・・』
　ここでいう頭脳循環は、3節の調査結果を参照すると、どちらかというと留学生の増加を担う新興国とではなく、研究先進国との国際交流を意味しているようである。

図9　「新成長戦略」の第15番目の戦略項目

図10　「新成長戦略」の第８番目の戦略項目

　図9において、他に留意すべきは、「リーディング大学院」により、博士課程の更なる強化を求めており、「学位プログラムに基づく」と記述していることである。また、イノベーションの言葉はあるが、起業家精神を育成するような記述は見当たらない。
　さらに、「新成長戦略」の８番目に、「グローバル人材の育成と高度人材の受け入れ拡大」がある（図10参照）。これより、留学生の増加への対応と、日本人学生特有の事情（内向き志向）への対応がカバーされていると言えるが、８大学工学系部局の使命に則した具体的議論が必要である。

【科学立国のつまづき】
　日本経済新聞社編集委員・滝順一氏は、「科学立国のつまづき」と題した記事[22]の中で、次のような指摘を行っている。
　『しかし先行する米国との研究開発力の差は縮まるどころか、生命科学や情報技術の分野では逆に広がった。背後では中国や韓国に脅かされる。・・・』
　『博士の育成策も破綻している。・・・「大学はお金をかけて育てた若手を即戦力としてだけ利用し、使い捨てにしている」・・・そんな大学院の実態は学生に知られ、博士課程の進学率は減り始めている。』
　『産学連携もうまくいっていない。・・・成果主義の波にさらされ・・・「産業貢献を重視すべき工学部でも論文数を競う風潮が定着した」・・・』
　大学毎に事情が違うので、これらを追認するものではないが、このような指摘があることには耳を傾けなければならない。特に、先述のＧ３０の目的にある高等教育とは、８大学工学系部局にとっては博士課程教育を指し、そのキャリアパスの設計が求められているといえる。これは学位プログラムの提案において考慮することもできる。
　また、同記事には、次の指摘も引用されている。
　『一部の政策だけ欧米流をまねても、社会の仕組みが変わらなくては矛盾が生じるばかりだ』（東京大学名誉教授・黒川清氏）
　『人や情報が自由に動くグローバル化が進み、中国などのキャッチアップに有利な環境が生まれてきた』（政策研究大学院・角南篤准教授）
　一番目の指摘はデリケートな論点（既得権に絡む）を含む。また、二番目の指摘に関連すると思われる、次の報告がある。

【中国成長の源泉＝起業家精神】
　東京大学・丸川知雄教授は、次の指摘を行っている[23]。
　『これまで中国は労働力の豊富さを武器に高い成長を続けてきたが、今後は起業家精神の旺盛さこそが、中国成長の源泉となるであろう。』
　ＧＤＰ２位の座を譲り渡した日本では、新成長戦略に起業家精神という言葉が見当たらない（正確にはＰＤＦファイルで検索をかけると１か所だけ）こともあって、当初は耳を疑わざるをえなかった。しかし、丸川教授の報告内容は、3節で調査した、グローバル社会のイノベーションの仕組みを地で行くものであった。

4.3 国際戦略の一つの枠組
　８大学工学系部局が対応すべき環境の変化として、次の３つを挙げた。
　・新興国からの留学生の急激な増加
　・イノベーションのオープン化・グロバール化
　・日本人学生の内向き志向
　これらへの対応の指針は、一応、次のように述べることができるであろう。
　「８大学工学系部局の高等教育（特に博士課程教育）において、日本人学生と留学生とが切磋琢磨する環境を整え、オープン・イノベーションに対応できる高度なグローバル人材の養成を図る。」
　ここで、高等教育とは、８大学工学系部局にとっては博士課程教育であることを是認せざるを得ない。したがって、修士課程教育と学士課程教育は、できるだけ博士課程教育に繋げることが必要となる。
　しかし、上の指針の中にどのような戦略的意図を組み込むのか、また世界的な人材獲得競争の相手と目される米国や将来の中国と、ＷＩＮ－ＷＩＮの関係を築けるかが問題となる。
　3節の調査結果から、次の仮説が立てられていた。
　「21世紀のイノベーションは、シリコンバレーを核にした分散的・補完的グローバル・ネットワークの中から、現代のアルゴノーツ（米国で教育を受けた外国生まれの起業家で頭脳循環の担い手）によって生み出される」
　以下では、これを受け入れて議論を進める。
　そうすると、まず日本の立ち位置を定める必要がある。そのためには、シリコンバレーのＩＴやバイオ分野での実績を考えて、補完的かつ競争優位性のある分野を特定することになる。これは、一般的には、環境・エネルギー分野と言えそうであるが、最終的には、８大学がそれぞれに、または連携した大学（外国の大学との連携も含めて）が特定すべきものである。以下では、特定された競争優位分野と呼ぶ。そして重要なのは、シリコンバレーを核にした分散的・補完的グローバル・ネットワークの中に入れるかどうかである（もちろん入らなくとも良いというスタンスもあり得る）。
　次に検討すべきは、博士課程修了生のキャリアパスの設計である。これは日本人学生と留学生とで異なる。また、どのような学位プログラムを準備するかにも関わる。まず、留学生は「現代の日本版アルゴノーツ」（日本で教育を受けた外国生まれの起業家で頭脳循環の担い手）になることが期待される。日本には起業家精神を育む土壌がないと言われているが、そのような土壌を準備できないとまでは言い切れない。大企業の中にも企業内ベンチャーの試みがあるし、また中小企業こそ起業家精神なしには生き残れない[24]とよく聞く。しかし現代の起業家の活動はもっと迅速・自由である。また起業するテーマは、将来母国に貢献できるものであることが望ましく、学位プログラムの編成に留意すべきである。一方、博士課程を修了した日本人学生の就職は、2/3が企業と大学で、残り1/3が「行方不明」と揶揄されている[22]。博士課程では留学生とともに切磋琢磨するのであるから、その繋がりを生かせることが望ましい。

図11 日本と途上国との間の頭脳循環

　以上から、留学生向けと日本人学生向けの戦略は表裏一体とし、お互いに切磋琢磨する高度職業人養成のための高等教育（博士課程）の整備が、８大学工学系部局の国際戦略の一つの枠組みとして考えられる（図11参照）。

【留学生向け国際戦略】
　「工学分野の基礎教育（学士課程、修士課程）に力を入れている途上国から、相当数の留学生を受け入れる。留学生に、特定された競争優位分野での学位プログラムによる高等教育（博士課程）を通して、博士号を取得させる。彼らには、日本での起業を促す。将来的には、日本との連携を保ちながら、母国での起業を促す。」
【日本人学生向け国際戦略】
　「オープン・イノベーションに対応できる高度グローバル人材養成を希望する日本人学生に、特定された競争優位分野での学位プログラムによる高等教育（博士課程）を通して、博士号を取得させる。（特に、シリコンバレーでのアントレプレナーシップの研修を受けさせる。）彼らには、同じ学位プログラムを修了した留学生とともに、日本での起業を促す。将来的には、帰国した留学生の日本における連携先となる。」
　これら表裏一体の国際戦略の意図は、現代の日本版アルゴノーツを生み出し、日本と途上国がイノベーションの成果を分かち合い、それを世界に広めることと言える。

4.4 検討すべき課題
　８大学工学系部局の国際戦略の一つの枠組みを述べたが、ここでは関連するいくつかの具体的課題を掘り下げてみる。

【博士課程学位プログラムの編成】
　㈱産業革新機構・産業革新委員長・吉川弘之氏からの次のメッセージがある[25]。
　『人類が遭遇している課題の中には気候変動だけでなく、生物多様性の喪失、資源の枯渇、水問題、食糧問題、新病の発生などの諸問題があり、それらが世界的な人口増加のもとで貧困を解消するという社会的問題と重なって、極めて複雑な難問を抱える時代となりました。また、これらの諸問題に立ち向かうための具体的な行動が全世界で求められる時代が近づいています。』
　工学には、人類が直面する諸問題の解決に有用な物作り（解析・設計・生産ツールの提供）を求められている。途上国では生活インフラの整備が将来的な課題であるが、状況に応じて多様な技術革新が必要となる。水問題を取り上げると、たとえば図12に示す文献のように当該分野の研究成果がまとまれば、競争優位性を前面に出した博士課程の学位プログラムの編成ができるであろう。
　重要なのは、キャリアパスを考慮した学位プログラムの編成である。そこでは、当該学位プログラムを修了する日本人学生と留学生が起業し、それが将来留学生の母国に将来母国に貢献できるように配慮することが望ましい。

【優秀な博士課程学生の確保】
　まず、日本人学生については、オープン・イノベーションに対応できる高度グローバル人材養成の意義とそのキャリアパス（アントレプレナーシップ）の魅力について、十分な説明を行うことが肝要である。これにはオープン・イノベーションの実態に触れさせるなど、それなりの工夫が必要である。
　次に、留学生については、工学分野の基礎教育（学士課程、修士課程）に力を入れている途上国から、相当数を受け入れることが肝要である。もし、基礎教育が不十分であれば、８大学工学系部局のカリキュラムの輸出やサテライト・キャンパスの設置が考えられる。たとえば、ミシガン大学は上海交通大学キャンパス内にサテライトキャンパス[26]を開設し、優秀な学生に確保に努めていると聞く（図13参照）。いずれしても米国との人材獲得競争になるので、これまで以上に８大学工学系部局の魅力を高め、広報することが必要がある[27], [28], [29]。
　相当数の留学生の確保が必要としているのは、起業支援のために、同胞のコミュニティまたは同窓会の果たす役割が大きい[2]とされているからである。
　そもそも８大学工学系部局でどの程度の博士課程留学生の受け入れが可能なのだろうか。専攻分野が学際的となっていることもあって、各大学の個別調査の機会を逃しているが、九州大学の工学系部局の場合、留学生数は、2010年11月時点で、学部=74人、修士=52人、博士=152人、研究生=42人である。

図12　水問題に関わる研究成果をまとめた文献の例

図13 上海交通大学キャンパス内のミシガン大学サテライトキャンパス[26]

　日本学生支援機構の2010年5月時点での集計[30]によると、８大学の留学生受け入れ数は、北大=1162人、東北大=1511人、東大=2772人、東工大=1247人、名大=1501人、京大=1530人、阪大=1662人、九大=1713人、合計13098人である。すべての留学生数141774人のうち、工学分野は22567人（構成比15.9%）であるので、単純にこの構成比を掛けると2082人である。国立大学の大学院と学部の留学生数はそれぞれ24355人と10362人であるので、この割合で８大学の大学院の留学生数を推定すると1460人である。博士課程留学生数の割合を3/4とすると(九大の実績から)、８大学工学系部局で受け入れが可能な博士課程留学生の推定数は1095人、１学年当たり365人となる。10年間で3650人である。これを少ないとみるか多いとみるかの判断は別にして、イノベーション人材として教育できるかが問われている。また、図8で示した台湾へ帰国したハイテク移民数の多さを実感せざるを得ない。

【シリコンバレーとのコネクション】
　８大学工学系部局にとって、シリコンバレーはどのような存在で、そのコネクションは強化すべきなのであろうか。言いかえれば、図7のＩＴ・バイオ分野の頭脳循環拠点シリコンバレーを核とするネットワークに、図11の特定された競争優位分野の頭脳循環拠点日本を核とするネットワークを融合すべきかである。
　おそらく二つの考え方があって、一つ目はどのような分野もＩＴには深く関わるので融合すべきという肯定的なもの、また二つ目はＩＴ分野とは競合しないので融合する必要ないという否定的なものである。
　3節の調査結果によれば、シリコンバレーと同様な（ＩＴ・バイオ分野での）頭脳循環拠点を築くことの難しさ、すなわち起業家精神の開放的な生態系（ecosystem）を築くことの難しさが強調されている。イノベーションがＩＴ分野からバイオ分野へ広がっているのであれば、一つの可能性は、特定された競争優位分野が、シリコンバレーを核にした分散的・補完的グローバル・ネットワークの中に入ることである。
　実は、８大学の一部は、次の組織を界して、シリコンバレーとのコネクションをもっている。
　JUNBA – サンフランシスコ・ベイエリア大学間連携ネットワーク[31]
・大阪大学サンフランシスコ教育研究センター
・九州大学カリフォルニアオフィス
・東北大学米国代表事務所
・東京工業大学バークレーオフィス
したがって、シリコンバレーの流儀を理解させるための海外経験を、日本人学生向け国際戦略の中に、組み込むことが考えられる。
　ちなみに、九大では、次の組織を立ち上げ、アントレプレナーシップの教育に乗り出そうとしている。
　九州大学ロバート・ファン／アントレプレナーシップ・センター[32]
ロバート・ファン氏は、九大で学位を取られ、シリコンバレーで起業されている。