本文介紹《Higher mode effects in hinged wall with BRBs in base-frame structures using distributed parameter models》中關於底部帶BRB的鉸支牆（HWBB）-鉸接框架的模態和內力求解的matlab程序。
共分為四個文件，僅需把這四個文件放到一個文件夾中，在matlab中打開文件夾的路徑，按F5運行即可。
第一個文件為：model.m 用於定義結構的基本情況
function F=spring(x)
E=3.7e+10; %混凝土彈性模量
I=2.86e-3; %混凝土截面慣性矩
h=0.8; %層高
H=4; %結構總高度
m=2.424e+3; %沿結構高度均勻分布的質量
ba=0.163; %牆抗彎剛度和轉動彈簧轉動剛度的相對剛度比
kb=ba*E*I/h; %轉角彈簧轉動剛度
q=kb/E/I; %中間量
l=4.0; %結構總高度
F=-q/x.*cos(x*l).*cosh(x*l)+sin(x*l).*cosh(x*l)-cos(x*l).*sinh(x*l)-q/x;
end
註：這裡最後一個函數的未知量為a，即原文中的公式(7).
第二個文件為：solveequation.m 是用於求解頻率方程，從而求得結構的周期和頻率。
x0=1.8;
x=fsolve(@spring,x0)
註：這裡求出的量實際為a,根據公式(7)可以求出圓頻率。因為結構的模態周期是從大到小，所以頻率是從小到大，即一階頻率最小。因為只有數值解，這裡通過試算得到前幾階頻率的值。建議：將x0的值從0.001增大，會找到前幾階頻率的值。
第三個文件為：reaction.m 求解模態和結構的位移及內力
syms A4 x;
A4=1;
E=3.7e+10;
I=2.86e-3;
ba=0.163;
h=0.8;
kb=ba*E*I/h;
q=kb/E/I;
l=4;
m=1.7e+5;
mt=8.16e+6;
a=0.7373;
w=a^2*(E*I/m)^0.5;
T=6.28/w;
A2=-A4;
A3=A4*(1+cosh(a*l)*cos(a*l)-sinh(a*l)*sin(a*l))/(cosh(a*l)*sin(a*l)-sinh(a*l)*cos(a*l));
A1=A4*(1+cos(a*l)*cosh(a*l)+sin(a*l)*sinh(a*l))/(sin(a*l)*cosh(a*l)-cos(a*l)*sinh(a*l));
%natural mode shape
y=A1*sin(a*x)+A2*cos(a*x)+A3*sinh(a*x)+A4*cosh(a*x);
%
y1=@(x) m*(A1*sin(a*x)+A2*cos(a*x)+A3*sinh(a*x)+A4*cosh(a*x));
lnh=integral(y1,0,l);
z=@(x) m*(A1*sin(a*x)+A2*cos(a*x)+A3*sinh(a*x)+A4*cosh(a*x)).*(A1*sin(a*x)+A2*cos(a*x)+A3*sinh(a*x)+A4*cosh(a*x));
mn=integral(z,0,l);
tn=lnh/mn; %modal particitation factor;
mnx=tn*lnh; %effective modal mass
mnxpercent=mnx/mt;
c=0;
x=0:0.4:4; %用於控制沿結構高度求解量的大小
y=A1*sin(a*x)+A2*cos(a*x)+A3*sinh(a*x)+A4*cosh(a*x);
ust=tn/w/w*(A1*sin(a*x)+A2*cos(a*x)+A3*sinh(a*x)+A4*cosh(a*x)); %靜位移
mw=E*I*tn/w/w*(-A1*a^2*sin(a*x)-A2*a^2*cos(a*x)+A3*a^2*sinh(a*x)+A4*a^2*cosh(a*x)); %牆體彎矩-靜力
vw=-E*I*tn/w/w*(-A1*a^3*cos(a*x)+A2*a^3*sin(a*x)+A3*a^3*cosh(a*x)+A4*a^3*sinh(a*x)) %牆體剪力-靜力
上述求得結構的靜位移、靜態內力之後，再乘以各階周期對應的偽加速度值，可以得到各階周期的動力反應；對於具有稀疏頻譜的結構，可以通過SRSS將各階振型的反應量疊加，即可獲得總的結構反應。