科学计算—理论、方法
及其基于MATLAB的实现与分析
解线性方程组的迭代法线性方程组的理论求解公式
x1-
=（１）
A
b
在应用于实际问题的计算时，通常面临两方面的问题
１、计算过程复杂，
２、不能保证算法的稳定性；
此外，当初始数据（可能）存在误差时，按公式（１）即使求出了“精确解”意义也不大，因此，对于存在初始数据误差、特别是大型的线性方程组求解，需要寻求能达到精度要求的、操作和计算过程相对简单的求解方法。下面将要介绍的迭代法就属于这类方法。
迭代法求解线性方程组的基本思想是
１）不追求“一下子”得到方程组的解，而是在逐步逼近方程组的精确解的迭代过程中获得满足精度要求的近似解，这一点与直接法不同；
２）通过对问题的转化，避免（困难的）矩阵求逆运算。
用迭代法求解线性方程组，首先要把线性方程组写成等价的形式
=
Ax+
?
=（２）
x
b
f
Mx
式（２）的右端称为迭代格式，由迭代格式（２）确定如下的迭代算法：
,2,10
1=???∈?+=+m R x f
Mx x n m m （３）
对于给定的线性方程组，可以写成不同的（无穷多）迭代格式，有意义的（可用的）迭代格式应具有收敛性―生成的解向量序列{}x n 收敛于方程组的解；而好的迭代法应具有较高的收敛速度。
关于迭代法收敛性的两个判别条件：
a 、充分必要条件是：矩阵M 的谱半径
(){}1,,2,1max <==n i M i
i
λρ
b 、充分条件是：矩阵M 的某个算子范数M <1。
设x 是方程组（２）的解，{}m x 是迭代法（３）生成的任一序列，因为
f Mx x +=，f Mx x m m +=+1
所以
()()()0221x x M x x M x x M x x m m m m -==-=-=--- （４）
设11--=?=TJT M J MT T ，其中矩阵J 是矩阵M 的Jordan 标准型，那么容易验证1-=T TJ M
m m
，并且
()[]
()
()()1
0lim ,2,10
lim lim 0
lim lim 010
→+∞
→-+∞
→+∞
→+∞
→M M n
i x x T J T x x M x x m
m m i n m m m m m m ρρλ （５）
此外，因为
(
)0
22
11x x M
x x M x x M x x M x x m m m m m -≤≤-≤-≤-=---- （６）
所以
x x x x M
M m m m m m
m =?=-?=?<+∞
→+∞
→+∞
→lim 0lim 0lim 1 （７）
注：迭代格式（２）所确定的迭代法收敛与否，完全由系数矩阵M 决定，而与常数项f 无关．
常用的迭代法
１、Jacobian 迭代法:
U L D a a a a a a a a a a a a a a a a a a A n n n nn n nn nn n n n n --=?
?
??????????----????????????----????????????=?
?
???
???????=--00000000
0000001112112122112
1
22221
11211
()()?
??=+=?++=???
?
--==----+b D f U L D M b D x U L D x U
L D A b Ax m m 1
1111 （８）
例1 解下面方程组（精确解为T x )1,1,1(*=）.
???
??=++-=+-=++.
14103,53102,
14310321
321321x x x x x x x x x
解 1） 改写成等价形式
???
?
?
?
???--=++=----=--=).314(101),325(101)325(101),314(10121321312321x x x x x x x x x x x
2） 构造迭代公式，即为雅可比迭代公式
???
?
??
???=--=++=--=+++.,2,1,0),314(101),
325(101),314(101)(2)(1)1(3)
(2)(1)1(2)(3)(2)1(1 k x x x x x x x x x k k k k k k k k k 3） 取初始向量T x )0,0,0()0(=，即,0)
0(3)0(2)0(1===x x x 代入上式，求出
4.110
14,5.0105,4.11014)
1(3)1(2)1(1=====
x x x . 再代回公式中，求出
11.1)4.15.0314(10
1
)2(1=-?-=
x ， 2.1)5.034.125(10
1
)2(2=?+?+=
x ， 11.1)2.1311.114(10
1
)2(3=?--=
x . 依次迭代，计算结果如表4-1.
Example intera_j.m
itera_j
2、Gauss-Seidel 迭代法:
()()()()???-=-=?-+-=???
?--==----+b
L D f U L D M b L D Ux L D x U
L D A b Ax n n 1
1
1
11 （9）
根据GS 迭代法（9），可进一步得到
()()()b
D Ux Lx D x b Ux x L D b
L D Ux L D x n n n n n n n 111111
1
1)(-+-++--+++=?+=-?-+-= （10）
即
???
?
??????????+??????? ??????????????---+
??????? ??????????????---=??????? ??-+++--+++b x x x a a a x x x a a a D x x x N n n n n n n N n n n nn n N n n n
2
1
1112
12
1111
121112
111
00
00
00000000 （11）
式（11）表明：Gauss-Seidel 迭代法在计算第k 个迭代值k
n x 1+时，及时地利用了在此步迭代中得到的新的迭代值：i n x 1+，1,,2,1-=k i ，由于第
1+n 步的迭代值通常比第n 步的迭代值更接近方程组的精确解，所以，
在Jacobian 迭代法和Gauss-Seidel 迭代法都收敛的情况下，Gauss-Seidel 迭代法的收敛速度比Jacobian 迭代法的收敛速度快。
例2 解下面方程组（与例1相同，精确解为T x )1,1,1(*=）.
???
??=++-=+-=++.
14103,53102,14310321
321321x x x x x x x x x 解 1） 原方程组改为等价方程组
???
?
?
?
???
--=++=--=).314(101),325(101),314(101213212321x x x x x x x x x .
2） 构造迭代公式，即为高斯-赛德尔迭代公式
???
?
??
???=--=++=--=++++++.,2,1,0),314(101),325(101),314(101)1(2)1(1)1(3)
(2)1(1
)1(2)(3)(2)1(1 k x x x x x x x x x k k k k k k k k k . 3） 取初始向量T x )0,0,0()0(=，即,0)
0(3)0(2)0(1===x x x 代入上式，求出
4.1)00314(10
1
)1(1=-?-=
x ， 78.0)034.125(10
1
)
1(2=?+?+=
x ， 026.1)78.034.114(10
1
)1(3=?--=
x . 迭代计算下去，得表4-2.
Example itera_gs.m
itera_gs
对Gauss-Seidel 迭代法做进一步的研究，式（９）还可以写成如下的形式：
()()()()
()n n J n
S G n n n n n n n n n n n n n n x x L D x x x x L D b D x U L D x b D Lx Ux Lx Lx D x b
Ux x L D b
L D Ux L D x -+=?-+++=?+-++=?+=-?-+-=+-++---+-+-++--+11
_111111111111
1
1)()( （12） 即Gauss-Seidel 迭代法是在Jacobian 迭代法的基础上增加了一个修正项：
()n n x x L D -+-11
并且这个修正项起到了加速的作用，受这一事实的启发，为获得更快的收敛速度，我们对Gauss-Seidel 迭代法再做进一步的研究
()()()()()()
b x D U Lx D x x b x D U Lx Dx Dx b
Ux x L D b
L D Ux L D x n n n n n n n n n n n n +-++=?+-++=?+=-?-+-=+-++++--+1111111
1
1 （13）
式（13）表明Gauss-Seidel 迭代法的第1+n 步迭代值是在第n 步迭代值的基础上增加了一个修正项：
()()b x D U Lx D n n +-++-11 （14）
并且这个修正项使第1+n 步迭代值更接近方程组的精确解，受这一事实的启发，我们希望通过用一个适当的因子乘修正项（14）的办法达到获得更快收敛速度的目的：
()()()()()()[]()
()[]()b
wL D w x wU D w wL D x wb x wU D w x wL D b x D U Lx w Dx Dx b x D U Lx wD x x n n n n n n n n n n n n 11
111111111--+++++-+-++--=?++-=-?+-++=?+-++= （15）
从而得到
3、SOR （Successive Over Relaxation Method ） 迭代法:
()
()[]()x D wL w D wU x w D wL b n n +--=--++-11
11 （16）
其中w 称为松弛因子，由于
()()[]
()()()???
??
?
?????
?------???
??
???????=+--=----nn n n n nn nn n a w wa a w wa wa a w a wa wa a wa a wU D w wL D M 10
101000112211211
1
1
1
2221111
（17） ()()()()∏∏∏====-=-=
?n
k k n
n
k kk
n
k kk
n M w a
a w M 1
1
111det λ （1８）
()()1120111>??>?
>∨<-?∏=M M w w w k
n
k k n
λλ （1９）
所以，为保证SOR 迭代法的收敛性，要求20< 例3 用超松驰迭代法解下面方程组，取松驰因子4.1=ω. ????? ???????=????????????????????????------010121001210012100 124321x x x x . 解 方程组的精确解为：T x )8.0,6.1,4.1,2.1(=.取初值T x )1,1,1,1()0(=，用高斯-赛德尔迭代10次，得T x )798216.0,5964336.1,3955917.1,1966324.1()10(=. 利用SOR 方法，构造迭代公式 ????? ??????-+ =+-++=+-+=+-+ =+++++++). 2(2 ),21(2 ), 2(2 ),21(2 )(4)1(3)(4)1(4)(4)(3)1(2)(3)1(3)(3)(2)1(1)(2)1(2) (2)(1)(1) 1(1k k k k k k k k k k k k k k k k k k x x x x x x x x x x x x x x x x x x ω ω ω ω 与高斯-赛德尔方法相同，初值为T x )1,1,1,1()0(=.迭代计算结果列于表4-3. 表4-3 SOR 迭代的数值结果 k )(1k x )(2k x ) (3k x )(1k x 1 2 3 4 5 1 1 1.343 1.19545 1.203472 1 1.49 1.4753 1.40236 1.4028735 1.7 1.616 1.65095 1.6019815 1.5939059 0.79 0.8152 0.829585 0.789553 0.7999129 由表4-3可知，用超松驰迭代法只迭代了5次，结果与G-S 法迭代10次的结果大体相同，可见，SOR 方法的松驰因子起到了加速收敛的重要作用. Example itera_sor.m itera_sor 关于迭代法收敛性的两个判别条件： a 、充分必要条件是：矩阵M 的谱半径 (){}1,,2,1max <==n i M i i λρ b 、充分条件是：矩阵M 的某个算子范数M <1。 特例：对于下面方程组，证明：雅可比迭代法收敛而高斯-赛德尔迭代法发散. ?? ?? ? ?????=????????????????????-111122*********X X X . 证明 1) 对于雅可比迭代矩阵 ?? ?? ? ?????-----=+=-022101220)(1U L D B J ， J B 的特征方程为 02 21122)det(=-=-λ λλλJ B I . J B 的特征多项式0)(3==λλf ，所以0=λ为J B 的特征根，显然10)(<=J B ρ,因此，由迭代收敛基本定理可知雅可比迭代法收敛. 2) 对于高斯-赛德尔迭代法，其迭代矩阵为 ?? ?? ? ?????---=-=--200320220)(1U L D B S G ， S G B -的特征方程为 0)2(2 003202 2)det(2=+=+-+-=--λλλλλλS G B I . 显然，12)(>=-S G B ρ，因此高斯-赛德尔迭代法发散. 关于收敛性重要结论： 设线性方程组b Ax =. （1）A 为严格对角占优矩阵，则解b Ax =的雅可比迭代过程和高斯 -赛德尔迭代过程均收敛. （2）若A 是对称正定方阵，则解b Ax =的高斯-赛德尔迭代过程收 敛. （3）A 为严格对角占优矩阵，10≤<ω，则解b Ax =的SOR 迭代过程收敛. 解线性方程组的迭代法 1.rs里查森迭代法求线性方程组Ax=b的解 function[x,n]=rs(A,b,x0,eps,M) if(nargin==3) eps=1.0e-6;%eps表示迭代精度 M=10000;%M表示迭代步数的限制值elseif(nargin==4) M=10000; end I=eye(size(A)); n=0; x=x0; tol=1; %迭代过程 while(tol>eps) x=(I-A)*x0+b; n=n+1;%n为最终求出解时的迭代步数tol=norm(x-x0); x0=x; if(n>=M) disp('Warning:迭代次数太多，可能不收敛！'); return; end end 2.crs里查森参数迭代法求线性方程组Ax=b的解 function[x,n]=crs(A,b,x0,w,eps,M) if(nargin==4) eps=1.0e-6;%eps表示迭代精度 M=10000;%M表示迭代步数的限制值 elseif(nargin==5) M=10000; end I=eye(size(A)); n=0; x=x0; tol=1; %迭代过程 while(tol>eps) x=(I-w*A)*x0+w*b; n=n+1;%n为最终求出解时的迭代步数tol=norm(x-x0); x0=x; if(n>=M) disp('Warning:迭代次数太多，可能不收敛！'); return; end end 3.grs里查森迭代法求线性方程组Ax=b的解 function[x,n]=grs(A,b,x0,W,eps,M) if(nargin==4) eps=1.0e-6;%eps表示迭代精度 M=10000;%M表示迭代步数的限制值 elseif(nargin==5) M=10000; end I=eye(size(A)); n=0; x=x0; tol=1;%前后两次迭代结果误差 %迭代过程 while(tol>eps) x=(I-W*A)*x0+W*b;%迭代公式 n=n+1;%n为最终求出解时的迭代步数tol=norm(x-x0); x0=x; if(n>=M) disp('Warning:迭代次数太多，可能不收敛！'); return; end end 4.jacobi雅可比迭代法求线性方程组Ax=b的解 function[x,n]=jacobi(A,b,x0,eps,varargin) if nargin==3 eps=1.0e-6; M=200; elseif nargin<3 error return elseif nargin==5 M=varargin{1}; end D=diag(diag(A));%求A的对角矩阵 L=-tril(A,-1);%求A的下三角阵 MatLab解线性方程组一文通 当齐次线性方程AX=0,rank(A)=r 一、Maple V 系统 Maple V是由Waterloo大学开发的数学系统软件，它不但具有精确的数值处理功能，而且具有无以伦比的符号计算功能。Maple V的符号计算能力还是MathCAD和MATLAB等软件的符号处理的核心。Maple提供了2000余种数学函数，涉及范围包括：普通数学、高等数学、线性代数、数论、离散数学、图形学。它还提供了一套内置的编程语言，用户可以开发自己的应用程序，而且Maple自身的2000多种函数，基本上是用此语言开发的。 Maple采用字符行输入方式，输入时需要按照规定的格式输入，虽然与一般常见的数学格式不同，但灵活方便，也很容易理解。输出则可以选择字符方式和图形方式，产生的图形结果可以很方便地剪贴到Windows应用程序内。 二、MATLAB 系统 MATLAB原是矩阵实验室（Matrix Laboratory）在70年代用来提供Linpack和Eispac k软件包的接口程序，采用C语言编写。从80年代出现3.0的DOS版本，逐渐成为科技计算、视图交互系统和程序语言。MATLAB可以运行在十几个操作平台上，比较常见的有基于W indows 9X/NT、OS/2、Macintosh、Sun、Unix、Linux等平台的系统。 MATLAB程序主要由主程序和各种工具包组成，其中主程序包含数百个内部核心函数，工具包则包括复杂系统仿真、信号处理工具包、系统识别工具包、优化工具包、神经网络工具包、控制系统工具包、μ分析和综合工具包、样条工具包、符号数学工具包、图像处理工具包、统计工具包等。而且5.x版本还包含一套几十个的PDF文件，从MATLAB的使用入门到其他专题应用均有详细的介绍。 MATLAB是数值计算的先锋，它以矩阵作为基本数据单位，在应用线性代数、数理统计、自动控制、数字信号处理、动态系统仿真方面已经成为首选工具，同时也是科研工作人员和大学生、研究生进行科学研究的得力工具。MATLAB在输入方面也很方便，可以使用内部的E ditor或者其他任何字符处理器，同时它还可以与Word6.0/7.0结合在一起，在Word的页面里直接调用MATLAB的大部分功能，使Word具有特殊的计算能力。 三、MathCAD 系统 MathCAD是美国Mathsoft公司推出的一个交互式的数学系统软件。从早期的DOS下的1. 0和Windows下的4.0版本，到今日的8.0版本，功能也从简单的数值计算，直至引用Map le强大的符号计算能力，使得它发生了一个质的飞跃。 MathCAD是集文本编辑、数学计算、程序编辑和仿真于一体的软件。MathCAD7.0 Profe ssional（专业版）运行在Win9X/NT下，它的主要特点是输入格式与人们习惯的数学书写格式很近似，采用WYSWYG（所见所得）界面，特别适合一般无须进行复杂编程或要求比较特殊的计算。MathCAD 7.0 Professional 还带有一个程序编辑器，对于一般比较短小，或者要求计算速度比较低时，采用它也是可以的。这个程序编辑器的优点是语法特别简单。 MathCAD可以看作是一个功能强大的计算器，没有很复杂的规则；同时它也可以和Wor d、Lotus、WPS2000等字处理软件很好地配合使用，可以把它当作一个出色的全屏幕数学公式编辑器。 四、Mathematica 系统 Mathematica是由美国物理学家Stephen Wolfram领导的Wolfram Research开发的数学系统软件。它拥有强大的数值计算和符号计算能力，在这一方面与Maple类似，但它的符 在这章中我们要学习线性方程组的直接法，特别是适合用数学软件在计算机上求解的方法. 3.1 方程组的逆矩阵解法及其MATLAB 程序 3.1.3 线性方程组有解的判定条件及其MATLAB 程序 判定线性方程组A n m ?b X =是否有解的MATLAB 程序 function [RA,RB,n]=jiepb(A,b) B=[A b];n=length(b); RA=rank(A); RB=rank(B);zhica=RB-RA; if zhica>0, disp('请注意：因为RA~=RB ，所以此方程组无解.') return end if RA==RB if RA==n disp('请注意：因为RA=RB=n ，所以此方程组有唯一解.') else disp('请注意：因为RA=RB 第四讲绘图功能 作为一个功能强大的工具软件，Matlab 具有很强的图形处理功能，提供了大量的二维、三维图形函数。由于系统采用面向对象的技术和丰富的矩阵运算，所以在图形处理方面即常方便又高效。 4.1 二维图形 一、plot函数 函数格式：plot(x,y)其中x和y为坐标向量函数功能：以向量x、y为轴，绘制曲线。【例1】在区间0≤X≤2 内，绘制正弦曲线Y=SIN（X），其程序为： x=0:pi/100:2*pi; y=sin(x); plot(x,y) 一、plot函数 【例2】同时绘制正、余弦两条曲线Y1=SIN（X）和Y2=COS（X），其程序为： x=0:pi/100:2*pi; y1=sin(x); y2=cos(x); plot(x,y1,x,y2) plot函数还可以为plot(x,y1,x,y2，x,y3，…)形式，其功能是以公共向量x为X轴，分别以y1，y2，y3，…为Y轴，在同一幅图内绘制出多条曲线。 一、plot函数 （一）线型与颜色 格式：plot(x,y1,’cs’,...) 其中c表示颜色，s表示线型。 【例3】用不同线型和颜色重新绘制例4.2图形，其程序为：x=0:pi/100:2*pi; y1=sin(x); y2=cos(x); plot(x,y1,'go',x,y2,'b-.') 其中参数'go'和'b-.'表示图形的颜色和线型。g表示绿色，o表示图形线型为圆圈；b表示蓝色，-.表示图形线型为点划线。 一、plot函数 （二）图形标记 在绘制图形的同时，可以对图形加上一些说明，如图形名称、图形某一部分的含义、坐标说明等，将这些操作称为添加图形标记。 title(‘加图形标题'); xlabel('加X轴标记'); ylabel('加Y轴标记'); text(X,Y,'添加文本'); 《MATLAB语言》课成论文 利用MATLAB求线性方程组 姓名：郭亚兰 学号：12010245331 专业：通信工程 班级：2010级通信工程一班 指导老师：汤全武 学院：物电学院 完成日期：2011年12月17日 利用MATLAB求解线性方程组 （郭亚兰 12010245331 2010 级通信一班） 【摘要】在高等数学及线性代数中涉及许多的数值问题，未知数的求解，微积分，不定积分，线性方程组的求解等对其手工求解都是比较复杂，而MATLAB语言正是处理线性方程组的求解的很好工具。线性代数是数学的一个分支，它的研究对象是向量，向量空间（或称线性空间），线性变换和有限维的线性方程组。因而，线性代数被广泛地应用于抽象代数和泛函分析中；由于科学研究中的非线性模型通常可以被近似为线性模型，使得线性代数被广泛地应用于自然科学和社会科学中。线性代数是数学的一个分支，它的研究对象是向量，向量空间（或称线性空间），线性变换和有限维的线性方程组。因而，线性代数被广泛地应用于抽象代数和泛函分析中；由于科学研究中的非线性模型通常可以被近似为线性模型，使得线性代数被广泛地应用于自然科学和社会科学中。线性代数是讨论矩阵理论、与矩阵结合的有限维向量空间及其线性变换理论的一门学科。 【关键字】线性代数MATLAB语言秩矩阵解 一、基本概念 1、N级行列式A：A等于所有取自不同性不同列的n个元素的积的代数和。 2、矩阵B：矩阵的概念是很直观的，可以说是一张表。 3、线性无关：一向量组（a1,a2,…,an）不线性相关，既没有不全为零的数 k1,k2,………kn使得：k1*a1+k2*a2+………+kn*an=0 4、秩：向量组的极在线性无关组所含向量的个数成为这个向量组的秩。 5、矩阵B的秩：行秩，指矩阵的行向量组的秩；列秩类似。记：R(B) 《Matlab与科学计算》作业 第一章MATLAB环境 1、MATLAB通用操作界面窗口包括哪些？命令窗口、历史命令窗口、当前目录窗口、工作空间窗口各有哪些功能？ 答：MATLAB通用操作界面窗口包括：命令窗口、历史命令窗口、当前目录浏览器窗口、工作空间窗口、变量编辑器窗口、M文件编辑/调试器窗口、程序性能剖析窗口、MATLAB帮助。 命令窗口是MATLAB命令操作的最主要窗口，可以把命令窗口当做高级的“草稿纸”。在命令窗口中可以输入各种MATLAB的命令、函数和表达式，并显示除图形外的所有运算结果。 历史命令窗口用来记录并显示已经运行过的命令、函数和表达式，并允许用户对它们进行选择、复制和重运行，用户可以方便地输入和修改命令，选择多行命令以产生M文件。 当前目录窗口用来设置当前目录，可以随时显示当前目录下的M、MKL等文件的信息，扬文件类型、文件名、最后个修改时间和文件的说明信息等，并可以复制、编辑和运行M文件及装载MAT数据文件。 工作空间窗口用来显示所有MATLAB工作空间中的变量名、数据结构、类型、大小和字节数。 2、熟悉课本中表格1.4、1.5、1.6、1.7、1.8的内容。 3、如何生成数据文件？如何把数据文件中的相关内容输入到工作空间中，用实例进行操作。 生成数据文件： 把数据文件中的相关内容输入到工作空间中： 结果： 4、在工作空间中可以通过哪些命令管理变量，写出每种语法的具体操作过程。答：（1）把工作空间中的数据存放到MAT数据文件。 语法：save filename 变量1 变量2 ……参数。 （2）从数据文件中取出变量存放到工作空间。 语法：load filename 变量1 变量2 ……。 盐城师范学院《MATLAB与科学计算》期末论文 2016－2017学年度第一学期 用MATLAB解决解析几何的图形问题 学生姓名吴梦成 学院数学与统计学院 专业信息与计算科学 班级数15(5)信计 学号 15213542 用MATLAB 解决解析几何的图形问题 摘 要 将 MATLAB 的图形和动画功能都用于解析几何教学，可使教学形象生动。以图形问题为例，详细给出了实例的程序编写和动画实现过程 。在解析几何教学中有一定的应用价值。 【关键词】： MATLAB ； 解析几何 ；图形 ； 动 画；编程 1 引 言 在解析几何的教学中，使用传统的教学方法。许多曲线及曲面的形成过程与变换过程只通过传统的教师讲授静态图示就很难形象生动地表示出来 。在解析几何教学中使用MATLAB 软件辅助教学，不仅可以很容易绘制出复杂的立体图形。把曲线、曲面的形成和变化过程准确地模拟出来 ，而且还能够对它们进行翻转 、旋转 ，甚 至还能够轻而易举地实现图形的动画效果!这对提高教学效率和培养学生的空间想象能力可起到事半功倍的效果。下面结合实例从几个方面说明MATLAB 在解析几何画图方面的应用。 2 利用 MATLAB 绘制三维曲线 在空间解析几何中，各种曲线和曲面方程的建立都离不开图形 ，而空间曲线和曲面图形既难画又费时。借助MATLAB 的绘图功能 ，可以快捷 、 准确地绘出图形，使教学变得形象 、生动 。有利于学生观察三维空间图形的形状 ， 掌握图形的性质 。 一 般地 ，MATLAB 可用plot3，ezplot3，comet3等函数来各种三维曲线 。 例如画螺旋曲线的图形，其参数方程设为 ：t at cos x =，t b sin t y -=，ct =z 。使用 plot3语句画螺旋曲线图形的方法如下( 设a =2 ,b=4,c=3)： );*3),sin(*.*4),cos(*.*2(3;*10:50/:0t t t t t plot pi pi t -= MATLAB 用两条简单的语句就可以画出螺旋 曲线(图1)，但上述方法是静态的 ，为了体 3.1 方程组的逆矩阵解法及其MATLAB 程序 3.1.3 线性方程组有解的判定条件及其MATLAB 程序 判定线性方程组A n m ?b X =是否有解的MATLAB 程序 function [RA,RB,n]=jiepb(A,b) B=[A b];n=length(b); RA=rank(A); RB=rank(B);zhica=RB-RA; if zhica>0, disp('请注意：因为RA~=RB ，所以此方程组无解.') return end if RA==RB if RA==n disp('请注意：因为RA=RB=n ，所以此方程组有唯一解.') else disp('请注意：因为RA=RB Matlab 与科学计算考试样题（客观题） 1 下面的MATLAB 语句中正确的有： a) 2a ＝pi 。 b) record_1=3+4i c) a=2.0, d) c=1+6j 2. 已知水的黏度随温度的变化公式如下，其中a=0.03368，b=0.000221，计算温度t 为20，30，40度时的粘度分别是： 2 1at bt μμ=++0μ为0℃水的黏度，值为31.78510-?；a 、b 为常数，分别为0.03368、0.000221。 3. 请补充语句以画出如图所示的图形： [x,y]=meshgrid(-2:0.1:2, -2:0.1:2)。 Z=x.*exp(-x.^2-y.^2)。 。 a) Plot3(x,y,Z) b) plot3(x,y,Z) c) mesh(x,y,Z) d) plot3(x,y,z) 2 a) 0.4900 1.2501 0.8560 b) 0.8560 1.2501 0.4900 c) -0.6341 3.8189 -3.7749 d) 3.8189 -3.7749 2.8533 解释说明： >> x=0.5:0.5:3.0。 >> y=[1.75,2.45,3.81,4.80,8.00,8.60]。 >> a=polyfit(x,y,2) a = 0.4900 1.2501 0.8560 >> x1=[0.5:0.25:3.0]。 >> y1=a(1)*x1.^2+a(2)*x1+a(3) >> plot(x,y,'*') >> hold on >> plot(x1,y1,'--r') 5. 求方程在 x=0.5附近的根. 21 x x += a) 0.6180 b) -1.1719e-25 c) －1 d) -1.6180 6. 用Newton-Cotes方法计算如下积分 1 5 x? （a）133.6625 (b)23.8600 (c) 87.9027 (d) -1.6180 7. y=ln(1+x)，求x=1时y" a) -0.25 b) 0.5 c) -0.6137 d) -1.6137 8.某公司用3台轧机来生产规格相同的铝合金薄板。取样测量薄板的 厚度，精确至‰厘M。得结果如下： 轧机1：0.236 0.238 0.248 0.245 0.243 轧机2：0.257 0.253 0.255 0.254 0.261 轧机3：0.258 0.264 0.259 0.267 0.262 计算方差分析结果，并判定各台轧机所生产的薄板的厚度有无显著的差异？ a) p＝1.3431e-005，没有显著差异。 实验1.1 用matlab 求解线性方程组 第一节 线性方程组的求解 一、齐次方程组的求解 rref （A ） %将矩阵A 化为阶梯形的最简式 null （A ） %求满足AX ＝0的解空间的一组基，即齐次线性方程组的基 础解系 【例１】 求下列齐次线性方程组的一个基础解系，并写出通解： 我们可以通过两种方法来解： 解法1： >> A=[1 -1 1 -1;1 -1 -1 1;1 -1 -2 2]; >> rref(A) 执行后可得结果： ans= 1 -1 0 0 0 0 -1 1 0 0 0 0 由最简行阶梯型矩阵，得化简后的方程 ??? ??=+--=+--=-+-0 22004321 43214321x x x x x x x x x x x x 取x2，x4为自由未知量，扩充方程组为 即 提取自由未知量系数形成的列向量为基础解系，记 所以齐次方程组的通解为 解法2： clear A=[1 -1 1 -1;1 -1 -1 1;1 -1 -2 2]; B=null(A, 'r') % help null 看看加个‘r’是什么作用, 若去掉r ，是什么结果？ 执行后可得结果： B= 1 0 1 0 0 1 0 1 ?? ?=-=-0 04321x x x x ?????? ?====4 4432221x x x x x x x x ??? ??? ??????+????????????=????? ???????1100001142 4321x x x x x x , 00111????? ? ??????=ε, 11002????? ???????=ε2 211εεk k x += Matlab线性方程组求解 1. Gauss消元法： function x=DelGauss(a,b) % Gauss消去法 [n,m]=size(a); nb=length(b); det=1; %存储行列式值 x=zeros(n,1); for k=1:n-1 for i=k+1:n if a(k,k)==0 return end m=a(i,k)/a(k,k); for j=k+1:n a(i,j)=a(i,j)-m*a(k,j); end b(i)=b(i)-m*b(k); end det=det*a(k,k); %计算行列式 end det=det*a(n,n); for k=n:-1:1 %回代求解 for j=k+1:n b(k)=b(k)-a(k,j)*x(j); end x(k)=b(k)/a(k,k); end Example: >> A=[1.0170 -0.0092 0.0095;-0.0092 0.9903 0.0136;0.0095 0.0136 0.9898]; >> b=[1 0 1]'; >> x=DelGauss(A,b) x = 0.9739 -0.0047 1.0010 2. 列主元Gauss消去法： function x=detGauss(a,b) % Gauss列主元消去法 [n,m]=size(a); nb=length(b); det=1; %存储行列式值 x=zeros(n,1); for k=1:n-1 amax=0; %选主元 for i=k:n if abs(a(i,k))>amax amax=abs(a(i,k));r=i; end end if amax<1e-10 return; end if r>k %交换两行 for j=k:n 求解线性方程组 solve，linsolve 例： A=[5 0 4 2;1 -1 2 1;4 1 2 0;1 1 1 1]; %矩阵的行之间用分号隔开，元素之间用逗号或空格 B=[3;1;1;0] X=zeros(4,1);%建立一个4元列向量 X=linsolve(A,B) diff（fun，var，n）：对表达式fun中的变量var求n阶导数。 例如：F=sym（'u(x,y)*v(x,y)'）; %sym（）用来定义一个符号表达式 diff(F); %matlab区分大小写 pretty(ans) %pretty（）：用习惯书写方式显示变量；ans是答案表达式 非线性方程求解 fsolve(fun,x0,options) 为待解方程或方程组的文件名；fun其中 x0位求解方程的初始向量或矩阵； option为设置命令参数 建立文件fun.m： function y=fun(x) y=[x(1)-0.5*sin(x(1))-0.3*cos(x(2)), ... x(2) - 0.5*cos(x(1))+0.3*sin(x(2))]; >>clear;x0=[0.1,0.1];fsolve(@fun,x0,optimset('fsolve')) 注： ...为续行符 m文件必须以function为文件头，调用符为@；文件名必须与定义的函数名相同；fsolve（）主要求解复杂非线性方程和方程组，求解过程是一个逼近过程。Matlab求解线性方程组 AX=B或XA=B 在MATLAB中，求解线性方程组时，主要采用前面章节介绍的除法运算符“/”和“\”。如： X=A\B表示求矩阵方程AX＝B的解； 的解。XA=B表示矩阵方程B/A＝X． 对方程组X＝A\B，要求A和B用相同的行数，X和B有相同的列数，它的行数等于矩阵A的列数，方程X＝B/A同理。 如果矩阵A不是方阵，其维数是m×n，则有： m＝n 恰定方程，求解精确解； m>n 超定方程，寻求最小二乘解； m 单元测验已完成成绩：100.0分 1 【单选题】 MATLAB一词来自（）的缩写。 ?A、 Mathematica Laboratory ?B、 Matrix Laboratory ?C、 MathWorks Lab ?D、 Matrices Lab 我的答案：B得分：50.0分 2 【单选题】 下列选项中能反应MATLAB特点的是（）。 ?A、 算法最优 ?B、 不需要写程序 ?C、 程序执行效率高 . . ?D、 编程效率高 我的答案：D得分：50.0分 单元测验已完成成绩：96.4分 1 【单选题】 当在命令行窗口执行命令时，如果不想立即在命令行窗口中输出结果，可以在命令后加上（）。 ?A、 冒号（:） ?B、 逗号（,） ?C、 分号（;） ?D、 百分号（%） 我的答案：C得分：7.1分 2 【单选题】 fix(264/100)+mod(264,10)*10的值是（）。 ?A、 86 . . ?B、 62 ?C、 423 ?D、 42 我的答案：D得分：7.1分 3 【单选题】 在命令行窗口输入下列命令后，x的值是（）。 >> clear >> x=i*j ?A、 不确定 ?B、 -1 ?C、 1 ?D、 i*j 我的答案：B得分：7.1分 . . 4 【单选题】 使用语句x=linspace(0,pi,6)生成的是（）个元素的向量。 ?A、 8 ?B、 7 ?C、 6 ?D、 5 我的答案：C得分：7.1分 5 【单选题】 ceil(-2.1)的结果为（）。 ?A、 -2 ?B、 -3 ?C、 1 . . ?D、 2 我的答案：A得分：7.1分 6 【单选题】 eval('sqrt(4)+2')的值是（）。 ?A、 sqrt(4)+2 ?B、 4 ?C、 2 ?D、 2+2 我的答案：B得分：7.1分 7 【单选题】 已知a为3×5矩阵，则执行完a(:,[2,4])=[]后（）。 ?A、 a变成行向量 ?B、 a变为3行2列 . . 线性方程组求解 1.直接法 Gauss消元法： function x=DelGauss(a,b) % Gauss消去法 [n,m]=size(a); nb=length(b); det=1;%存储行列式值 x=zeros(n,1); for k=1:n-1 for i=k+1:n if a(k,k)==0 return end m=a(i,k)/a(k,k); for j=k+1:n a(i,j)=a(i,j)-m*a(k,j); end b(i)=b(i)-m*b(k); end det=det*a(k,k); end det=det*a(n,n); for k=n:-1:1 %回代 for j=k+1:n b(k)=b(k)-a(k,j)*x(j); end x(k)=b(k)/a(k,k); end Example: >> A=[1.0170 -0.0092 0.0095;-0.0092 0.9903 0.0136;0.0095 0.0136 0.9898]; >> b=[1 0 1]'; >> x=DelGauss(A,b) x = 0.9739 -0.0047 1.0010 列主元Gauss消去法： function x=detGauss(a,b) % Gauss列主元消去法 [n,m]=size(a); nb=length(b); det=1;%存储行列式值 x=zeros(n,1); for k=1:n-1 amax=0;% 选主元 for i=k:n if abs(a(i,k))>amax amax=abs(a(i,k));r=i; end end if amax<1e-10 return; end if r>k %交换两行 for j=k:n z=a(k,j);a(k,j)=a(r,j);a(r,j)=z; end z=b(k);b(k)=b(r);b(r)=z;det=-det; end 用matlab解线性方程组 2008-04-12 17:00 一。高斯消去法 1.顺序高斯消去法 直接编写命令文件 a=[] d=[]' [n,n]=size(a); c=n+1 a(:,c)=d; for k=1:n-1 a(k+1:n, k:c)=a(k+1:n, k:c)-(a(k+1:n,k)/ a(k,k))*a(k, k:c); %消去 end x=[0 0 0 0]' %回带 x(n)=a(n,c)/a(n,n); for g=n-1:-1:1 x(g)=(a(g,c)-a(g,g+1:n)*x(g+1:n))/a(g,g) end 2.列主高斯消去法 *由于“[r,m]=max(abs(a(k:n,k)))”返回的行是“k:n,k”内的第几行，所以要通过修正来把m 改成真正的行的值。该程序只是演示程序，真正机器计算不需要算主元素所在列以下各行应为零的值。 直接编写命令文件 a=[] d=[] ' [n,n]=size(a); c=n+1 a(:,c)=d; %(增广) for k=1:n-1 [r,m]=max(abs(a(k:n,k))); %选主 m=m+k-1; %(修正操作行的值) if(a(m,k)~=0) if(m~=k) a([k m],:)=a([m k],:); %换行 end a(k+1:n, k:c)=a(k+1:n, k:c)-(a(k+1:n,k)/ a(k,k))*a(k, k:c); %消去end end x=[0 0 0 0]' %回带 x(n)=a(n,c)/a(n,n); for g=n-1:-1:1 x(g)=(a(g,c)-a(g,g+1:n)*x(g+1:n))/a(g,g) end 雅可比迭代 使用雅可比迭代法求解线性方程组的步骤 步骤1：输入系数矩阵A和方程组右端向量B； 步骤2：将矩阵A分解为下三角阵L对角阵D和上三角阵U 可分解为（D+L+U）X=B for o=1:n d(o,o)=a(o,o); u(o,o+1:n)=-a(o,o+1:n); end for p=2:n l(p,1:p-1)=-a(p,1:p-1); end； 步骤3：将上式化简为x=B0x+f,其中B0=-D-1(L+U),f=D-1B for i=1:n b0(i,i+1:n)=u(i,i+1:n)/a(i,i); f(i,:)=b(i,:)/a(i,i); end for p=2:n b0(p,1:p-1)=l(p,1:p-1)/a(p,p);； 步骤4：采用迭代公式在允许误差范围e=1e-7内求得解向量x x0=x; x=b0*x+f 雅可比迭代法matlab程序： function [x,k]=jacobi(a,b) n=length(a); e=1e-7; m=100; x0=zeros(n,1); x=x0; k=0; d=zeros(n); l=zeros(n); u=zeros(n); b0=zeros(n); f=zeros(n,1); x0=x+2*e; for o=1:n d(o,o)=a(o,o); u(o,o+1:n)=-a(o,o+1:n); end for p=2:n l(p,1:p-1)=-a(p,1:p-1); end for i=1:n b0(i,i+1:n)=u(i,i+1:n)/a(i,i); f(i,:)=b(i,:)/a(i,i); end for p=2:n b0(p,1:p-1)=l(p,1:p-1)/a(p,p); end while max(abs(x0-x))>e&k Matlab与科学计算考试样题（客观题） 1 下面的MATLAB语句中正确的有： a) 2a＝pi; b) record_1=3+4i c) a=2.0, d) c=1+6j 2. 已知水的黏度随温度的变化公式如下，其中a=0.03368，b=0.000221，计算温度t为20，30，40度时的粘度分别是： 为0℃水的黏度，值为；a、b为常数，分别为0.03368、0.000221。 （a）0.0018 0.0010 0.0007 (b) 0.0010 0.0007 0.0005 (0.0010 0.0008 0.0007) (c) 1.7850e-003 1.0131e-003 6.6092e-004 (d) 1.0131e-003 6.6092e-004 4.6772e-004 (1.0131e-003 8.0795e-004 6.6092e-004) a=0.03368;b=0.000221;u0=1.785e-3; t=[20 30 40];u=u0./(1+a*t+b*t.^2) >>format short %format short e >>u 3.请补充语句以画出如图所示的图形： [x,y]=meshgrid(-2:0.1:2, -2:0.1:2); Z=x.*exp(-x.^2-y.^2); ; a) Plot3(x,y,Z) b) plot3(x,y,Z) c) mesh(x,y,Z) d) plot3(x,y,z) 2 x 0.5 1.0 1.5 2.0 2.5 3.0 y 1.75 2.45 3.81 4.80 8.00 8.60 a) 0.4900 1.2501 0.8560 b) 0.8560 1.2501 0.4900 c) -0.6341 3.8189 -3.7749 d) 3.8189 -3.7749 2.8533 解释说明： >> x=0.5:0.5:3.0; Matlab求解线性方程组、非线性方程组 姓名：罗宝晶学号：15 专业：材料学院高分子系 第一部分数值计算 Matlab求解线性方程组AX=B或XA=B 在MATLAB中，求解线性方程组时，主要采用除法运算符“/”和“\”。如：X=A\B表示求矩阵方程AX＝B的解； X＝B/A表示矩阵方程XA=B的解。 对方程组X＝A\B，要求A和B用相同的行数，X和B有相同的列数，它的行数等于矩阵A的列数，方程X＝B/A同理。 如果矩阵A不是方阵，其维数是m×n，则有： m＝n 恰定方程，求解精确解； m>n 超定方程，寻求最小二乘解； m MATLAB与科学计算复习题题库 Chapter1 1.MA TLAB的主界面是一个高度集成的工作环境，有四个不同职责分工的窗口，分别为、、、和窗口. 2. 欲将MA TLAB的数值显示格式设置为近似有理数格式，可在命令提示符后输入． 3. 欲将MA TLAB的数值显示格式设置为长格式，可在命令提示符后输入． 4. 欲将MA TLAB的数值显示格式设置为短格式，可在命令提示符后输入． 5. 删除MA TLAB命令窗口的所有内容用命令． 6．欲退出MA TLAB环境，可在命令窗口中输入或命令． （提示：实验一下exit、quit函数） 1.请叙述MA TLAB中逗号、分号、冒号、方括号的作用． 1．下面工作空间管理的命令中，用于从磁盘中调入数据变量的命令为( ) A、clear B、load C、whos D、pack 1．下面的MATLAB特殊变量中，表示函数输入变量的数目的为( ) A、ans B、eps C、whos D、nargin 1．检查指定名字的变量或函数文件的存在性用________命令． A、clear B、load C、exist D、whos 2．欲退出MATLAB环境，可在命令窗口中输入命令． A、interp1 B、quit C、polyfit D、feval Chapter2 1.MA TLAB软件有一些常用的系统预定义的变量，如无穷大、圆周率π分别用表示、． 2.MA TLAB计算中， 应在命令窗口命令提示符后输入，的值，可在命令窗口命令提示符后输入． 3. MA TLAB语言中，合法的变量名须以开头，后可跟、、． 4.MA TLAB软件中，要输入矩阵 12 34 ?? ?? ?? ，应在命令窗口命令提示符后输入. 5．MA TLAB软件中，生成23 ?阶的全零、全一、单位矩阵时，应分别输入、、．6．MA TLAB软件中，用于求可逆矩阵A的行列式和逆矩阵的函数分别为、．7．语句A=linspace(2,18,9)，B=reshape(A,3,3)的执行结果为B=． 8. MA TLAB计算中，语句A=[1 2 3]；b=[-1 2 6]；c=dot(a,b)运行结果为c=． (提示：先用help dot在命令窗口中在线查询函数dot的用法) 9. 已知A=[1 1;2 4];B=[1 1;3 4];运行A.*B,A*B的结果分别为、． 10.A=zeros(2,4),A(:)=1:8;s=[2 3 7];则A(s)=． 11.A=[1 2 3;3 4 5;5 6 7]; 删除矩阵A的第三行的语句为． 12.A=[1 2 3;3 4 5;5 6 7]; 删除矩阵A的第三行的语句为． 13. A=[1 2 3;3 4 5;5 6 7]; C=[A,A+1];C(2,2:3)=． 14.MATLAB . 15.MATLAB软件中，要计算行列式12 34 ，应在命令窗口命令提示符后输入.MATLAB代码 解线性方程组的迭代法
MatLab求解线性方程组
MATLAB精通科学计算_偏微分方程求解
MATLAB解线性方程组的直接方法
科学计算与MATLAB语言(第四课)
利用MATLAB求线性方程组
《Matlab与科学计算》作业 2010010099
《MATLAB与科学计算》期末论文
线性方程组求解matlab实现
Matlab与科学计算样题(加主观题答案)
实验一用matlab求解线性方程组
Matlab线性方程组求解(Gauss消去法)
Matlab求解线性方程组非线性方程组
科学计算与MATLAB
线性方程组求解Matlab程序(精.选)
用matlab解线性方程组
雅可比解线性方程组matlab
Matlab与科学计算样题 (加主观题答案)
Matlab求解线性方程组、非线性方程组
MATLAB与科学计算复习题题库
相关主题
文本预览