?190?4T1?T2?T3??T4?140??T1?4T2. ??T?4T3?T4?1401??T2?T3?4T4?100??在Matlab命令窗口输入以下命令
>> A = [4,-1,-1,0;-1,4,0,-1;-1,0,4,-1;0,-1,-1,4]; b = [190;140;140;100]; >> x = A\\b; x’
Matlab执行后得 ans =
82.9167 70.8333 70.8333 60.4167
可见T1 = 82.9167, T2 = 70.8333, T3 = 70.8333, T4 = 60.4167.
参考文献
陈怀琛, 高淑萍, 杨威, 工程线性代数, 北京: 电子工业出版社, 2007. 页码: 15-16.
Matlab实验题
假定下图中的平板代表一条金属梁的截面, 并忽略垂直于该截面方向上的热传导. 已知平板内部有30个节点, 每个节点的温度近似等于与它相邻的四个节点温度的平均值. 设4条边界上的温度分别等于每位同学学号的后四位的5倍, 例如学号为16308209的同学计算本题时, 选择Tl = 40, Tu = 10, Tr = 0, Td = 45.
Tu Tu Tu Tl T1 T6 T26 Tr Tr Tl T2 T7 T27 Tl T5 T10 T30 Tr Td Td Td 图10 一块平板的温度分布图
(1) 建立可以确定平板内节点温度的线性方程组. (2) 用Matlab软件求解该线性方程组.
(3) 用Matlab中的函数mesh绘制三维平板温度分布图.
9
案例五. CT图像的代数重建问题
X射线透视可以得到3维对象在2维平面上的投影, CT则通过不同角度的X射线得到3维对象的多个2维投影, 并以此重建对象内部的3维图像. 代数重建方法就是从这些2维投影出发, 通过求解超定线性方程组, 获得对象内部3维图像的方法.
图11双层螺旋CT 图12 CT图像
这里我们考虑一个更简单的模型, 从2维图像的1维投影重建原先的2维图像. 一个长方形图像可以用一个横竖均匀划分的离散网格来覆盖, 每个网格对应一个像素, 它是该网格上各点像素的均值. 这样一个图像就可以用一个矩阵表示,其元素就是图像在一点的灰度值(黑白图像). 下面我们以3?3图像为例来说明.
表4 消耗与产出情况 3?3图像 水平方向上
各点的灰度值 的叠加值
x1 = 1 x2 = 0 x3 = 0 x1 + x2 + x3 = 1 x4 = 0 x5 = 0.5 x6 = 0.5 x4 + x5 + x6 = 1 x7 = 0.5 x8 = 0 x9 = 1 x7 + x8 + x9 = 1.5 竖直方向上x1 + x4 + x7 x2 + x5 + x8 x3 + x6 + x9
= 1.5 = 0.5 = 1.5 的叠加值
每个网格中的数字xi代表其灰度值, 范围在[0, 1]内. 0表示白色, 1表示黑色, 0.5表示灰色. 如果我们不知道网格中的数值, 只知道沿竖直方向和水平方向的叠加值, 为了确定网格中的灰度值, 可以建立线性方程组(含有6个方程, 9个未知数)
?x1?x2?x3?1??x4?x5?x6?1 ????x3?x6?x9?1显然该方程组的解是不唯一的, 为了重建图像, 必须增加叠加值. 如我们增加从右上方到左下方的叠加值, 则方程组将增加5个方程
x1 = 1,
x2 + x4 = 0, x3 + x5 + x7 = 1, x6 + x8 = 0.5, x9 = 1,
和上面的6个方程放在一起构成一个含有11个方程, 9个未知数的线性方程组.
10
【模型准备】设3?3图像中第一行3个点的灰度值依次为x1, x2, x3, 第二行3个点的灰度值依次为x4, x5, x6, 第三行3个点的灰度值依次为x7, x8, x9. 沿竖直方向的叠加值依次为1.5, 0.5, 1.5, 沿水平方向的叠加值依次为1, 1, 1.5, 沿右上方到左下方的叠加值依次为1, 0, 1, 0.5, 1. 确定x1, x2, …, x9的值.
【模型建立】由已知条件可得(含有11个方程, 9个未知数的)线性方程组
?x1?x2?x3?1??x4?x5?x6?1 ????x9?1【模型求解】在Matlab命令窗口输入以下命令
>> A = [1,1,1,0,0,0,0,0,0;0,0,0,1,1,1,0,0,0;0,0,0,0,0,0,1,1,1;
1,0,0,1,0,0,1,0,0;0,1,0,0,1,0,0,1,0;0,0,1,0,0,1,0,0,1; 1,0,0,0,0,0,0,0,0;0,1,0,1,0,0,0,0,0;0,0,1,0,1,0,1,0,0; 0,0,0,0,0,1,0,1,0;0,0,0,0,0,0,0,0,1];
>> b = [1;1;1.5;1.5;0.5;1.5;1;0;1;0.5;1]; >> x = A\\b; x’
Matlab执行后得
Warning: Rank deficient, rank = 8 tol = 4.2305e-015. ans =
1.0000 0.0000 0 -0.0000 0.5000 0.5000 0.5000 -0.0000 1.0000 可见上述方程组的解不唯一. 其中的一个特解为
x1 = 1, x2 = 0, x3 = 0, x4 = 0, x5 = 0.5, x6 = 0.5, x7 = 0.5, x8 = 0, x9 = 1.
【模型分析】上述结果表明, 仅有三个方向上的叠加值还不够.可以再增加从左上方到右下方的叠加值. 在实际情况下, 由于测量误差, 上述线性方程组可能是超定的. 这时可以将超定方程组的近似解作为重建的图像数据.
Matlab实验题
给定一个3?3图像的2个方向上的灰度叠加值: 沿左上方到右下方的灰度叠加值依次为0.8, 1.2, 1.7, 0.2, 0.3; 沿右上方到左下方的灰度叠加值依次为0.6, 0.2, 1.6, 1.2, 0.6.
(1) 建立可以确定网格数据的线性方程组, 并用Matlab求解. (2) 将网格数据乘以256, 再取整, 用Matlab绘制该灰度图像.
11
案例六. 平衡结构的梁受力计算
在桥梁、房顶、铁塔等建筑结构中, 涉及到各种各样的梁. 对这些梁进行受力分析是设计师、工程师经常做的事情.
图13埃菲尔铁塔全景 图14 埃菲尔铁塔局部
下面以双杆系统的受力分析为例, 说明如何研究梁上各铰接点处的受力情况. 【模型准备】在图15所示的双杆系统中, 已知杆1重G1 = 200牛顿, 长L1 = 2米, 与水平方向的夹角为?1 = ?/6, 杆2重G2 = 100牛顿, 长L2 = 2米, 与水平方向的夹角为?2 = ?/4. 三个铰接点A, B, C所在平面垂直于水平面. 求杆1, 杆2在铰接点处所受到的力.
C 杆1 A 杆2 ?/6 ?/4 B
图15双杆系统
【模型假设】假设两杆都是均匀的. 在铰接点处的受力情况如图16所示. 【模型建立】对于杆1:
水平方向受到的合力为零, 故N1 = N3,
竖直方向受到的合力为零, 故N2 + N4 = G1,
1以点A为支点的合力矩为零, 故(L1sin?1)N3 + (L1cos?1)N4 = (L1cos?1)G1.
2N4 N3 N2 A 杆1 N1 G1 图16 两杆受力情况
对于杆2类似地有
C C N8 G2 N5 N7 杆2 N6 B
12
1N5 = N7, N6 = N8 + G2, (L2sin?2)N7 = (L2cos?2)N8 + (L2cos?2)G2.
2此外还有N3 = N7, N4 = N8. 于是将上述8个等式联立起来得到关于N1, N2, …, N8的线性方程组:
?N1?N3?0?N?N?G?241 ????N4?N8?0【模型求解】在Matlab命令窗口输入以下命令
>> G1=200; L1=2; theta1=pi/6; G2=100; L2=sqrt(2); theta2=pi/4; >> A = [1,0,-1,0,0,0,0,0;0,1,0,1,0,0,0,0;
0,0,L1*sin(theta1),L1*cos(theta1),0,0,0,0;0,0,0,0,1,0,-1,0; 0,0,0,0,0,1,0,-1;0,0,0,0,0,0,L2*sin(theta2),-L2*cos(theta2); 0,0,1,0,0,0,-1,0;0,0,0,1,0,0,0,-1];
>> b = [0;G1;0.5*L1*cos(theta1)*G1;0;G2;0.5*L2*cos(theta2)*G2;0;0]; >> x = A\\b; x’ Matlab执行后得 ans =
95.0962 154.9038 95.0962 45.0962 95.0962 145.0962 95.0962 45.0962
【模型分析】最后的结果没有出现负值, 说明图16中假设的各个力的方向与事实一致. 如果结果中出现负值, 则说明该力的方向与假设的方向相反.
参考文献
陈怀琛, 高淑萍, 杨威, 工程线性代数, 北京: 电子工业出版社, 2007. 页码: 157- 158.
Matlab实验题
有一个平面结构如下所示, 有13条梁(图中标号的线段)和8个铰接点(图中标号的圈)联结在一起. 其中1号铰接点完全固定, 8号铰接点竖直方向固定, 并在2号, 5号和6号铰接点上, 分别有图示的10吨, 15吨和20吨的负载. 在静平衡的条件下,任何一个铰接点上水平和竖直方向受力都是平衡的. 已知每条斜梁的角度都是45o.
(1) 列出由各铰接点处受力平衡方程构成的线性方程组. (2) 用Matlab软件求解该线性方程组, 确定每条梁受力情况.
13
数学建模案例分析--线性代数建模案例(20例)
