matlab操作说明

MATLAB实验指导书

一、基础知识

1.1 常见数学函数

ceil(x)= -4 -2 0 2 5 7

fix(x) = -4 -2 0 1 4 6

floor(x) = -5 -3 -1 1 4 6

round(x) = -5 -2 0 1 5 7

1.2 系统的在线帮助

1 help 命令：

1.当不知系统有何帮助内容时，可直接输入help以寻求帮助:

>>help（回车）

2.当想了解某一主题的内容时，如输入：

>> help syntax（了解Matlab的语法规定）

3.当想了解某一具体的函数或命令的帮助信息时，如输入：

>> help sqrt （了解函数sqrt的相关信息）

2 lookfor命令

现需要完成某一具体操作，不知有何命令或函数可以完成，如输入：

>> lookfor line （查找与直线、线性问题有关的函数）

1.3 常量与变量

系统的变量命名规则：变量名区分字母大小写；变量名必须以字母打头，其后可以是任意字母，数字，或下划线的组合。此外，系统内部预先定义了几个有特殊意

1 数值型向量（矩阵）的输入

1．任何矩阵（向量），可以直接按行方式

．．．输入每个元素：同一行中的元素用逗号（，）或者用空格符来分隔；行与行之间用分号（；）分隔。所有元素处于一方括号（[ ]）内；

例1：

>> Time = [11 12 1 2 3 4 5 6 7 8 9 10]

>> X_Data = [2.32 3.43；4.37 5.98]

2

上面函数的具体用法，可以用帮助命令help得到。如：meshgrid(x,y)

输入x=[1 2 3 4]; y=[1 0 5]; [X,Y]=meshgrid(x, y)，则

X = Y =

1 2 3 4 1 1 1 1

1 2 3 4 0 0 0 0

1 2 3 4 5 5 5 5

目的是将原始数据x，y转化为矩阵数据X，Y。

2 符号向量（矩阵）的输入

1．用函数sym定义符号矩阵：

函数sym实际是在定义一个符号表达式，这时的符号矩阵中的元素可以是任何的符号或者是表达式，而且长度没有限制。只需将方括号置于单引号中。

例2：

>> sym_matrix = sym（'[a b c；Jack Help_Me NO_WAY]'）

sym_matrix =

[ a, b, c]

[Jack, Help_Me, NO_WAY]

2．用函数syms定义符号矩阵

先定义矩阵中的每一个元素为一个符号变量，而后像普通矩阵一样输入符号矩阵。

例3：

>> syms a b c ；

>> M1 = sym（'Classical'）；

>> M2 = sym（' Jazz'）；

>> M3 = sym（'Blues'）；

>> A = [a b c；M1，M2，M3；sym（[2 3 5]）]

A =

[ a, b, c]

[Classical, Jazz, Blues]

[ 2, 3, 5]

1.4 数组（矩阵）的点运算

运算符：+（加）、-（减）、./（右除）、.\（左除）、.^（乘方），

例4：

>> g = [1 2 3 4]；h = [4 3 2 1]；

>> s1 = g + h, s2 = g.*h, s3 = g.^h, s4 = g.^2, s5 = 2.^h

1.5 矩阵的运算

运算符：+（加）、-（减）、*（乘）、/（右除）、\（左除）、^（乘方）、’（转置）等；

常用函数：det（行列式）、inv（逆矩阵）、rank（秩）、eig（特征值、特征向量）、rref （化矩阵为行最简形）

例5：

>> A=[2 0 –1;1 3 2]; B=[1 7 –1;4 2 3;2 0 1];

>> M = A*B % 矩阵A与B按矩阵运算相乘

>> det_B = det(B) % 矩阵A的行列式

>> rank_A = rank（A）% 矩阵A的秩

>> inv_B = inv（B）% 矩阵B的逆矩阵

>> [V,D] = eig(B)% 矩阵B的特征值矩阵V与特征向量构成的矩阵D

>> X = A/B % A/B = A*B-1，即XB=A，求X

>> Y = B\A % B\A = B-1*A，即BY=A，求Y

上机练习（一）：

1．练习数据和符号的输入方式，将前面的命令在命令窗口中执行通过；

2．输入A=[7 1 5；2 5 6；3 1 5]，B=[1 1 1; 2 2 2; 3 3 3]，在命令窗口中执行下列表达式，掌握其含义：

A(2, 3) A(:,2) A(3,:) A(:,1:2:3) A(:,3).*B(:,2)

A(:,3)*B(2,:) A*B A.*B A^2 A.^2 B/A B./A 3．输入C=1:2:20，则C（i）表示什么？其中i=1,2,3, (10)

4．查找已创建变量的信息，删除无用的变量；

5．欲通过系统做一平面图，请查找相关的命令与函数，获取函数的帮助信息。

二、编程

2.1 无条件循环

当需要无条件重复执行某些命令时，可以使用for循环：

for 循环变量t=表达式1 : 达式2 : 表达式3

语句体

end

说明：表达式1为循环初值，表达式2为步长，表达式3为循环终值；当表达式2省略时则默认步长为1；for语句允许嵌套。

例6：如：矩阵输入程序

生成3×4阶的Hiltber矩阵。m=input(‘矩阵行数：m=’);

for i=1 : 3n= input(‘矩阵列数：n=’);

for j=1 : 4 for i=1:m

H（i，j）=1/（i+j-1）；for j=1:n

end disp([‘输入第’,num2str(i),’行，第’, num2str(j),’列元素’]) end A(i, j) = input (‘’)

end end

2.2 条件循环

1) if-else-then 语句

if-else-then 语句的常使用三种形式为：

(1) if 逻辑表达式 (3) if 逻辑表达式1 语句体 语句体1

end elseif 逻辑表达式2 语句体2 (2) if 逻辑表达式1 elseif 逻辑表达式3 语句体1 … else else

语句体2 语句体n

end end 2) while 循环语句

while 循环的一般使用形式为：

while 表达式 语句体 end 例7：

用二分法计算多项式方程523

--x x = 0在[0，3]内的一个根。 解：

a = 0；fa = -inf ；

b = 3；fb = inf ； while b-a > eps*b x =（a+b ）/2； fx = x^3-2*x-5； if sign(fx)== sign(fa) a =x ；fa = fx ； else

b = x ；fb = fx ； end end x

运行结果为：x = 2.0945515148154233

2.3 分支结构

若需要对不同的情形执行不同的操作，可用switch 分支语句： switch 表达式（标量或字符串）

case 值1

语句体1

case 值2

语句体2

…

otherwise

语句体n

end

说明：当表达式不是“case”所列值时，执行otherwise语句体。

2.4 建立M文件

将多个可执行

．．．的系统命令，用文本编辑器编辑后并存放在后缀为.m 的文件中，若在MATLAB命令窗口中输入该m-文件的文件名（不跟后缀.m!），即可依次执行该文件中的多个命令。这个后缀为.m的文件，也称为Matlab的脚本文件(Script File)。

注意：文件存放路径必须在Matlab能搜索的范围内。

2.5 建立函数文件

对于一些特殊用户函数，系统提供了一个用于创建用户函数的命令function，以备用户随时调用。

1．格式：

function [输出变量列表]=fun_name(输入变量列表)

用户自定义的函数体

2．函数文件名为：fun_name，注意：保存时文件名

．．．最好相同；

．．．与函数名

3．存储路径：最好在系统的搜索路径上。

4. 调用方法：输出参量=fun_name (输入变量)

例8：

计算s = n!，在文本编辑器中输入：

function s=pp(n);

s=1;

for i=1:n

s=s*i;

end

s;

在MA TLAB命令窗口中输入：s=pp(5)

结果为s = 120

上机练习（二）：

1．编写程序，计算1+3+5+7+…+(2n+1)的值（用input 语句输入n 值）。

2．编写分段函数 ??

?

??≤≤-<≤=其它021210)(x x x x

x f 的函数文件，存放于文件ff.m 中，计算出

)3(-f ，)2(f ，)(∞f 的值。

三、矩阵及其运算

3.1 矩阵的创建

1.加、减运算

运算符：“＋”和“－”分别为加、减运算符。

运算规则：对应元素相加、减，即按线性代数中矩阵的“十”、“一”运算进行。 例3-1 在Matlab 编辑器中建立m 文件：LX0701.m

A=[1, 1, 1; 1, 2, 3; 1, 3, 6] B=[8, 1, 6; 3, 5, 7; 4, 9, 2] A ＋B=A+B A-Ｂ=A-B

在Matlab 命令窗口建入LX0701，则 结果显示：A+B=

9 2 7 4 7 10 5 12 8 A －B=

-7 0 -5 -2 -3 -4 -3 -6 4 2.乘法

运算符：*

运算规则：按线性代数中矩阵乘法运算进行，即放在前面的矩阵的各行元素，分别与放在后面的矩阵的各列元素对应相乘并相加。

（1）两个矩阵相乘

例3-2 在Mtalab 编辑器中建立M 文件：LX0702.m

X= [2 3 4 5 1 2 2 1]；

Y=[0 1 1

1 1 0

0 0 1

1 0 0]；

Z=X*Y

存盘

在命令行中建入LX0702，回车后显示：

Z=

8 5 6

3 3 3

（2）矩阵的数乘：数乘矩阵

上例中：a=2*X

则显示：a =

4 6 8 10

2 4 4 2

（3）向量的点乘（内积）：维数相同的两个向量的点乘。

命令：dot向量点乘函数

例：X=[-1 0 2]；

Y=[-2 -1 1];

Z=dot(X, Y)

则显示：Z =

4

还可用另一种算法：

sum(X.*Y)

ans=

4

（4）向量叉乘

在数学上，两向量的叉乘是一个过两相交向量的交点且垂直于两向量所在平面的向量。在Matlab中，用函数cross实现。

命令cross 向量叉乘函数

例3-3 计算垂直于向量(1, 2, 3)和(4, 5, 6)的向量。

在Mtalab编辑器中建立M文件：LX0703.m

a=[1 2 3];

b=[4 5 6];

c=cross(a,b)

结果显示：

c=

-3 6 -3

可得垂直于向量(1, 2, 3)和(4, 5, 6)的向量为±(-3, 6, -3)

（5）混合积

混合积由以上两函数实现：

例3-4 计算向量a=(1, 2, 3)、b=(4, 5, 6)和c=(-3, 6, -3) 的混合积)c b (a ?? 在Matlab 编辑器中建立M 文件：LX0704.m a=[1 2 3]; b=[4 5 6]; c=[-3 6 -3]; x=dot(a, cross(b, c)) 结果显示：x = 54

注意：先叉乘后点乘，顺序不可颠倒。 3.矩阵的除法

Matlab 提供了两种除法运算：左除（\）和右除（/）。一般情况下，x=a\b 是方程a*x =b 的解，而x=b/a 是方程x*a=b 的解

例：a=[1 2 3; 4 2 6; 7 4 9]

b=[4; 1; 2]; x=a\b 则显示：x=

-1.5000

2.0000

0.5000

如果a 为非奇异矩阵，则a\b 和b/a 可通过a 的逆矩阵与b 阵得到： a\b = inv(a)*b b/a = b*inv(a)

4.矩阵乘方

运算符：^ 运算规则：

（1）当A 为方阵，p 为大于0的整数时，A^P 表示A 的P 次方，即A 自乘P 次；p 为小于0的整数时，A^P 表示A -1的P 次方。

（2）当A 为方阵，p 为非整数时，则1p nn p 11V d d V P ^A -????

?

???

??= 其中V 为A 的特征向量，?

???

???

?nn 11d d

为特征值矩阵 5.矩阵的转置

运算符：′

运算规则：与线性代数中矩阵的转置相同。

6.矩阵的逆矩阵

例3-5 求???

?

? ??=343122321A 的逆矩阵

方法一：在Matlab 编辑器中建立M 文件：LX07051.m A=[1 2 3; 2 2 1; 3 4 3]; inv(A)或A^(-1) 则结果显示为 ans =

1.0000 3.0000 -

2.0000 -1.5000 -

3.0000 2.5000 1.0000 1.0000 -1.0000

方法二：由增广矩阵???

?

? ??=100343010122001321B 进行初等行变换

在Matlab 编辑器中建立M 文件：LX07052.m

B=[1, 2, 3, 1, 0, 0; 2, 2, 1, 0, 1, 0; 3, 4, 3, 0, 0, 1]； C=rref(B) %化行最简形 X=C(:, 4:6)

在Matlab 命令窗口建入LX07052，则显示结果如下： C =

1.0000 0 0 1.0000 3.0000 -

2.0000 0 1.0000 0 -1.5000 -

3.0000 2.5000 0 0 1.0000 1.0000 1.0000 -1.0000 X =

1.0000 3.0000 -

2.0000 -1.5000 -

3.0000 2.5000 1.0000 1.0000 -1.0000 这就是A 的逆矩阵。 7.方阵的行列式

命令： det 计算行列式的值 例3-6 计算上例中A 的行列式的值

在Matlab 编辑器中建立M 文件：LX0706.m A=[1 2 3; 2 2 1; 3 4 3]; D=det(A)

则结果显示为 D = 2

3.2 符号矩阵的运算

1.符号矩阵的四则运算

Matlab 5.x 抛弃了在4.2版中为符号矩阵设计的复杂函数形式，把符号矩阵的四则运算简化为与数值矩阵完全相同的运算方式，其运算符为：加（＋），减（－）、乘（×）、除（/、\）等或：符号矩阵的和（symadd ），差（symsub ），乘 (symmul)。

例3-7 );])3x /(1),2x /(1);1x /(1,x /1[(sym A '+++'= )]0,2x ;1,x [(sym B '+'=；

C=B-A D=a\b 则显示：

C=

x-1/x 1-1/(x+1)

x+2-1/(x+2) -1/(x+3) D=

-6*x-2*x^3-7*x^2 1/2*x^3+x+3/2*x^2

6+2*x^3+10*x^2+14*x -2*x^2-3/2*x-1/2*x^3 2.其他基本运算

符号矩阵的其他一些基本运算包括转置（'）、行列式（det ）、逆（inv ）、秩（rank ）、幂（^）和指数（exp 和expm ）等都与数值矩阵相同

3.符号矩阵的简化

符号工具箱中提供了符号矩阵因式分解、展开、合并、简化及通分等符号操作函数。 （1）因式分解

命令：factor 符号表达式因式分解函数 格式：factor(s)

说明：s 为符号矩阵或符号表达式。常用于多项式的因式分解 例3-8 将x 9-1分解因式 在Matlab 命令窗口建入

syms x factor(x^9-1) 则显示：ans =

(x-1)*(x^2+x+1)*(x 6+x^3+1)

例3-9 问入取何值时，齐次方程组

???

??=λ-++=+λ-+=+-λ-0

)1(0)3(2042)1(321321321x x x x x x x x x 有非0解

解：在Matlab编辑器中建立M文件：LX0709.m

syms k

A=[1-k -2 4;2 3-k 1;1 1 1-k];

D=det(A)

factor(D)

其结果显示如下：

D =

-6*k+5*k^2-k^3

ans =

-k*(k-2)*(-3+k)

从而得到：当k=0、k=2或k=3时，原方程组有非0解。

（2）符号矩阵的展开

命令expand符号表达式展开函数

格式：expand(s)

说明：s为符号矩阵或表达式。常用在多项式的因式分解中，也常用于三角函数，指数函数和对数函数的展开中

例3-10 将(x+1)3、sin(x+y)展开

在Matlab编辑器中建立M文件：LX0710.m

syms x y

p=expand((x+1)^3)

q=expand(sin(x+y))

则结果显示为

p =

x^3+3*x^2+3*x+1

q =

sin(x)*cos(y)+cos(x)*sin(y)

（3）同类式合并

命令：Collect合并系数函数

格式：Collect(s,v) 将s中的变量v的同幂项系数合并。

Collect(s) s —矩阵或表达式，此命令对由命令findsym函数返回的默认变量进行同类项合并。

（4）符号简化

命令：simple或simplify寻找符号矩阵或符号表达式的最简型

格式：Simple(s) s —矩阵或表达式

说明：Simple(s)将表达式s的长度化到最短。若还想让表达式更加精美，可使用函数Pretty。

格式：Pretty(s) 使表达式s更加精美

例3-11 计算行列式

4

44422221111d c b a d c b a d c b a 的值。

在Matlab 编辑器中建立M 文件：LX0711.m syms a b c d

A=[1 1 1 1;a b c d;a^2 b^2 c^2 d^2;a^4 b^4 c^4 d^4]; d1=det(A)

d2=simple(d1) %化简表达式d1

pretty(d2) %让表达式d2符合人们的书写习惯 则显示结果如下： d1 =

b*c^2*d^4-b*d^2*c^4-b^2*c*d^4+b^2*d*c^4+b^4*c*d^2-b^4*d*c^2-a*c^2*d^4+a*d^2*c ^4+a*b^2*d^4-a*b^2*c^4-a*b^4*d^2+a*b^4*c^2+a^2*c*d^4-a^2*d*c^4-a^2*b*d^4+a^2*b*c^4+a^2*b^4*d-a^2*b^4*c-a^4*c*d^2+a^4*d*c^2+a^4*b*d^2-a^4*b*c^2-a^4*b^2*d+a^4*b^2*c

d2 =

(-d+c)*(b-d)*(b-c)*(-d+a)*(a-c)*(a-b)*(a+c+d+b) (-d+c)(b-d)(b-c)(-d+a)(a-c)(a-b)(a+c+d+b)

例3-12 设????? ??=343122321A ，??? ??=3512B ，?

??

?

? ??=130231C

求矩阵X ，使满足：AXB = C

在Matlab 编辑器中建立M 文件：LX0712.m A=[1 2 3;2 2 1;3 4 3]; B=[2,1;5 3]; C=[1 3;2 0;3 1]; X=A\C/B

则结果显示如下： X =

-2.0000 1.0000 10.0000 -4.0000 -10.0000 4.0000

例3-13 计算5

)t cos()t sin()t sin()t cos(??

?

??-

在Matlab 编辑器中建立M 文件：LX0713.m

syms t

A =[cos(t) -sin(t); sin(t), cos(t)];

B=sympow(A, 5) %计算A 的5次幂

C=simple(B) %化简 pretty(C)

则显示结果如下 B =

[ cos(t)*(cos(t)*(cos(t)*(cos(t)^2-sin(t)^2)-2*sin(t)^2*cos(t))-sin(t)*(sin(t)*(cos(t)^2-sin(t)^2)+2*cos(t)^2*sin(t)))-sin(t)*(sin(t)*(cos(t)*(cos(t)^2-sin(t)^2)-2*sin(t)^2*cos(t))+cos(t)*(sin(t)*(cos(t)^2-sin(t)^2)+2*cos(t)^2*sin(t))),

cos(t)*(cos(t)*(-2*cos(t)^2*sin(t)-sin(t)*(cos(t)^2-sin(t)^2))-sin(t)*(cos(t)*(cos(t)^2-sin(t)^2)-2*s in(t)^2*cos(t)))-sin(t)*(sin(t)*(-2*cos(t)^2*sin(t)-sin(t)*(cos(t)^2-sin(t)^2))+cos(t)*(cos(t)*(cos(t )^2-sin(t)^2)-2*sin(t)^2*cos(t)))]

[ sin(t)*(cos(t)*(cos(t)*(cos(t)^2-sin(t)^2)-2*sin(t)^2*cos(t))-sin(t)*(sin(t)*(cos(t)^2-sin(t)^2)+2*cos(t)^2*sin(t)))+cos(t)*(sin(t)*(cos(t)*(cos(t)^2-sin(t)^2)-2*sin(t)^2*cos(t))+cos(t)*(sin(t )*(cos(t)^2-sin(t)^2)+2*cos(t)^2*sin(t))),

sin(t)*(cos(t)*(-2*cos(t)^2*sin(t)-sin(t)*(cos(t)^2-sin(t)^2))-sin(t)*(cos(t)*(cos(t)^2-sin(t)^2)-2*sin(t)^2*cos(t)))+cos(t)*(sin(t)*(-2*cos(t)^2*sin(t)-sin(t)*(cos(t)^2-sin(t)^2))+cos(t)*(cos(t)*(cos(t)^2-sin(t)^2)-2*sin(t)^2*cos(t)))]

C =

[ cos(5*t), -sin(5*t)] [ sin(5*t), cos(5*t)]

[cos(5 t) -sin(5 t)]

[ ] [sin(5 t) cos(5 t) ]

四、秩与线性相关性

4.1 矩阵和向量组的秩以及向量组的线性相关性。

矩阵A 的秩是矩阵A 中最高阶非零子式的阶数；向量组的秩通常由该向量组构成的矩阵来计算。

命令：rank

格式：rank(A) A 为矩阵式向量组构成的矩阵

例4-1 求向量组)3221(1-=α，)3142(2--=α，)3021(3-=α，

)3260(4=α，)4362(5-=α的秩，并判断其线性相关性。

在Matlab 编辑器中建立M 文件：LX0714.m A=[1 -2 2 3;-2 4 -1 3;-1 2 0 3;0 6 2 3;2 -6 3 4]; B=rank(A)

运行后结果如下： B =

3

由于秩为3 < 向量个数，因此向量组线性相关。

4.2 向量组的最大无关组

矩阵的初等行变换有三条：

1.交换两行 j i r r ? （第i 、第j 两行交换）

2.第i 行的K 倍 i r k

3.第i 行的K 倍加到第j 行上去 i j r k r +

通过这三条变换可以将矩阵化成行最简形，从而找出列向量组的一个最大无关组，Matlab 将矩阵化成行最简形的命令是：

命令：rref

格式：rref(A) A 为矩阵 例4-2 求向量组a1=(1,-2,2,3)，a2=(-2,4,-1,3)，a3=(-1,2,0,3)，a4=(0,6,2,3)，a5=(2,-6,3,4)的一个最大无关组。

在Matlab 编辑器中建立M 文件：LX0715.m a1=[1 -2 2 3]'; a2=[-2 4 -1 3]'; a3=[-1 2 0 3]'; a4=[0 6 2 3]'; a5=[2 -6 3 4]';

A=[a1 a2 a3 a4 a5]

format rat %以有理格式输出 B=rref(A) %求A 的行最简形 运行后的结果为 A =

1 -

2 -1 0 2 -2 4 2 6 -6 2 -1 0 2

3 3 3 3 3

4 B =

1 0 1/3 0 16/9 0 1 2/3 0 -1/9 0 0 0 1 -1/3 0 0 0 0 0 从B 中可以得到：向量a1 a

2 a4为其中一个最大无关组

五、线性方程的组的求解

我们将线性方程的求解分为两类：一类是方程组求唯一解或求特解，另一类是方程组求无穷解即通解。可以通过系数矩阵的秩来判断：

若系数矩阵的秩r=n （n 为方程组中未知变量的个数），则有唯一解 若系数矩阵的秩r

线性方程组的无穷解 = 对应齐次方程组的通解+非齐次方程组的一个特解；其特解的求法属于解的第一类问题，通解部分属第二类问题。

5.1 求线性方程组的唯一解或特解（第一类问题）

这类问题的求法分为两类：一类主要用于解低阶稠密矩阵 —— 直接法；另一类是解大型稀疏矩阵 —— 迭代法。

5.1.1 利用矩阵除法求线性方程组的特解（或一个解） 方程：AX=b 解法：X=A\b

例5-1 求方程组?

??

?

???=+=++=++=++=+1

x 5x 0x 6x 5x 0x 6x 5x 0

x 6x 5x 1

x 6x 5545434323212

1的解 解：在Matlab 编辑器中建立M 文件：LX0716.m A=[5 6 0 0 0 1 5 6 0 0 0 1 5 6 0 0 0 1 5 6 0 0 0 1 5]; B=[1 0 0 0 1]';

R_A=rank(A) %求秩 X=A\B %求解 运行后结果如下 R_A = 5 X =

2.2662 -1.7218 1.0571 -0.5940 0.3188

这就是方程组的解。 例5-2 求方程组

???

??=--+=+--=--+0

x 8x 9x 5x 4x 4x 3x x 31x x 3x x 432143214321 的一个特解

解：在Matlab 编辑器中建立M 文件：LX0717.m

A=[1 1 -3 -1;3 -1 -3 4;1 5 -9 -8];

B=[1 4 0]'; X=A\B X = 0 0

-0.5333 0.6000

5.1.2 利用矩阵的LU 、QR 和cholesky 分解求方程组的解

1.LU 分解：

LU 分解又称Gauss 消去分解，可把任意方阵分解为下三角矩阵的基本变换形式（行交换）和上三角矩阵的乘积。即A=LU ，L 为下三角阵，U 为上三角阵。

则：A*X=b 变成L*U*X=b

∴X=U\(L\b) 这样可以大提高运算速度。 命令 [L ，U]=lu (A) 例5-3 求方程组

???

??=+=+-=-+8x 3x 1110x 2x x 32x x 2x 42

1321321 的一个特解

解

]8,10,2[0311213124'=????

? ??--=b A

在Matlab 编辑器中建立M 文件：LX0718.m A=[4 2 -1;3 -1 2;11 3 0]; B=[2 10 8]'; D=det(A) [L,U]=lu(A) X=U\(L\B)

在Matlab 命令窗口建入LX0718，回车后显示结果如下： D =

0 L =

0.3636 -0.5000 1.0000 0.2727 1.0000 0 1.0000 0 0 U =

11.0000 3.0000 0 0 -1.8182 2.0000 0 0 0.0000

Warning: Matrix is close to singular or badly scaled.

Results may be inaccurate. RCOND = 2.018587e-017. > In D:\Matlab\pujun\lx0720.m at line 4 X =

1.0e+016 * -0.4053 1.4862 1.3511

说明：结果中的警告是由于系数行列式为零产生的。可以通过A*X 验证其正确性。 2.Cholesky 分解

若A 为对称正定矩阵，则Cholesky 分解可将矩阵A 分解成上三角矩阵和其转置的乘积，即：R R A *'= 其中R 为上三角阵。

方程 A*X=b 变成 b X *R R =*'

所以 )b \R (\R X '=

命令 R=chol(A) 3.QR 分解

对于任何长方矩阵A ，都可以进行QR 分解，其中Q 为正交矩阵，R 为上三角矩阵的初等变换形式，即：A=QR

方程 A*X=b 变形成 QRX=b

所以 X=R\(Q\b) 命令 [Q, R]=qr(A) 上例中 [Q, R]=qr(A)

X=R\(Q\B) 说明：这三种分解，在求解大型方程组时很有用。其优点是运算速度快、可以节省磁盘空间、节省内存。

5.2 求线性齐次方程组的通解

在Matlab 中，函数null 用来求解零空间，即满足A ·X=0的解空间，实际上是求出解

空间的一组基（基础解系）。

格式：z = null % z 的列向量为方程组的正交规范基，满足I Z Z =?'

)r ,A (n u l l z ''= % z 的列向量是方程AX=0的有理基 例5-4 求解方程组的通解：

???

??=---=--+=+++0

x 3x 4x x 0x 2x 2x x 20

x x 2x 2x 432143214321 解：在Matlab 编辑器中建立M 文件：LX0719.m A=[1 2 2 1;2 1 -2 -2;1 -1 -4 -3]; format rat %指定有理式格式输出 B=null(A,'r') %求解空间的有理基 运行后显示结果如下：

B =

2 5/

3 -2 -4/3 1 0 0 1 写出通解： syms k1 k2

X=k1*B(:,1)+k2*B(:,2) %写出方程组的通解 pretty(X) %让通解表达式更加精美 运行后结果如下： X =

[ 2*k1+5/3*k2] [ -2*k1-4/3*k2] [ k1] [ k2]

% 下面是其简化形式 [2k1 + 5/3k2 ] [ ] [-2k1 - 4/3k2]

[ ] [ k1 ] [ ] [ k2 ]

5.3 求非齐次线性方程组的通解

非齐次线性方程组需要先判断方程组是否有解，若有解，再去求通解。因此， 步骤为：第一步：判断AX=b 是否有解，若有解则进行第二步

matlab电力系统潮流计算

华中科技大学 信息工程学院课程设计报告书题目: 电力系统潮流计算 专业：电气工程及其自动化 班级： 学号： 学生姓名： 指导教师： 2015年 11 月 10 日

2015年11月12日

信息工程学院课程设计成绩评定表

摘要 电力系统稳态分析包括潮流计算和静态安全分析。本文主要运用的事潮流计算，潮流计算是电力网络设计与运行中最基本的运算，对电力网络的各种设计方案及各种运行方式进行潮流计算，可以得到各种电网各节点的电压，并求得网络的潮流及网络中的各元件的电力损耗，进而求得电能损耗。本位就是运用潮流计算具体分析，并有MATLAB仿真。 关键词：电力系统潮流计算 MATLAB仿真

Abstract Electric power system steady flow calculation and analysis of the static safety analysis. This paper, by means of the calculation, flow calculation is the trend of the power network design and operation of the most basic operations of electric power network, various design scheme and the operation ways to tide computation, can get all kinds of each node of the power grid voltage and seek the trend of the network and the network of the components of the power loss, and getting electric power. The standard is to use the power flow calculation and analysis, the specific have MATLAB simulation. Key words: Power system; Flow calculation; MATLAB simulation

潮流计算(matlab)实例计算

潮流例题：根据给定的参数或工程具体要求（如图），收集和查阅资料；学习相关软件（软件自选：本设计选择Matlab进行设计）。 2.在给定的电力网络上画出等值电路图。 3.运用计算机进行潮流计算。 4.编写设计说明书。 一、设计原理 1．牛顿-拉夫逊原理 牛顿迭代法是取x0 之后，在这个基础上，找到比x0 更接近的方程的跟，一步一步迭代，从而找到更接近方程根的近似跟。牛顿迭代法是求方程根的重要方法之一，其最大优点是在方程f(x) = 0 的单根附近具有平方收敛，而且该法还可以用来求方程的重根、复根。电力系统潮流计算，一般来说，各个母线所供负荷的功率是已知的，各个节点电压是未知的（平衡节点外）可以根据网络结构形成节点导纳矩阵，然后由节点导纳矩阵列写功率方程，由于功率方程里功率是已知的，电压的幅值和相角是未知的，这样潮流计算的问题就转化为求解非线性方程组的问题了。为了便于用迭代法解方程组，需要将上述功率方程改写成功率平衡方程，并对功率平衡方程求偏导，得出对应的雅可比矩阵，给未知节点赋电压初值，一般为额定电压，将初值带入功率平衡方程，得到功率不平衡量，这样由功率不平衡量、雅可比矩阵、节点电压不平衡量（未知的）构成了误差方程，解误差方程，得到节点电压不平衡量，节点电压加上节点电压不平衡量构成新的节点电压初值，将新的初值带入原来的功率平衡方程，并重新形成雅可比矩阵，然后计算新

的电压不平衡量，这样不断迭代，不断修正，一般迭代三到五次就能收敛。 牛顿—拉夫逊迭代法的一般步骤： （1）形成各节点导纳矩阵Y。 （2）设个节点电压的初始值U和相角初始值e 还有迭代次数初值为0。 （3）计算各个节点的功率不平衡量。 （4）根据收敛条件判断是否满足，若不满足则向下进行。 （5）计算雅可比矩阵中的各元素。 （6）修正方程式个节点电压 （7）利用新值自第（3）步开始进入下一次迭代，直至达到精度退出循环。 （8）计算平衡节点输出功率和各线路功率 2．网络节点的优化 １）静态地按最少出线支路数编号 这种方法由称为静态优化法。在编号以前。首先统计电力网络个节点的出线支路数，然后，按出线支路数有少到多的节点顺序编号。当由n 个节点的出线支路相同时，则可以按任意次序对这n 个节点进行编号。这种编号方法的根据是导纳矩阵中，出线支路数最少的节点所对应的行中非零元素也２）动态地按增加出线支路数最少编号在上述的方法中，各节点的出线支路数是按原始网络统计出来的，在编号过程中认为固定不变的，事实上，在节点消去过程中，每消去一个节点以后，与该节点相连的各节点的出线支路数将发生变化(增加，减少或保持不变)。因此，如果每消去一个节点后，立即修正尚未编号节点的出线支路数，然后选其中支路数最少的一个节点进行编号，就可以预期得到更好的效果，动态按最少出线支路数编号方法的特点就是按出线最少原则编号时考虑了消去过程中各节点出线支路数目的变动情况。 3．MATLAB编程应用 Matlab 是“Matrix Laboratory”的缩写，主要包括：一般数值分析，矩阵运算、数字信号处理、建模、系统控制、优化和图形显示等应用程序。由于使用Matlab 编程运算与人进行科学计算的思路和表达方式完全一致，所以不像学习高级语言那样难于掌握，而且编程效率和计算效率极高，还可在计算机上直接输出结果和精美的图形拷贝，所以它的确为一高效的科研助手。 二、设计内容 1.设计流程图

matlab中GUI设计

MATLAB的GUI 程序设计 Chapter 8: Design of MATLAB of GUI program GUI(Graphical User Interfaces):由各种图形对象组成的用户界面，在这种用户界面下，用户的命令和对程序的控制是通过“选择”各种图形对象来实现的。 基本图形对象分为控件对象和用户界面菜单对象，简称控件和菜单。 一. 控件对象及属性(Object and its attributes of controller)) 1. GUI控件对象类型(The mode of controller object) 控件对象是事件响应的图形界面对象。当某一事件发生时，应用程序会做出响应并执行某些预定的功能子程序（Callback）. 控件对象及其功能：（表7—1） 2. 控件对象的描述(Description of controller object) MATLAB中的控件大致可分为两种，一种为动作控件，鼠标点击这些控件时会产生相应的响应。一种为静态控件，是一种不产生响应的控件，如文本框等。

每种控件都有一些可以设置的参数，用于表现控件的外形、功能及效果，既属性。属性由两部分组成：属性名和属性值，它们必须是成对出现的。 （1）按钮(Push Buttons)：执行某种预定的功能或操作； （2）开关按钮(Toggle Button)：产生一个动作并指示一个二进制状态（开或关），当鼠点击它时按钮将下陷，并执行callback（回调函数）中指定的内容，再次点击，按钮复原，并再次执行callback 中的内容； （3）单选框(Radio Button)：单个的单选框用来在两种状态之间切换，多个单选框组成一个单选框组时，用户只能在一组状态中选择单一的状态，或称为单选项； （4）复选框(Check Boxes)：单个的复选框用来在两种状态之间切换，多个复选框组成一个复选框组时，可使用户在一组状态中作组合式的选择，或称为多选项； （5）文本编辑器(Editable Texts)：用来使用键盘输入字符串的值，可以对编辑框中的内容进行编辑、删除和替换等操作； （6）静态文本框(Static Texts)：仅仅用于显示单行的说明文字； （7）滚动条(Slider)：可输入指定范围的数量值；

matlab潮流计算

附录1 使用牛顿拉夫逊法进行潮流计算的Matlab程序代码 % 牛拉法计算潮流程序 %----------------------------------------------------------------------- % B1矩阵:1、支路首端号;2、末端号;3、支路阻抗;4、支路对地电纳 % 5、支路的变比;6、支路首端处于K侧为1,1侧为0 % B2矩阵:1、该节点发电机功率;2、该节点负荷功率;3、节点电压初始值 % 4、PV节点电压V的给定值;5、节点所接的无功补偿设备的容量 % 6、节点分类标号:1为平衡节点(应为1号节点);2为PQ节点;3为PV 节点; %------------------------------------------------------------------------ clear all; format long; n=input('请输入节点数:nodes='); nl=input('请输入支路数:lines='); isb=input('请输入平衡母线节点号:balance='); pr=input('请输入误差精度:precision='); B1=input('请输入由各支路参数形成的矩阵:B1='); B2=input('请输入各节点参数形成的矩阵:B2='); Y=zeros(n);e=zeros(1,n);f=zeros(1,n);V=zeros(1,n);sida=zeros(1,n);S1=zeros(nl); %------------------------------------------------------------------ for i=1:nl %支路数 if B1(i,6)==0 %左节点处于1侧 p=B1(i,1);q=B1(i,2); else %左节点处于K侧 p=B1(i,2);q=B1(i,1); end Y(p,q)=Y(p,q)-1、/(B1(i,3)*B1(i,5)); %非对角元 Y(q,p)=Y(p,q); %非对角元 Y(q,q)=Y(q,q)+1、/(B1(i,3)*B1(i,5)^2)+B1(i,4); %对角元K侧 Y(p,p)=Y(p,p)+1、/B1(i,3)+B1(i,4); %对角元1侧 end %求导纳矩阵 disp('导纳矩阵Y='); disp(Y) %------------------------------------------------------------------- G=real(Y);B=imag(Y); %分解出导纳阵的实部与虚部 for i=1:n %给定各节点初始电压的实部与虚部 e(i)=real(B2(i,3)); f(i)=imag(B2(i,3));

matlab实现数值计算功能源程序(个人整理)

matlab数值计算功能 1,基础运算 （1）多项式的创建与表达 将多项式(x-6)(x-3)(x-8)表示为系数形式 a=[6 3 8] % 写成根矢量 pa=poly(a)% 求出系数矢量 ppa=poly2sym(pa,'x') % 表示成符号形式 ezplot(ppa,[-50,50]) 求3介方阵A的特征多项式 a=[6 2 4;7 5 6;1 3 6 ]; pa=poly(a)% 写出系数矢量 ppa=poly2sym(pa) %表示成符号形式 ezplot(ppa,[-50,50]) % 绘图 求x^3-6x^2-72x-27的根。 a=[1,-6,-72,-85]; % 写出多项式系数矢量 r=roots(a) % 求多项式的根 （2）多项式的乘除运算 c=conv(a,b) %乘法 [q,r]=deconv(c,a)% 除法 求a(s)=s^2+2s+3乘以b(s)=4s^2+5s+6的乘积 a=[1 2 3] b=[4 5 6] % 写出系数矢量 c=conv(a,b) c=poly2sym(c,'s') % 写成符号形式的多项式 展开(s^2+2s+2)(s+4)(s+1)并验证结果被(s+4),(s+3)除后的结果。c=conv([1,2,2],conv([1,4],[1,1])) cs=poly2sym(c,'s') c=[1 7 16 18 8] [q1,r1]=deconv(c,[1,4]) [q2,r2]=deconv(c,[1,3]) cc=conv(q2,[1,3]) test=((c-r2)==cc)

其他常用的多项式运算命令 pa=polyval(p,s) % 按数组规则计算给定s时多项式的值 pm=polyvalm(p,s)% 按矩阵规则计算给定s时多项式的值 [r,p,k]=residue(b,a) % 部分分式展开，b,a分别是分子，分母多项式系数矢量。r,p,k分别是留数，极点和值项矢量。 p=poly(x,y,n) % 用n介多项式拟合给定的数据 polyder(p) %多项式微分 x=[1 2 3 4 5]; y=[5.5 43.1 128 290.7 498.4]; p=polyfit(x,y,3) x2=1:0.1:5; y2=polyval(p,x2); plot(x,y,'o',x2,y2) 2，线性代数 1，求解方程的根 t=[0 0.3 0.8 1.1 1.6 2.3]'; y=[0.82 0.72 0.63 0.60 0.55 0.50]'; e=[ones(size(t)) exp(-t)] c=e\y t1=[0:0.1:2.5]'; y1=[ones(size(t1)),exp(-t1)]*c; plot(t1,y1,'b',t,y,'ro') 2，逆矩阵及行列式 inv(a) det(a) pinv(a) 3，矩阵分解（略） 4，数据分析 max(x)求x各列的最大元素 mean(x)求x各列的平均值 median(x)找出x各列的中位元素 min(x)求出x各列的最小元素 S=cumsum(x)求x各列元素累计和 sort(x)使x的各列元素按递增排序。 std(x)求x各列的标准差。 sum(x)求x各列元素之和 prod(x)求x各列元素之积

基于MATLAB的电力系统潮流计算

基于MATLAB的电力系统潮流计算 %简单潮流计算的小程序，相关的原始数据数据数据输入格式如下： %B1是支路参数矩阵，第一列和第二列是节点编号。节点编号由小到大编写%对于含有变压器的支路，第一列为低压侧节点编号，第二列为高压侧节点%编号，将变压器的串联阻抗置于低压侧处理。 %第三列为支路的串列阻抗参数。 %第四列为支路的对地导纳参数。 %第五烈为含变压器支路的变压器的变比 %第六列为变压器是否是否含有变压器的参数，其中“1”为含有变压器，%“0”为不含有变压器。 %B2为节点参数矩阵，其中第一列为节点注入发电功率参数；第二列为节点%负荷功率参数；第三列为节点电压参数；第六列为节点类型参数，其中 %“1”为平衡节点，“2”为PQ节点，“3”为PV节点参数。 %X为节点号和对地参数矩阵。其中第一列为节点编号，第二列为节点对地%参数。 n=input('请输入节点数:n='); n1=input('请输入支路数:n1='); isb=input('请输入平衡节点号:isb='); pr=input('请输入误差精度:pr='); B1=input('请输入支路参数:B1='); B2=input('请输入节点参数:B2='); X=input('节点号和对地参数:X='); Y=zeros(n); Times=1; %置迭代次数为初始值 %创建节点导纳矩阵 for i=1:n1 if B1(i,6)==0 %不含变压器的支路 p=B1(i,1); q=B1(i,2); Y(p,q)=Y(p,q)-1/B1(i,3); Y(q,p)=Y(p,q); Y(p,p)=Y(p,p)+1/B1(i,3)+0.5*B1(i,4); Y(q,q)=Y(q,q)+1/B1(i,3)+0.5*B1(i,4); else %含有变压器的支路 p=B1(i,1); q=B1(i,2); Y(p,q)=Y(p,q)-1/(B1(i,3)*B1(i,5)); Y(q,p)=Y(p,q); Y(p,p)=Y(p,p)+1/B1(i,3);

基于MATLAB的潮流计算源程序代码(优.选)

%*************************电力系统直角坐标系下的牛顿拉夫逊法潮流计算********** clear clc load E:\data\IEEE014_Node.txt Node=IEEE014_Node; weishu=size(Node); nnum=weishu(1,1); %节点总数 load E:\data\IEEE014_Branch.txt branch=IEEE014_Branch; bwei=size(branch); bnum=bwei(1,1); %支路总数 Y=(zeros(nnum)); Sj=100; %********************************节点导纳矩阵******************************* for m=1:bnum; s=branch(m,1); %首节点 e=branch(m,2); %末节点 R=branch(m,3); %支路电阻 X=branch(m,4); %支路电抗 B=branch(m,5); %支路对地电纳 k=branch(m,6); if k==0 %无变压器支路情形 Y(s,e)=-1/(R+j*X); %互导纳 Y(e,s)=Y(s,e); end if k~=0 %有变压器支路情形 Y(s,e)=-(1/((R+j*X)*k)); Y(e,s)=Y(s,e); Y(s,s)=-(1-k)/((R+j*X)*k^2); Y(e,e)=-(k-1)/((R+j*X)*k); %对地导纳 end Y(s,s)=Y(s,s)-j*B/2; Y(e,e)=Y(e,e)-j*B/2; %自导纳的计算情形 end for t=1:nnum; Y(t,t)=-sum(Y(t,:))+Node(t,12)+j*Node(t,13); %求支路自导纳 end G=real(Y); %电导 B=imag(Y); %电纳 %******************节点分类************************************* * pq=0; pv=0; blancenode=0; pqnode=zeros(1,nnum); pvnode=zeros(1,nnum); for m=1:nnum; if Node(m,2)==3 blancenode=m; %平衡节点编号 else if Node(m,2)==0 pq=pq+1; pqnode(1,pq)=m; %PQ 节点编号 else if Node(m,2)==2 pv=pv+1; pvnode(1,pv)=m; %PV 节点编号 end end end end %*****************************设置电压初值********************************** Uoriginal=zeros(1,nnum); %对各节点电压矩阵初始化 for n=1:nnum Uoriginal(1,n)=Node(n,9); %对各点电压赋初值 if Node(n,9)==0;

基于matlab--psat软件的电力系统潮流计算课程设计

东北电力大学课程设计改革试用任务书: 电力系统潮流计算课程设计任务书 设计名称：电力系统潮流计算课程设计 设计性质：理论计算，计算机仿真与验证 计划学时：两周 一、设计目的 1.培养学生独立分析问题、解决问题的能力； 2.培养学生的工程意识，灵活运用所学知识分析工程问题的能力 3.编制程序或利用电力系统分析计算软件进行电力系统潮流分析。 二、原始资料 1、系统图：IEEE14节点。 2、原始资料：见IEEE14节点标准数据库 三、课程设计基本内容： 1.采用PSAT仿真工具中的潮流计算软件计算系统潮流； 1)熟悉PSAT仿真工具的功能； 2)掌握IEEE标准数据格式内容； 3)将IEEE标准数据转化为PSAT计算数据； 2.分别采用NR法和PQ分解法计算潮流，观察NR法计算潮流中雅可比矩阵的变化情况， 分析两种方法计算潮流的优缺点； 3.分析系统潮流情况，包括电压幅值、相角，线路过载情况以及全网有功损耗情况。

4.选择以下内容之一进行分析： 1)找出系统中有功损耗最大的一条线路，给出减小该线路损耗的措施，比较各种措施 的特点，并仿真验证； 2)找出系统中电压最低的节点，给出调压措施，比较各种措施的特点，并仿真验证； 3)找出系统中流过有功功率最大的一条线路，给出减小该线路有功功率的措施，比较 各种措施的特点，并仿真验证； 5.任选以下内容之一作为深入研究：（不做要求） 1)找出系统中有功功率损耗最大的一条线路，改变发电机有功出力，分析对该线路有 功功率损耗灵敏度最大的发电机有功功率，并进行有效调整，减小该线路的损耗； 2)找出系统中有功功率损耗最大的一条线路，进行无功功率补偿，分析对该线路有功 功率损耗灵敏度最大的负荷无功功率，并进行有效调整，减小该线路的损耗； 3)找出系统中电压最低的节点，分析对该节点电压幅值灵敏度最大的发电机端电压， 并有效调整发电机端电压，提高该节点电压水平； 四、课程设计成品基本要求： 1.绘制系统潮流图，潮流图应包括： 1)系统网络参数 2)节点电压幅值及相角 3)线路和变压器的首末端有功功率和无功功率 2.撰写设计报告，报告内容应包括以下几点： 1)本次设计的目的和设计的任务； 2)电力系统潮流计算的计算机方法原理，分析NR法和PQ分解法计算潮流的特点； 3)对潮流计算结果进行分析，评价该潮流断面的运行方式安全性和经济性； 4)找出系统中运行的薄弱环节，如电压较低点或负载较大线路，给出调整措施； 5)分析各种调整措施的特点并比较它们之间的差异； 6)结论部分以及设计心得； 五、考核形式 1.纪律考核：学生组织出勤情况和工作态度等； 2.书面考核：设计成品的完成质量、撰写水平等； 3.答辩考核：参照设计成品，对计算机方法进行电力系统潮流计算的相关问题等进行答辩； 4.采用五级评分制：优、良、中、及格、不及格五个等级。

第06章_MATLAB数值计算_例题源程序汇总

第6章 MATLAB 数值计算 例6.1 求矩阵A 的每行及每列的最大和最小元素，并求整个矩阵的最大和最小元素。 1356 78256323578255631 01-???? -? ?=???? -??A A=[13,-56,78;25,63,-235;78,25,563;1,0,-1]; max(A,[],2) %求每行最大元素 min(A,[],2) %求每行最小元素 max(A) %求每列最大元素 min(A) %求每列最小元素 max(max(A)) %求整个矩阵的最大元素。也可使用命令：max(A(:)) min(min(A)) %求整个矩阵的最小元素。也可使用命令：min(A(:)) 例6.2 求矩阵A 的每行元素的乘积和全部元素的乘积。 A=[1,2,3,4;5,6,7,8;9,10,11,12]; S=prod(A,2) prod(S) %求A 的全部元素的乘积。也可以使用命令prod(A(:)) 例6.3 求向量X =(1！，2！，3！，…，10！)。 X=cumprod(1:10) 例6.4 对二维矩阵x ，从不同维方向求出其标准方差。 x=[4,5,6;1,4,8] %产生一个二维矩阵x y1=std(x,0,1) y2=std(x,1,1) y3=std(x,0,2) y4=std(x,1,2) 例6.5 生成满足正态分布的10000×5随机矩阵，然后求各列元素的均值和标准方差，再求这5列随机数据的相关系数矩阵。 X=randn(10000,5); M=mean(X) D=std(X) R=corrcoef(X)

例6.6 对下列矩阵做各种排序。 185412613713-?? ??=?? ??-?? A A=[1,-8,5;4,12,6;13,7,-13]; sort(A) %对A 的每列按升序排序 -sort(-A,2) %对A 的每行按降序排序 [X,I]=sort(A) %对A 按列排序，并将每个元素所在行号送矩阵I 例6.7 给出概率积分 2 (d x x f x x -? e 的数据表如表6.1所示，用不同的插值方法计算f (0.472)。 x=0.46:0.01:0.49; %给出x ，f(x) f=[0.4846555,0.4937542,0.5027498,0.5116683]; format long interp1(x,f,0.472) %用默认方法，即线性插值方法计算f(x) interp1(x,f,0.472,'nearest') %用最近点插值方法计算f(x) interp1(x,f,0.472,'spline') %用3次样条插值方法计算f(x) interp1(x,f,0.472,'cubic') %用3次多项式插值方法计算f(x) format short 例6.8 某检测参数f 随时间t 的采样结果如表6.2，用数据插值法计算t =2，7，12，17，22，17，32，37，42，47，52，57时的f 值。 T=0:5:65; X=2:5:57;

MATLAB下的潮流计算实现-稀疏技术毕业设计

毕业设计（论文）MATLAB下的潮流计算实现-稀疏技术

毕业设计（论文）原创性声明和使用授权说明 原创性声明 本人郑重承诺：所呈交的毕业设计（论文），是我个人在指导教师的指导下进行的研究工作及取得的成果。尽我所知，除文中特别加以标注和致谢的地方外，不包含其他人或组织已经发表或公布过的研究成果，也不包含我为获得及其它教育机构的学位或学历而使用过的材料。对本研究提供过帮助和做出过贡献的个人或集体，均已在文中作了明确的说明并表示了谢意。 作者签名：日期： 指导教师签名：日期： 使用授权说明 本人完全了解大学关于收集、保存、使用毕业设计（论文）的规定，即：按照学校要求提交毕业设计（论文）的印刷本和电子版本；学校有权保存毕业设计（论文）的印刷本和电子版，并提供目录检索与阅览服务；学校可以采用影印、缩印、数字化或其它复制手段保存论文；在不以赢利为目的前提下，学校可以公布论文的部分或全部内容。 作者签名：日期：

摘要 电力系统潮流计算是研究电力系统稳态运行情况的一种计算，它根据给定的运行条件及系统接线情况确定整个电力系统各部分的运行状态：各母线的电压，各元件中流过的功率，系统的功率损耗等等。在电力系统规划的设计和现有电力系统运行方式的研究中，都需要利用潮流计算来定量地分析比较供电方案或运行方式的合理性、可靠性和经济性。因此潮流计算是研究电力系统的一种很重要和很基础的计算。由于电力系统结构及参数的一些特点，并且随着电力系统不断扩大，潮流问题的方程式阶数越来越高，对这样的方程式并不是任何数学方法都能保证给出正确答案的。这种情况成为促使电力系统计算人员不断寻求新的更可靠方法的重要因素。 本文旨在于研究潮流计算的牛顿—拉夫逊法的基本原理，在Matlab环境中实现牛顿—拉夫逊法潮流计算的数学模型，程序流程以及编制相应程序，并在程序中融合了节点优化编号和稀疏技术，以提高计算效率。最后用IEEE-3O节点标准测试系统验证所编程序。 关键词：潮流计算Newtom-Raphson法节点优化稀疏技术Matlab ABSTRACT Power flow calculation is fundanmental of analysis. Network reconfiguration,fault management,state estimator etc also need the data of electrial system power flow.There is important significance to develop power flow calculation in allusion to traits of distribution network. This paper introduces the principle of Newtom-Raphson algorithm, which is developed for calculation of power flow calculation ,where zero sequence network is open.With this algorithm,the three-phase load is resolved into positive/negative sequence power and coupling power,thus,decoupling three phase power flow into sequencet component power flow.The power flow can be obtained by just finding the positive sequence power flow and then finding the negative sequent component from the coupling https://www.doczj.com/doc/6710715877.html,pared with the existing methods,the jacobian matrix with the proposed algorithm is of much lower order,thus substantially reducing the computation burden.The proposed algorithm,together with a reference algorithm,has been simulated on an actual IEEE-30 system using statistic load date.And then it will

数值分析的MATLAB程序

列主元法 function lianzhuyuan(A,b) n=input('请输入n:') %选择阶数A=zeros(n,n); %系数矩阵A b=zeros(n,1); %矩阵b X=zeros(n,1); %解X for i=1:n for j=1:n A(i,j)=(1/(i+j-1)); %生成hilbert矩阵A end b(i,1)=sum(A(i,:)); %生成矩阵b end for i=1:n-1 j=i; top=max(abs(A(i:n,j))); %列主元 k=j; while abs(A(k,j))~=top %列主元所在行 k=k+1; end for z=1:n %交换主元所在行a1=A(i,z); A(i,z)=A(k,z); A(k,z)=a1; end a2=b(i,1); b(i,1)=b(k,1); b(k,1)=a2; for s=i+1:n %消去算法开始m=A(s,j)/A(i,j); %化简为上三角矩阵 A(s,j)=0; for p=i+1:n A(s,p)=A(s,p)-m*A(i,p); end b(s,1)=b(s,1)-m*b(i,1); end end X(n,1)=b(n,1)/A(n,n); %回代开始 for i=n-1:-1:1 s=0; %初始化s for j=i+1:n s=s+A(i,j)*X(j,1);

end X(i,1)=(b(i,1)-s)/A(i,i); end X 欧拉法 clc clear % 欧拉法 p=10; %贝塔的取值 T=10; %t取值的上限 y1=1; %y1的初值 r1=1; %y2的初值 %输入步长h的值 h=input('欧拉法please input number(h=1 0.5 0.25 0.125 0.0625):h=') ; if h>1 or h<0 break end S1=0:T/h; S2=0:T/h; S3=0:T/h; S4=0:T/h; i=1; % 迭代过程 for t=0:h:T Y=(exp(-t)); R=(1/(p-1))*exp(-t)+((p-2)/(p-1))*exp(-p*t); y=y1+h*(-y1); y1=y; r=r1+h*(y1-p*r1); r1=r; S1(i)=Y; S2(i)=R; S3(i)=y; S4(i)=r; i=i+1; end t=[0:h:T]; % 红线为解析解，'x'为数值解 plot(t,S1,'r',t,S3,'x')

Matlab关于数值计算的实现

Matlab关于数值计算的实现 摘要：数值计算（numerical computation computation），主要研究更好的利用计算机更好的进行数值计算，解决各种数学问题。数值分析包括离散傅里叶变换，考虑截断误差，计算误差，函数的敛散性与稳定性等。在数学方面，数值计算的主要研究数值微分与积分，数据的处理与多项式计算，最优化问题，线性方程与非线性方程的求解，常微分方程的数值求解等。同时，数值计算在物理，化学，经济等方面也有研究，本文暂且不表。M atlab软件历经二十多年来的发展，已成为风靡世界的数学三大软件（matlb，Mathematica l，Maple）之一，在数学类科技应用软件中在数值计算方面首屈一指。Matlab以矩阵为数据操作的基本单位，使得矩阵运算十分便捷快速，同时Matlab还提供了海量的计算函数，而且使用可靠地算法进行计算，能使用户在繁复的数学运算中解脱，Matlab还具有方便且完善的图形处理功能，方便绘制二维和三维图形并修饰。

目录 1.数值分析(离散傅里叶变换，考虑截断误差，计算误差，函数 的敛散性与稳定性) 2.数值计算(数值微分与积分，数据的处理与多项式计算， 最优化问题，线性方程与非线性方程的求解，常微分方程的数值求解) 3.图形处理功能(方便绘制二维和三维图形并修饰) 4.总结

1.数据统计与分析 Matlab 可以进行求矩阵的最大最小元素，平均值与中值，关于矩阵元素的求和与求积，累加和与累乘积，标准方程，相关系数，元素排序。现在以求标准方差举例说明Matlab 的实现。 在Matlab 中，实现标准方差计算的函数为std 。对于向量(Y ），std （Y ）实现返回一个标准方差，而对于矩阵（A ），std （A ）返回一个行向量，该行向量的每个元素对应着矩阵A 各行或各列的标准方差。一般调用std 函数的格式为std （A ，flag ，dim ） Dim 取1或者2分别对应求各列或各行的标准方差，flag 取1时，按照标准方差的计算公式 ∑-=-=N i x x S i N 1 2 1)(11来计算。若flag 取2，则用公式 ∑-==N i x x S i N 1 2 2) (1 进行计算。默认的flag 取值为0，dim 取值为1。课本page143 2. 离散傅里叶变换 离散傅里叶变换广泛应用于信号的分析，光谱和声谱分析、全息技术等各个领域。但直接计算dft 的运算量与变化的长度N 的平方成正比，当N 较大时，计算量太大。随着计算机技术的迅速发展，在计算机上进行离散傅里叶变换计算成为可能。特别是快速傅里叶变换算法的出现，为傅里叶变换的应创造了条件。 （1）：傅里叶变换算法的简述。 傅立叶变换是一种分析信号的方法，它可分析信号的成分，也可用这些成分合成信号。许多波形可作为信号的成分，比如正弦波、方波、锯齿波等，傅立叶变换用正弦波作为信号的成分. f(t)是t 的周期函数，如果t 满足狄里赫莱条件：在一个以2T 为周期内f(X)连续或只有有限个第一类间断点，附f （x ）单调或可划分成有限个单调区间，则F （x ）以2T 为周期的傅里叶级数收敛，和函数S （x ）也是以2T 为周期的周期函数，且在这些间断点上，函数是有限值；在一个周期内具有有限个极值点；绝对可积。则有下图①式成立。称为积分运算f(t)的傅立叶变换， ②式的积分运算叫做F(ω)的傅立叶逆变换。F(ω)叫做f(t)的像函数，f(t)叫做 F(ω)的像原函数。F(ω)是f(t)的像。f(t)是F(ω)原像。 ①傅立叶变换 ②傅立叶逆变换

Matlab的gui界面设计实例练习

一个不错的Matlab的gui界面设计实例 %非常漂亮的日历， function CalendarTable; % calendar 日历 % Example: % CalendarTable; S=datestr(now); [y,m,d]=datevec(S); % d is day % m is month % y is year DD={'Sun','Mon','Tue','Wed','Thu','Fri','Sat'}; close all figure; for k=1:7; uicontrol(gcf,'style','text',... 'unit','normalized','position',[0.02+k*0.1,0.55,0.08,0.06],... 'BackgroundColor',0.6*[1,1,1],'ForegroundColor','b',... 'String',DD(k),'fontsize',16,'fontname','times new roman'); end h=1; ss='b'; qq=eomday(y,m); for k=1:qq; n=datenum(y,m,k); [da,w] = weekday(n); if k==d; ss='r'; end uicontrol(gcf,'style','push',... 'unit','normalized','position',[0.02+da*0.1,0.55-h*0.08,0.08,0.06],... 'BackgroundColor',0.6*[1,1,1],'ForegroundColor',ss,... 'String',num2str(k)); ss='b'; if da==7; h=h+1;

Matlab牛拉法潮流计算程序

%本程序的功能是用牛顿——拉夫逊法进行潮流计算 % B1矩阵：1、支路首端号；2、末端号；3、支路阻抗；4、支路对地电纳 % 5、支路的变比；6、支路首端处于K侧为1，1侧为0 % B2矩阵：1、该节点发电机功率；2、该节点负荷功率；3、节点电压初始值 % 4、PV节点电压V的给定值；5、节点所接的无功补偿设备的容量 % 6、节点分类标号：1为平衡节点（应为1号节点）；2为PQ节点； % 3为PV节点； clear; n=input('请输入节点数:n='); nl=input('请输入支路数:nl='); isb=input('请输入平衡母线节点号:isb='); pr=input('请输入误差精度:pr='); B1=input('请输入由各支路参数形成的矩阵:B1='); B2=input('请输入各节点参数形成的矩阵:B2='); Y=zeros(n);e=zeros(1,n);f=zeros(1,n);V=zeros(1,n);sida=zeros(1,n);S1=zeros(nl); % % %--------------------------------------------------- for i=1:nl %支路数 if B1(i,6)==0 %左节点处于1侧 p=B1(i,1);q=B1(i,2); else %左节点处于K侧 p=B1(i,2);q=B1(i,1); end Y(p,q)=Y(p,q)-1./(B1(i,3)*B1(i,5)); %非对角元 Y(q,p)=Y(p,q); %非对角元 Y(q,q)=Y(q,q)+1./(B1(i,3)*B1(i,5)^2)+B1(i,4)./2; %对角元K侧 Y(p,p)=Y(p,p)+1./B1(i,3)+B1(i,4)./2; %对角元1侧 end %求导纳矩阵 disp('导纳矩阵Y='); disp(Y) %---------------------------------------------------------- G=real(Y);B=imag(Y); %分解出导纳阵的实部和虚部 for i=1:n %给定各节点初始电压的实部和虚部 e(i)=real(B2(i,3)); f(i)=imag(B2(i,3)); V(i)=B2(i,4); %PV节点电压给定模值 end for i=1:n %给定各节点注入功率 S(i)=B2(i,1)-B2(i,2); %i节点注入功率SG-SL B(i,i)=B(i,i)+B2(i,5); %i节点无功补偿量 end %=================================================================== P=real(S);Q=imag(S); %分解出各节点注入的有功和无功功率 ICT1=0;IT2=1;N0=2*n;N=N0+1;a=0; %迭代次数ICT1、a；不满足收敛要求的节点数IT2

MATLAB实验二MATLAB的数值运算和程序

课程名称：Matlab语言 开设时间：2016—2017学年第 2 学期 专业班级：学生学号：学生姓名： 实验名称：实验二、MATLAB的数值运算和程序实验成绩： 指导教师：批改时间： 一、实验目的和要求 1）掌握基本的矩阵运算及常用的函数。 2）掌握MATLAB函数的编写及调试方法。 3）掌握MATLAB常用的数值运算函数。 二、实验仪器和设备 计算机一台 三、实验过程 1、一维数组在命令窗口执行下面指令，观察输出结果，体味数组创建和寻访方法，%号后面的为注释，不用输入。 rand('state',0) % 把均匀分布伪随机发生器置为0 状态 x=rand(1,5) % 产生（1*5）的均布随机数组 x(3) % 寻访数组x 的第三个元素。 x([1 2 5]) % 寻访数组x 的第一、二、五个元素组成的子数组。 x(1:3) % 寻访前三个元素组成的子数组 x(3:end) % 寻访除前2 个元素外的全部其他元素。end 是最后一个元素的下标。 x(3:-1:1) % 由前三个元素倒排构成的子数组 x(find(x>0.5)) % 由大于0.5 的元素构成的子数组 x([1 2 3 4 4 3 2 1]) % 对元素可以重复寻访，使所得数组长度允许大于原数组。 x(3) = 0 % 把上例中的第三个元素重新赋值为0。 x([1 4])=[1 1] % 把当前x 数组的第一、四个元素都赋值为1。 x[3]=[] % 空数组的赋值操作

2、在命令窗口执行下面指令，观察输出结果 a=2.7358; b=33/79; % 这两条指令分别给变量 a ， b 赋值。 C=[1,2*a+i*b,b*sqrt(a);sin(pi/4),a+5*b,3.5+i] % 这指令用于创建二维组C M_r=[1,2,3;4,5,6],M_i=[11,12,13;14,15,16] % 创建复数数组的另一种方法 CN=M_r+i*M_i % 由实部、虚部数组构成复数数组 3. 记录下面题目的程序和运行后的结果。 1??????=654321a ??????-=531142b ???? ? ?????-=201c ??????????=063258741d 下列运算是否合法，为什么？如合法，结果是多少？