第三章线性系统状态方程的解

第三章线性系统状态方程的解
第三章线性系统状态方程的解

第三章 线性系统的运动分析

§3-1线性连续定常齐次方程求解

一、齐次方程和状态转移矩阵的定义

1、齐次方程

状态方程的齐次方程部分反映系统自由运动的状况(即没有输入作用的状况),设系统的状态方程的齐次部分为:)()(t Ax t x =

线性定常连续系统:Ax x =

2、状态转移矩阵的定义

齐次状态方程Ax x = 有两种常见解法:(1)幂级数法;(2)拉氏变换法。其解为)0()(x e

t x At

?=。

其中At

e 称为状态转移矩阵(或矩阵指数函数、矩阵指数),记为:At

e t =)(φ。

若初始条件为)(0t x ,则状态转移矩阵记为:)(00

)(t t A e t t -=-Φ

对于线性时变系统,状态转移矩阵写为),(0t t φ,它是时刻t ,t 0的函数。但它一般不能写成指数形式。

(1)幂级数法

设Ax x

= 的解是t 的向量幂级数 +++++=k

k t b t b t b b t x 2210)(

式中 ,,,

,,k b b b b 210都是n 维向量,则

+++++=-1232132)(k k t kb t b t b b t x

)(2210 +++++=k

k t b t b t b b A

故而有:

?????

??

????==

====003

230

2

12

01!1!

3131

2

121b A k b b A Ab b b A Ab b Ab b K K

且有0)0(b x =。 故

+++++=k k t b t b t b b t x 2210)( ++

+++=k

k t

b A k t b A t Ab b 02

02

00!

1!

21

)0()!

1!21(22

x t A k t A At I k

k

++

++

+=

定义:∑

==

++

+++=0

2

2!

1!

1!21K k k k

k At

t A k t A k t A At I e

则)0()(x e t x At ?=。

(2)拉氏变换解法

将Ax x

= 两端取拉氏变换,有 )()0()(s Ax x s sx =- )0()()(x s x A sI =- )0()()(1x A sI s x ?-=- 拉氏反变换,有

)0(])[()(11x A sI L t x ?-=-- 则

])[()(11---==A sI L e t At φ

【例3.1.1】 已知系统的状态方程为x x

??

?

???=00

10 ,初始条件为)0(x ,试求状态转移矩阵和状态方程的解。 解:(1)求状态转移矩阵 ++

++

+==k

k At

t

A k t A At I e t !

1!

21)(2

此题中:

????

??=00

10A , ??

?

???====00

0032n

A A A 所以

??

????=??????+??????=+==1010001001)(t t At I e t At φ (2)状态方程的解 )0(101

)0()(x t x e t x At

???

?

?

??=?=

【例 3.1.2】 已知系统状态方程为x x

??

?

???--=3210 ,初始条件为)0(x ,试求状态方程的解。

解:)0()(x e t x At ?=

??

?

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

??=-32

132

100

0s s s s A sI ??

??

?

?????++

+-++

+-+-

++-+=??

????-+++=--22112

21221

112

1122

13

)2)(1(1

)

(1

s s s s s s s s s s s s A sI ∴?

?

????+-+---=-==----------t t

t

t t

t t

t At

e e

e

e e e e e A sI L e t 22221

1

2222])[()(φ 故而

)0(2222)0()(2222x e e

e

e e

e e e x e

t x t t

t

t t

t t

t At

?

?

?

???+-+---=?=-------- 二、状态转移矩阵At e 的性质

++

++

+==k

k At

t A k t A At I e

t !

1!

21)(2

(1)I =)0(φ

(2)A t t A t )()()(φφφ== A =)0(φ

(3))()()()()(122121t t t t t t φφφφφ±=±=± 证明:)()()()()(1221)

()

()

(212121t t t t e

e

e

t t t A t A t t A φφφφφ±=±=?==±±±

(4))()(1t t -=-φφ,)()(1t t φφ=--

证明:)()()()()()0(1t t I t t t t -=?=-=-=-φφφφφφ (5))()()(00t x t t t x -=φ 证明:)0()()(x t t x φ=

)()()0()0()()(00100t x t x x t t x ?=?=-φφ,代入上式 ∴)()()()()()(00001t x t t t x t t t x -=?=-φφφ 证毕。

(6))()()(011202t t t t t t --=-φφφ

证明:)()()(0022t x t t t x -=φ………………………. …………………(1) )()()(0011t x t t t x -=φ……………………………………………(2) )()()()()()(001121122t x t t t t t x t t t x --=-=φφφ…………….(3) 比较(1)、(3)式,有)()()(011202t t t t t t --=-φφφ成立。证毕。

(7)[])()(kt t k

φφ=

证明:[])(][)()(kt e e e t kt A kAt k At k

φφ====

(8)若BA AB =,则At

Bt

Bt

At

t

B A e

e

e

e

e

?=?=+)(

若BA AB ≠,则At

Bt

Bt

At

t

B A e e e e e

?≠?≠+)(

(9)设)(t φ为Ax x

= 的状态转移矩阵,引入非奇异变换x P x =后的状态转移矩阵为: P e P t At

1)(-=φ

证明:将x P x =代入Ax x

= 中,有 x AP P x

1-= APt

P e t 1

)(-=φ

++

++

+=----k

k APt

P

t

AP P

k t AP P

APt P

I e

)(!

1)(!

211

2

2

1

1

1

++

+++=----k

k t

AP P

k t AP P

APt P P P )(!

1)(!

211

2

21

11

P t

A k t A At I P k

k

)!1!

21(2

2

1 ++

+++=-

P e P At 1-=

∴P e P t At 1)(-=φ。证毕。

(10)两种常见的状态转移矩阵

①设],,,[21n diag A λλλ =,即A 为对角阵,且具有互异元素。则

??

??

?

???

?

?

?

?=t t n e e t λλφ

00)(1

②设A 为m m ?约当阵

m

m A ??????????

??

?=λλ

λ11

,则???

???

?

??

??????????

?--=--t t m t

t t

m t

t

t

e e

t m te e e

t

m e

t te e

t λλλλλλλλφ

0)!

2(10)!

1(1

!21)(21

2

【例3.1.3】 已知状态转移矩阵为 ?

?

????+-+---=--------t t

t

t t

t t

t At

e e

e

e e

e e e e

22222222 试求)(1

t -φ

和A 。

解:(1)根据状态转移矩阵的性质4,可知 ?

?

????+-+---=-=-t t t

t t

t t

t e e e

e e

e e e t t 22221

2222)()(φφ

(2)根据状态转移矩阵的性质2,可知 ??

?

???--==??????--+-+-==--------3210

0442222)0(2222t e

e

e

e e

e e e A t

t

t

t

t

t

t

t φ

【例3.1.4】 已知

4

40

11

01

??????????

?

??=λλ

λλA 试求状态转移矩阵At e 。 解:根据状态转移矩阵的性质10,可知

??

??????????????==t

t t

t t

t

t

t

t

t At

e te e e t te e e

t e

t te e e t λλλλλλλλλλφ0

002106

1

21)(23

2

【例3.1.5】 验证如下矩阵是否为状态转移矩阵。 ????

?

???

??-t t

t t sin cos 0

cos sin 0

001

解:利用性质(1)I =)0(φ I t t t

t t ≠???

?

?

?????-==??????????-01

100

001

0sin cos 0

cos sin 0

001

,所以该矩阵不是状态转移矩阵。

【例3.1.6】 已知系统状态方程为Ax x

= , 当???

???-=11)0(x 时,?

?

????-=--t t e e t x 22)( 当???

???-=12)0(x 时,??

????-=--t t e e t x 2)( 试求系统矩阵A 和状态转移矩阵At

e 。

解:由性质(2)可知:)0(φ =A 由已知,有 )0()(x e

t x At

?= ??

????--=?

?????--?----1121222At t t

t

t

e e e

e

e

∴??

?

???--??????--=??

????--??????--=---------1121211212221

22t t

t

t t t

t

t At

e e

e

e e e

e

e e

?

?

????+-+---=--------t t

t

t t

t t

t e e

e

e e

e e e 22222222

∴??

?

???--=??????--+-+-===--------=3120

424222)(0

22220

t t

t t

t

t

t

t t t e

e

e

e e

e e e t A φ

§3-2 线性连续定常非齐次状态方程的解

线性定常非齐次状态方程:Bu Ax x

+= ,求)(t x 。 1、直接积分法

Bu Ax x

+= 左乘At e -,有 Bu e Ax x

e At At ?=---)( 由于)()(Ax x

e x e dt

d At

At

-=?-- 所以Bu e

x e

dt

d At

At ?=?--)(,两端同时积分,有

τττ

d Bu e

x t x e

A t

At

)()0()(0

?=

---?

∴τττd Bu e

x e t x t A t

At

)()0()()

(0

?+

=-?

τττφφd Bu t x t t

)()()0()(0

?-+

=?

注意:若取0t 作为初始时刻,积分可得: τττ

d Bu e

t x e t x e A t

t At At )()()(0

0?=

----?

τττd Bu e

t x e t x t A t

t t t A )()()()

(0)

(0

0?+

=--?

2、拉氏变换法

Bu Ax x

+= ,两边同时取拉氏变换 )()()0()(s Bu s Ax x s sx +=- )()0()()(s Bu x s x A sI +=-

则 )()()0()()(11s Bu A sI x A sI s x ---+-=

)]()[()0(])[()(1111s Bu A sI L x A sI L t x -----+-= 由拉氏变换卷积定理: ?

-=

-t

d f t f s F s F L 0

21211)().()]().([τττ

在此1)(--A sI 视为)(1s F ,)(s Bu 视为)(2s F 。则 τττd Bu e

x e t x t A t

At )()0()()

(0

?+=-?

【例3.2.1】 已知系统状态方程为u x x

??

?

???+??????--=1032

10 ,输入)(1)(t t u =, 初始条件为??

?

???=)0()0()0(21x x x ,试求解此非齐次状态方程。

解:由已知有

τττd Bu e

x e t x t A t

At

)()0()()

(0

?+

=-?

(1)先求At e ,由前面例题可知 ?

?

????+-+---=--------t t

t

t t

t t t At

e e

e

e e

e e e e

22222222 (2)求τττd Bu e

t A t )()

(0

?-?

ττττττττττττd e e

e e

e

e e e d Bu e

t

t t t t t t t t t A t

?

?

?

??????????+-+---=???-----------------102222)(0

)(2)

()(2)()

(2)()

(2)()

(0

τττττd e e e e t t t t t

?

?

?

???+--=

--------?

)(2)()(2)(0

2 )(2)(2)()(2)(0

τττττ-?

?

?

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

t d e e e e t t t t t

???

?????-+-=????????-+--=------------t

t t t t t t t e e e e t e e e e 22)(2)()(2)(2121021ττττ 故而

???

?????-+-+?

????????

??

??+-+---=------------t

t t t

t t

t

t t

t t

t e e e e x x e e

e

e e

e e e t x 222122222121)0()0(2222)(

特别说明:若?

?????=??????=00)0()0()0(21x x x ,则???

????

?-+-=----t t t t e

e e e t x 222121)( 其状态轨迹图可以MABLAB 绘出:

%Example 3.2.1 matlab program: grid;

xlabel('时间轴');

ylabel('x 代表x1,----*代表x2'); t=0:0.1:10;

x1=0.5-exp(-t)+0.5*exp(-2*t); x2=exp(-t)-exp(-2*t); plot(t,x1,'x',t,x2,'*')

end

§3-3 状态转移矩阵At e 的计算

1、直接幂级数法

==

++

++

+=0

2

2!

1!

1!

21K k k k

k At

t A k t A k t A At I e

2、拉氏变换法

])[(11---=A sI L e At

3、利用性质,采用对角化的方法

【例3.3.1】 已知系统状态方程为x x

??

????--=32

10 ,试利用对角化的方法求At

e 。

解:0)2)(1()det(=++=-λλλA I ,解出特征值11-=λ,22-=λ。 选用变换阵P ,使AP P 1-对角化。由于A 为友矩阵,故P 可选为: ?

?????--=????

??=21

1111

21

λλP , ??

?

???--=-11

121

P 根据P e

P e

At

APt

P

1

1

-=-可推出:1

1

--=P

Pe

e APt

P

At

而??

?

???==--??

????---t t t APt

P

e e e

e

220

010

01

∴?

?

????+-+---==----------t t

t

t t

t t t APt

P

At

e e

e

e e

e e e P

Pe

e 22221

22221

4、利用Caylay-Hamilton 定理计算(待定系数法)

(1)Caylay-Hamilton 定理

设n 阶矩阵A 的特征多项式为: 011

1)(a a a A I f n n n ++++=-=--λλλλλ

则A 满足其特征方程,即

0)(011

1=++++=--I a A a A

a A A f n n n (2)推论1

矩阵A 的k (n k ≥)次幂,可表示为A 的)1(-n 阶多项式

∑-==

1

n m m

m

k

A

a

A ,n k ≥

【例如】 已知??

??

??=1021A ,求?100

=A

解:A 的特征多项式为:

12)(2

+-=-=λλλλA I f

根据Caylay-Hamilton 定理,有

02)(2=+-=I A A A f , ∴I A A -=22

故I A A I A A A I A A AA A 23)2(22)2(223-=--=-=-== I A A I A A A I A A AA A 342)2(323)23(234-=--=-=-== 依次归纳,有:

I k kA A k )1(--=

所以有:??

?

???=??????-????

??=-=102001990099

1000200100

99100100

I A A

(3)推论2

状态转移矩At e 可表示为A 的)1(-n 阶多项式

∑-==

1

)(n m m

m

At

A

t a

e

式中,)(,),(),(110t a t a t a n - 均为幂函数。

【例3.3.2】 已知系统状态方程为x x

??

?

???--=32

10 , 试利用Caylay-Hamilton 定理求At

e 。

解:(1)求系统矩阵A 的特征值

0)d e t (=-A I λ 0)2)(1(=++?λλ, 解出11-=λ,22-=λ

(2)一般情况下,对于n 个互异的特征值n λλλ,,

, 21,写出如下方程组: ??

????

?=++++=++++=++++------t

n n n n n t n n t n n n e a a a a e a a a a e a a a a λλλλλλλλλλλλ112210

1

212222101112121102

1

并解出n a a a ,,, 10即可。对于本例:

???=+=+t t e

a a e a a 21210110λλλλ ???=-=-?--t

t e a a e a a 210102 解出t t e e a 202---=,t t e e a 21---=

(3)对于系统具有n 个互异的特征值n λλλ,,, 21的情况,按下式计算At

e

112210--++++=n n At A a A a A a I a e 对于本例有: ?

?

?

??

?+-+---=+=--------t t

t

t t

t t

t At

e e

e e e

e e e A a I a e 2222102222

§3-4 离散系统状态方程的解

一、由差分方程建立动态方程

线性离散系统的动态方程可以充分利用差分方程建立,也可以利用线性连续动态方程的离散化得到。

SISO 线性定常离散系统的差分方程一般形式为:

)

()1()1()()()1()1()(011011k u b k u b n k u b n k u b k y a k y a n k y a n k y n n n ++++-+++=++++-+++--

式中,k 表示kT 时刻;T 为采样周期;y(k)、u (k)分别为kT 时刻的输出量和输入量;i a 、i

b (n i ,,,,

210=, 且1=n a )为表征系统特征的常数。 考虑初始条件为零时的Z 变换关系有:

)()]([z y k y Z =, )()]([z y z n k y Z n

=+ 对上边式子两边取Z 变换,并整理为:

11

10

11

1)

()()(a z a z

a z

b z b z b z b z u z y z G n n n n n n

n ++++++++=

=

----

11

10

11

1a z a z a z z z

b n n n

n n n ++++++++

=---- βββ

)

()(z D z N b n +

=

按连续系统的方法,对)(/)(z D z N 做串联分解,最后可得到离散系统状态空间表达式的一种形式:

)(100

0)()()()(100001000

010)1()1()1()1(12112

1

0121k u k x k x k x k x a a a a

k x k x k x k x n n n n n ??

???

???

????????+?????????????????????????????

??

?----=????????????

????++++---

[])()()(12

10k u b k x k y n n +=-ββββ

简记为: ???+=+=+)

()()()()()1(k du k cx k y k hu k Gx k x

MIMO 线性定常离散系统的动态方程为: ??

?+=+=+)()()()()()1(k Du k Cx k y k Hu k Gx k x

离散系统的一般结构图

【例3.4.1】设某线性离散系统的差分方程为:

)(2)1()(16.0)1()2(k u k u k y k y k y ++=++++ 试写出系统的状态空间表达式。 解:离散系统的状态空间表达式为: ??

?+=+=+)

()()()()()1(k du k cx k y k hu k Gx k x

其中:??????--=116

.010

G , ??

?

???=10h , []12=c , 0=d

二、线性定常连续系统动态方程的离散化

线性定常非齐次状态方程Bu Ax x

+= 在)(0t x 及)(t u 作用下的解为: τττd Bu e

t x e t x t A t

t t t A )()()()

(0)(0

?+

=--?

或τττφφd Bu t t x t t t x t

t )()()()()(0

00?-+

-=?

令kT t =0,则)()()(0k x kT x t x ==

T k t )1(+=,则)1(])1[()(+=+=k x T k x t x 常数=+=)1()(k u k u ,于是 )(])1[()(])1[()1()1(k u Bd T k k x kT T k k x T

k kT

??-++-+=+?

+ττφφ

)(])1[()()()1(k u Bd T k k x T T

k kT

??-++=?

+ττφφ

记ττφBd T k T H T

k kT

?-+=

?

+)1(])1[()(

令ττ'=-+T k )1(,则代换后有 ττφττφBd Bd T H T

T

?

?

=

''=

)()()(

故离散化状态方程为:)()()()()1(k u T H k x T G k x +=+ 输出方程为:)()()(k Du k Cx k y +=

【例3.4.2】试写出连续时间系统

u x x

??

?

???+??????-=1020

10 采样周期为T 的离散化状态方程。 解:先求At e

1

11120

1])[()(----??

?

?

??

+-=-==s s L A sI L e

t At

φ

???

?

????

-=?????

?

?????

?++=---t

t

e

e s s s s L 2210

)1(2

11

210)2(1

1

T

t t T T G ===)

()()(φφ???

?????

-=--T

T

e

e 220)1(2

11

τττφτ

τ

d e

e Bd T H T

T

???

??????

?????-=

=

?

?

--100

)1(211)()(0

220

ττ

τ

d e e T

?

???

?

???

?-=--0

22)1(21

021412122T

e

e ???????

??

?-+=--τ

ττ????

?

????

?-+-=--T T e

e T 2221214

14121 所以:

)()()()()1(k u T H k x T G k x +=+ )(21214

14121

)()(0

)1(2

11)1()1(2221

2221k u e e

T k x k x e

e k x k x T T T

T

???

?

?

????

?-+-+??

???????

????

?

-=?

??

???++----

例3.4.2连续系统离散化MATLAB 程序:

%Example3.4.2 : Continuous to discrete system A=sym('[0,1;0,-2]') B=sym('[0;1]') T='T'

[G,H]=c2d(A,B,T)

%example3.4.2的另一种MA TLAB 程序: syms s t T;

A=sym('[0,1;0,-2]'); B=sym('[0;1]'); I=eye(2); L=inv(s*I-A) lap=ilaplace(L)

G=subs(lap,'T')

H=int(symmul(lap,B),0,T)

三、离散系统状态方程的解

两种解法:递推法和Z 变换法。

递 推 法:又称迭代法,对于定常和时变系统都适用。 Z 变换法:只适用于定常系统。

1、递推法

)()()()()1(k u T H k x T G k x +=+ 依次令 ,2,1,0=k ,从而有

0=k )0()()0()()1(u T H x T G x += 1=k )1()()1()()2(u T H x T G x +=

)1()()0()()()0()(2u T H u T H T G x T G ++=

2=k )2()()2()()3(u T H x T G x +=

)2()()1()()()0()()()0()(2

3u T H u T H T G u T H T G x T G +++= ………… 依此类推。

递推公式为:

∑-=--+

=1

1)()()()0()()(k i i

k k

i u T H T G

x T G k x

其中)(T G k

称为线性定常离散系统的状态转移矩阵,记为)(k φ。 )()()(k kT T G k

φφ==

()(k φ满足:)()1(k G k φφ=+; I =)0(φ)

【例3.4.3】已知某离散系统的状态方程是: )()()()()1(k u T H k x T G k x +=+

????

??

--=116

.010

G ,??????=11H ,初始状态??

?

???-=11)0(x ,1)(=k u , 试用递推法求解)(k x 。

解:)0()()0()()1(u T H x T G x +=??????+??????-????

??

--=1111116

.010

???

???=84.10 )1()()1()()2(u T H x T G x +=??????+??????????

??

--=1184.10116.010

??????-=84.084.2 )2()()2()()3(u T H x T G x +=??????+??????-????

??

--=1184.084.2116

.010

??

?

?

??=386.116.0

显然,用递推法求解所得到的不是一个封闭的解析形式,而是一个解序列。

采用MATLAB 语言,求解例3.4.3:

%Example 3.4.3 G=[0,1;-0.16,-1]; H=[1;1]; U=1; X1=[1;-1]; hold on; for k=1:400 X1=G*X1+H*U plot(X1(1),X1(2),'*'); end

2、Z 变换法

设定常离散系统的状态方程是: )()()1(k Hu k Gx k x +=+ 两边取Z 变换:

)()()0()(z Hu z Gx zx z zx +=-,整理有 )()0()()(z Hu zx z x G zI +=-

∴)()()0()()(11z Hu G zI zx G zI z x ---+-= 两边取Z 反变换: )]()

[()]0()

[()(1

1

1

1

z Hu G zI Z

zx G zI Z k x -----+-=

【例3.4.4】已知某离散系统的状态方程是: )()()()()1(k u T H k x T G k x +=+

????

??

--=116

.010

G ,??????=11H ,初始状态??

?

???-=11)0(x ,1)(=k u ,

试用Z 变换法求解)(k x 。 解:

1

1

116.01)

(--??

?

???+-=-z z

G zI ??

??

?

?

??

?

???

++++-+++++=)8.0)(2.0()8.0)(2.0(16.0)8.0)(2.0(1

)8.0)(2.0(1z z z z z z z z z z

?????

?????????++

+-+++-+-+

++-+

+=8.0342.031

8.01542

.0154

8.0352.035

8

.031

2.034z z z z z z z z

而 )]()0([)()(1

z Hu zx G zI z x +-=- 1

)(-=

z z z u

????

?

???????

-+--=??????

?

???

--+??????-=+121

11)()0(2

2

z z z z z z

z z

z

z z z Hu zx ∴?????

???????-+++--+++=)1)(8.0)(2.0()84.1()

1)(8.0)(2.0()2()(2

2

z z z z z z z z z z z z x

???

?

??????-+

+-+-++++-=)1(187)8.0(96.17)2

.0(64.3)1(

1825)8.0(922)2.0(617

z z z z z z z z

z z z z

对)(z x 取z 反变换,有

??????????+

-+-+-+--=187)8.0(96.17)2.0(6

4.31825)8.0(922)2.0(617)(k k k

k k x

线性方程组有解的判别定理

非齐次线性方程组同解的讨论 摘要 本文主要讨论两个非齐次线性方程组有相同解的条件,即如何判定这两个非齐次线性方程组有相同的解. 关键词 非齐次线性方程组 同解 陪集 零空间 引言 无论是解齐次线性方程组,还是解非齐次线性方程组.所用的方法都是消元法,即对其系数矩阵或增广矩阵施以行的初等变换,而得到比较简单的同解方程组.用矩阵理论来说,就是系数矩阵或增广矩阵左乘以可逆矩阵后所得线性方程组与原线性方程组据有相同的解.这仅为问题的一面,而问题的反面是,如果两个非齐次线性方程组同解,则它们的系数矩阵或增广矩阵之间是否存在一个可逆矩阵?答案是肯定的,此即是本文主要解决的问题。 下面是一个非齐次线性方程组,我们用矩阵的形式写出 11121121222212n n m m mn m a x a x a x b a x a x a x b a x a x a x b +++=??+++=????+++=? 令 A= 111212122212n n m m mn a a a a a a a a a ???????????? ,b= 12m b b b ???????????? 。 即非齐次线性方程组可写成Ax b =。 一 、线性方程组同解的性质 引理 1 如果非齐次线性方程组Ax b =与Bx d =同解,则矩阵[]A b 与[]B d 的秩相等. 证明 设非齐次线性方程组Ax b =的导出组的基础解系为111,,,r ξξξ ,其中1 r 为矩阵[]A b 的秩,再设非齐次线性方程组Bx=d 的导出组的基础解系为 2 12,,,r ηηη ,其中2r 为矩阵[]B d 的秩,如果*η是非齐次线性方程组Ax=b 与Bx=d 特解,由于这两个方程组同解,所以向量组1*11,,,,r ξξξη 与向量组2*12,,,,r ηηηη 等价。从而这两个线性无关的向量组所含的向量个数相等,于是有12,r r =则矩阵[]A b 与[]B d 的秩相等. 引理[1]2 设A 、B 为m n ?矩阵,则齐次线性方程组0Ax =与0Bx =同解的充

控制系统状态方程求解

第2章 控制系统的状态方程求解 要点: ① 线性定常状态方程的解 ② 状态转移矩阵的求法 ③ 离散系统状态方程的解 难点: ① 状态转移矩阵的求法 ② 非齐次状态方程的解 一 线性定常系统状态方程的解 1 齐次状态方程的解 考虑n 阶线性定常齐次方程 ? ? ?==0)0()()(x x t Ax t x & (2-1) 的解。 先复习标量微分方程的解。设标量微分方程为 ? ??==0)0(x x ax x & (2-2) 对式(2-2)取拉氏变换得 )()(0s aX X s sX =- 移项 0)()(x s X a s =- 则 a s x s X -= )(

取拉氏反变换,得 00 0!)()(x k at x e t x k k at ∑∞ === 标量微分方程可以认为是矩阵微分方程当n=1时的特征,因此矩阵微分方程的解与标量微分方程应具有形式的不变性,由此得如下定理: 定理2-1 n 阶线性定常齐次状态方程(2-1)的解为 00 0!)()(x k At x e t x k k At ∑∞ === (2-3) 式中,∑∞ ==0 !)(k k At k At e 推论2-1 n 阶线性定常齐次状态方程 ???==00 )()()(x t x t Ax t x & (2-4) 的解为 0)(0 )(x e t x t t A -= (2-5) 齐次状态方程解的物理意义是)(0 t t A e -将系统从初始时刻0t 的初始 状态0x 转移到t 时刻的状态)(t x 。故)(0 t t A e -又称为定常系统的状态转移 矩阵。 (状态转移矩阵有四种求法:即定义(矩阵指数定义)法、拉氏反变换法、特征向量法和凯来-哈密顿(Cayly-Hamilton )法) 从上面得到两个等式 ∑∞ ==0 !)(k k At k At e ])[(11---=A sI L e At 其中,第一式为矩阵指数定义式,第二式可为At e 的频域求法或拉氏反变换法

列主元高斯消去法和列主元三角分解法解线性方程

计算方法实验报告1 【课题名称】 用列主元高斯消去法和列主元三角分解法解线性方程 【目的和意义】 高斯消去法是一个古老的求解线性方程组的方法,但由它改进得到的选主元的高斯消去法则是目前计算机上常用的解低阶稠密矩阵方程组的有效方法。 用高斯消去法解线性方程组的基本思想时用矩阵行的初等变换将系数矩阵A 约化为具有简单形式的矩阵(上三角矩阵、单位矩阵等),而三角形方程组则可以直接回带求解 用高斯消去法解线性方程组b Ax =(其中A ∈Rn ×n )的计算量为:乘除法运算步骤为 32(1)(1)(21)(1)(1)262233n n n n n n n n n n n MD n ----+= +++=+-,加减运算步骤为 (1)(21)(1)(1)(1)(25) 6226 n n n n n n n n n n AS -----+= ++= 。相比之下,传统的克莱姆 法则则较为繁琐,如求解20阶线性方程组,克莱姆法则大约要19 510?次乘法,而用高斯消去法只需要3060次乘除法。 在高斯消去法运算的过程中,如果出现abs(A(i,i))等于零或过小的情况,则会导致矩阵元素数量级严重增长和舍入误差的扩散,使得最后的计算结果不可靠,所以目前计算机上常用的解低阶稠密矩阵方程的快速有效的方法时列主元高斯消去法,从而使计算结果更加精确。 2、列主元三角分解法 高斯消去法的消去过程,实质上是将A 分解为两个三角矩阵的乘积A=LU ,并求解Ly=b 的过程。回带过程就是求解上三角方程组Ux=y 。所以在实际的运算中,矩阵L 和U 可以直接计算出,而不需要任何中间步骤,从而在计算过程中将高斯消去法的步骤进行了进一步的简略,大大提高了运算速度,这就是三角分解法 采用选主元的方式与列主元高斯消去法一样,也是为了避免除数过小,从而保证了计算的精确度 【计算公式】 1、 列主元高斯消去法 设有线性方程组Ax=b ,其中设A 为非奇异矩阵。方程组的增广矩阵为 第1步(k=1):首先在A 的第一列中选取绝对值最大的元素 1l a ,作为第一步的主元素: 111211212222112[,]n n n l n nn n a a a a b a a a b a a a b ?? ???? ?? =?????? ?? ????a b

三角分解法解线性方程组

三角分解法解线性方程组 #include #include #include //----------------------------------------------全局变量定义区 const int Number=15; //方程最大个数 double a[Number][Number],b[Number],copy_a[Number][Number],copy_b[Number]; // 系数行列式 int A_y[Number]; //a[][]中随着横坐标增加列坐标的排列顺序,如 a[0][0],a[1][2],a[2][1]...则A_y[]={0,2,1...}; int lenth,copy_lenth; //方程的个数 char * x; //未知量a,b,c的载体 int i,j; //----------------------------------------------函数声明区 void input(); //输入方程组 void print_menu(); //打印主菜单 int Doolittle_check(double a[][Number],double b[Number]); //判断是否 行列式>0,若是, 调整为顺序主子式全>0 void xiaoqu_u_l(); //将行列式Doolittle分解 void calculate_u_l(); // 计算Doolittle结果 void exchange(int m,int i); //交换A_y[m],A_y[i] void exchange_lie(int j); //交换a[][j]与b[]; void exchange_hang(int m,int n);

C++实现 牛顿迭代 解非线性方程组

C++实现牛顿迭代解非线性方程组(二元二次为例) 求解{0=x*x-2*x-y+0.5; 0=x*x+4*y*y-4; }的方程 #include #include #define N 2 // 非线性方程组中方程个数、未知量个数#define Epsilon 0.0001 // 差向量1范数的上限 #define Max 100 // 最大迭代次数 using namespace std; const int N2=2*N; int main() { void ff(float xx[N],float yy[N]); //计算向量函数的因变量向量yy[N] void ffjacobian(float xx[N],float yy[N][N]); //计算雅克比矩阵yy[N][N] void inv_jacobian(float yy[N][N],float inv[N][N]); //计算雅克比矩阵的逆矩阵inv void newdundiedai(float x0[N], float inv[N][N],float y0[N],float x1[N]); //由近似解向量x0 计算近似解向量x1 float x0[N]={2.0,0.25},y0[N],jacobian[N][N],invjacobian[N][N],x1[N],errornorm; int i,j,iter=0; //如果取消对x0的初始化,撤销下面两行的注释符,就可以由键盘x读入初始近似解向量for( i=0;i>x0[i]; cout<<"初始近似解向量:"<

第三章线性系统状态方程的解

第三章 系统的分析——状态方程的解 §3-1线性连续定常齐次方程求解 一、齐次方程和状态转移矩阵的定义 1、齐次方程 状态方程的齐次方程部分反映系统自由运动的状况(即没有输入作用的状况),设系统的状态方程的齐次部分为: )()(t Ax t x =& 线性定常连续系统: Ax x =& 初始条件:00x x t == 2、状态转移矩阵的定义 齐次状态方程Ax x =&有两种常见解法:(1)幂级数法;(2)拉氏变换法。其解为 )0()(x e t x At ?=。其中At e 称为状态转移矩阵(或矩阵指数函数、矩阵指数),记为: At e t =)(φ。 若初始条件为)(0t x ,则状态转移矩阵记为:) (0 0)(t t A e t t -=-Φ 对于线性时变系统,状态转移矩阵写为),(0t t φ,它是时刻t ,t 0的函数。但它一般不能写成指数形式。 (1)幂级数法——直接求解 设Ax x =&的解是t 的向量幂级数 Λ ΛΛΛ+++++=k k t b t b t b b t x 2210)( 式中ΛΛ,,, ,,k b b b b 210都是n 维向量,是待定系数。则当0=t 时, 000b x x t === 为了求其余各系数,将)(t x 求导,并代入)()(t Ax t x =&,得: Λ ΛΛΛ&+++++=-1232132)(k k t kb t b t b b t x )(2210ΛΛΛΛ+++++=k k t b t b t b b A

上式对于所有的t 都成立,故而有: ????? ??????======00 3 230 21201!1!31312121b A k b b A Ab b b A Ab b Ab b K K M 且有:00x b = 故以上系数完全确定,所以有: Λ ΛΛΛ+++++=k k t b t b t b b t x 2210)( ΛΛ++++ +=k k t b A k t b A t Ab b 020200! 1 !21 )0()! 1!21(22x t A k t A At I k k ΛΛ+++++= 定义(矩阵指数或矩阵函数): ∑∞==+++++=022! 1!1!21K k k k k At t A k t A k t A At I e ΛΛ 则 )0()(x e t x At ?=。 (2)拉氏变换解法 将Ax x =&两端取拉氏变换,有 )()0()(s AX X s sX =- )0()()(X s X A sI =- )0()()(1X A sI s X ?-=- 拉氏反变换,有 )0(])[()(1 1x A sI L t x ?-=--

线性方程组解的判定

1 / 3 第四节 线性方程组解的判定 从本节开始,讨论含有n 个未知量、m 个方程的线性方程组的解. 11112211211222221122n n n n m m mn n m a x a x a x b a x a x a x b a x a x a x b +++=??+++=????+++= ? (13—2) 主要问题是要判断出方程组(13-2)何时有解?何时无解?有解时解有多少?如何求出方程组的解。 线性方程组有没有解,以及有怎样的解,完全决定于方程组的系数和常数项。因此,将线性方程组写成矩阵形式或向量形式,以矩阵或向量作为讨论线性方程组的工具,将带来极大的方便。 方程组(13-2)中各未知量的系数组成的矩阵111212122212n n m m mn a a a a a a A a a a ??????=?????? 称为方程组(13-2)的系数矩阵.由各系数与常数项组成的矩阵,称为增广矩阵,记作A ,即 11121121 222212n n m m mn m a a a b a a a b A a a a b ??????=?????? 方程组(13-2)中的未知量组成一个n 行、1列的矩阵(或列向量),记作X ;常数项组成一个m 行、1列 的矩阵(或列向量),记作b ,即12n x x X x ??????=??????,12m b b b b ??????=?????? 由矩阵运算,方程组(13—2)实际上是如下关系111212122212 n n m m mn a a a a a a a a a ????????????12n x x x ????????????=12m b b b ???????????? 即 AX=b

系统的状态方程

第2章 系统的状态空间描述 输入输出:可测量,欠全面 §2.1 基本概念 例2.1 密封水箱 1 ()(),y t x t μ = 1 d [()()]d [()()]d c x u t y t t u t x t t μ ?=-?=-? 即 μ 2 (m ) c 3 ()(m /s)u t 3 ()(m /s)y t ()(m) x t

11 ()()()x t x t u t c c μ'=-+. 解 t t c c x t x u c 001()e ()e d τμμττ- ??=+ ? ??? ?. 若()u t r ≡, 则 0()e 1e ,()t t c c x t x r r t μμμμ--??=+-?→∞ ? ? ??, 若想()x h ∞=, 只要()h u t μ =.

例2.2 LRC 123()()();i t i t i t =+ ()()()()()L R L C u t v t v t v t v t =+=+ 选1()()C i t v t 和; 则: 1 1()()()1()()()C C C Li t v t u t Cv t i t v t R '=-+???'?=-? 其余 2()()/, C i t v t R = ()()(),()(). L C R C v t u t v t v t v t =-=)(t v C ) (t v L L R C )(1t i )(t u )(2t i )(3t i 2.2 图

1. 系统的状态变量 状态变量: 完全表征系统,个数最少的一组变量 未来()x t :由0()x t 和0t t ≥的()u t 完全确定. 对定常, 常取00t =. 2. 状态向量和状态空间 状态向量:12()(),(),()T n x t x t x t x t =???? 状态空间:()x t 取值范围 状态轨线:()x t 的轨迹(无时间轴) 3.几点说明

列主元高斯消去法和列主元三角分解法解线性方程

列主元高斯消去法和列主元三角分解法解线性方程

计算方法实验报告1 【课题名称】 用列主元高斯消去法和列主元三角分解法解线性方程 【目的和意义】 高斯消去法是一个古老的求解线性方程组的方法,但由它改进得到的选主元的高斯消去法则是目前计算机上常用的解低阶稠密矩阵方程组的有效方法。 用高斯消去法解线性方程组的基本思想时用矩阵行的初等变换将系数矩阵A 约化为具有简单形式的矩阵(上三角矩阵、单位矩阵等),而三角形方程组则可以直接回带求解 用高斯消去法解线性方程组b Ax =(其中A ∈Rn ×n )的计算量为:乘除法运算步骤为 32(1)(1)(21)(1)(1)262233 n n n n n n n n n n n MD n ----+= +++=+-,加减运算步骤为 (1)(21)(1)(1)(1)(25) 6226 n n n n n n n n n n AS -----+= ++= 。相比之 下,传统的克莱姆法则则较为繁琐,如求解20阶线性方程组,克莱姆法则大约要19 510?次乘法, 而用高斯消去法只需要3060次乘除法。

在高斯消去法运算的过程中,如果出现abs(A(i,i))等于零或过小的情况,则会导致矩阵元素数量级严重增长和舍入误差的扩散,使得最后的计算结果不可靠,所以目前计算机上常用的解低阶稠密矩阵方程的快速有效的方法时列主元高斯消去法,从而使计算结果更加精确。 2、列主元三角分解法 高斯消去法的消去过程,实质上是将A 分解为两个三角矩阵的乘积A=LU ,并求解Ly=b 的过程。回带过程就是求解上三角方程组Ux=y 。所以在实际的运算中,矩阵L 和U 可以直接计算出,而不需要任何中间步骤,从而在计算过程中将高斯消去法的步骤进行了进一步的简略,大大提高了运算速度,这就是三角分解法 采用选主元的方式与列主元高斯消去法一样,也是为了避免除数过小,从而保证了计算的精确度 【计算公式】 1、 列主元高斯消去法 设有线性方程组Ax=b ,其中设A 为非奇异矩阵。方程组的增广矩阵为 1112112122221[,]n n l a a a a b a a a b ?????? ??=? ? ?? ??? ?a b L L M L L M M M M M M M M M

基于Matlab的牛顿迭代法解非线性方程组

基于Matlab 实现牛顿迭代法解非线性方程组 已知非线性方程组如下 2211221212 10801080x x x x x x x ?-++=??+-+=?? 给定初值0(0,0)T x =,要求求解精度达到0.00001 首先建立函数F(x),方程组编程如下,将F.m 保存到工作路径中: function f=F(x) f(1)=x(1)^2-10*x(1)+x(2)^2+8; f(2)=x(1)*x(2)^2+x(1)-10*x(2)+8; f=[f(1) f(2)]; 建立函数DF(x),用于求方程组的Jacobi 矩阵,将DF.m 保存到工作路径中: function df=DF(x) df=[2*x(1)-10,2*x(2);x(2)^2+1,2*x(1)*x(2)-10]; 编程牛顿迭代法解非线性方程组,将newton.m 保存到工作路径中: clear; clc x=[0,0]'; f=F(x); df=DF(x); fprintf('%d %.7f %.7f\n',0,x(1),x(2)); N=4; for i=1:N y=df\f'; x=x-y; f=F(x); df=DF(x); fprintf('%d %.7f %.7f\n',i,x(1),x(2)); if norm(y)<0.0000001 break ; else end end

运行结果如下: 0 0.0000000 0.0000000 1 0.8000000 0.8800000 2 0.9917872 0.9917117 3 0.9999752 0.9999685 4 1.0000000 1.0000000

求解系统的状态方程

求解系统的状态方程 一、实验设备 PC计算机,MATLAB软件,控制理论实验台 二、实验目的 (1)掌握状态转移矩阵的概念。学会用MATLAB求解状态转移矩阵 (2)学习系统齐次、非齐次状态方程求解的方法,计算矩阵指数,求状态响应; (3)通过编程、上机调试,掌握求解系统状态方程的方法,学会绘制输出响应和状态响应曲线; (4)掌握利用MATLAB导出连续状态空间模型的离散化模型的方法。 三、实验原理及相关基础 (1)参考教材P99~101“3.8利用MATLAB求解系统的状态方程” (2)MATLAB现代控制理论仿真实验基础 (3)控制理论实验台使用指导 四、实验内容 (1)求下列系统矩阵A对应的状态转移矩阵 (a)

(b) 代码: syms lambda A=[lambda 0 0;0 lambda 0;0 0 lambda];syms t;f=expm(A*t) (c) 代码: syms t;syms lambda;A=[lambda 0 0 0;0 lambda 1 0;0 0 lambda 1;0 0 0 lambda];f=expm(A*t) (2) 已知系统

a) 用MATLAB求状态方程的解析解。选择时间向量t,绘制系统的状态响应曲线。观察并记录这些曲线。 (1) 代码: A=[0 1; -2 -3]; B=[3;0]; C=[1 1]; D=[0]; u=1; syms t; f=expm(A*t);%状态转移矩阵 x0=0; s1=f*B*u; s2=int(s1,t,0,t)%状态方程解析解 状态曲线: (2)A=[0 1;-2 -3]; syms t; f=expm(A*t); X0=[1;0]; t=[0:0.5:10]; for i=1:length(t); g(i)=double(subs(f(1),t(i))); end plot(t,g)

第四章线性方程组直接法,矩阵三角分解

第四章 习题答案 1。用Gauss 消去法解方程组 1231231 2323463525433032 x x x x x x x x x ++=?? ++=??++=? 解:方程组写成矩阵形式为12323463525433032x x x ?????? ? ? ? = ? ? ? ? ? ????? ?? 对其进行Gauss 消去得12323441 4726002x x x ?? ???? ? ? ? ?-= ? ? ? ? ? ?????-?? 得方程组12312323 32346 131 44 822 24 x x x x x x x x x ++=?=-???? -=-?=????=?-=-?? 2。用Gauss 列主元素消去法解方程组 1233264107075156x x x -?????? ? ? ?-= ? ? ? ? ? ?-???? ?? 解:因为第一列中10最大,因此把10作为列主元素 1233264107075156x x x -?????? ? ? ?-= ? ? ? ? ? ?-??????12r r ????→1231070732645156x x x -?????? ? ? ?-= ? ? ? ? ? ?-???? ?? 21 3113 10122 31070716106101055052 2r r r r x x x +-? ??? ? ?-?? ? ? ? ? ????→-= ? ? ? ? ? ??? ? ?????23 r r ????→123107075505221 61061010x x x ? ??? ? ?-?? ? ? ? ? ?= ? ? ? ? ? ? ?? ? ?-????

线性方程组解的判定

第四节 线性方程组解的判定 从本节开始,讨论含有n 个未知量、m 个方程的线性方程组的解。 11112211211222 22 11 22n n n n m m mn n m a x a x a x b a x a x a x b a x a x a x b +++=??+ ++= ????+++=? (13—2) 主要问题是要判断出方程组(13-2)何时有解?何时无解?有解时解有多少?如何求出方程组的解。 线性方程组有没有解,以及有怎样的解,完全决定于方程组的系数和常数项。因此,将线性方程组写成矩阵形式或向量形式,以矩阵或向量作为讨论线性方程组的工具,将带来极大的方便。 方程组(13-2)中各未知量的系数组成的矩阵11121212221 2 n n m m mn a a a a a a A a a a ? ?? ? ? ?=?? ?? ? ? 称为方程组(13-2)的系数矩阵。由各系数与常数项组成的矩阵,称为增广矩阵,记作A ,即 11121121 222212 n n m m mn m a a a b a a a b A a a a b ?? ????=??? ??? 方程组(13-2)中的未知量组成一个n 行、1列的矩阵(或列向量),记作X;常数项组成一个m 行、1 列的矩阵(或列向量),记作b ,即12n x x X x ??????=?????? ,12 m b b b b ?? ????=?????? 由矩阵运算,方程组(13-2)实际上是如下关系111212122212 n n m m mn a a a a a a a a a ? ?? ? ? ? ?? ?? ? ? 12n x x x ???????????? =12m b b b ???????????? 即 AX=b

非线性方程的简单迭代法和Steffensen迭代法

《数值计算方法》实验报告 实验名称:实验1 非线性方程的简单迭代法和Steffensen 迭代法 实验题目:分别用简单迭代法和Steffensen 迭代法求方程 010423=-+x x 在 [1, 2] 内的一个实根. 实验目的:理解并掌握简单迭代法和Steffensen 迭代法 基础理论:简单迭代法和Steffensen 迭代法 1).简单迭代法的原理:将一元非线性方程:0)(=x f 改写成等价方程:)(x x ρ= ,对此,从某个初始值x0开始,对应式)(x x ρ= 构成迭代公式 ,...1,0),(1==+k x x k k ρ ,这样就可以确定序列 {}k x (k=0,1,2…)。如果 {}k x 有极限 *lim x x k k =∞→ ,由式 ,...1,0),(1==+k x x k k ρ 两边取极限可得 )(**x x ρ= ,可知 * x 为方程0)(=x f 的近似解。 2)Steffensen 迭代法的原理: 通过把改进的Aitken 方法应用于根据不动点迭代所得到的线性收敛序列,将收敛速度加速到二阶。

()???? ?????+---===+k k k k k k k k k k k x y z x y x x y z x y 2) ()(21ρρ []x x x x x x x +---=)(2)(()()(2ρρρρψ 实验环境:操作系统:Windows 7; 实验平台:Turbo C++ 实验过程:写出算法→编写程序→调试运行程序→计算结果 1)简单迭代法的算法: Input:初始近似值x0,精度要求del,最大迭代次数N Output:近似解x 或失败信息 1. n ←1 2. While n ≤N do; 3. x ←f(x0); 4. if | x-x0|

不动点迭代法求解非线性方程组

不动点迭代法求解非线性方程组 摘要:一般非线性方程组可以写成()0F x =的形式,其中:n m F R R →是定义在区域 n D R ?上的向量函数。把方程组()0F x =改写成与之等价的形式:(x G x =)。因此,求 方程组()0F x =的解就转化为求函数的(G x )的不动点。本文首先介绍了多变量函数()F x 的微积分性质,接着介绍了用不动点迭代法求解非线性方程组。 关键词:多变量函数;微积分;不动点 Fixed Point Iteration Method For Solving Nonlinear Equations Abstract :General nonlinear equations can be written in the form of ()F x θ=, where the vector function :n m F R R →is defined on the region n D R ?. Transform the equations ()0F x = into its equivalent form: (x G x =).Therefore, we can get the solution of ()0F x = by finding the fixed point of (G x ).In this paper, we first introduce some knowledge about multivariable calculus, then introduce the fixed point iteration method for solving nonlinear equations. Key words: multi-variable function; calculus; fixed point 1 引言 一般非线性方程组及其向量表示法:含有n 个方程的n 元非线性方程组的一般形式为 11221212(,...,)0(,...,)0......(,...,)0 n n m n f x x x f x x x f x x x =??=?? ? ?=? (1) 其中,()1,2,...,i f i n =是定义在n D R ?上的n 元实值函数,且i f 中至少有一个是非线性 函数。

线性方程组解的判定与解的结构

***学院数学分析课程论文 线性方程组解的判定与解的结构 院系数学与统计学院 专业数学与应用数学(师范) 姓名******* 年级 2009级 学号200906034*** 指导教师 ** 2011年6月

线性方程组解的判定与解的结构 姓名****** (重庆三峡学院数学与计算机科学学院09级数本?班) 摘 要:线性方程组是否有解,用系数矩阵和增广矩阵的秩来刻画.在方程组有解且有 多个解的情况下,解的结构就是了解解与解之间的关系. 关键词:矩阵; 秩; 线性方程组; 解 引言 通过系数矩阵和增广矩阵的秩是否相同来给出判定线性方程组的解的判别条件.在了解了线性方程组的判别条件之后,我们进一步讨论解的结构.对于齐次线性方程组,解的线性组合还是方程组的解.在线性方程组有无穷个解时可用有限多个解表示出来.另外以下还涉及到线性方程组通解的表达方式. 1 基本性质 下面我们分析一个线性方程组的问题,导出线性方程组有解的判别条件. 对于线性方程组 1111221121122222 1122n n n n s s sn n s a x a x a x b a x a x a x b a x a x a x b ++???+=??++???+=???????++???+=? (1) 引入向量 112111s αααα??????=?????????,122222s αααα??????=?????????,…12n n n sn αααα??????=????????? ,12s b b b β?? ?? ??=??????? ?? 方程(1)可以表示为 1122n n x x x αααβ++???+= 性质 线性方程组⑴有解的充分必要条件为向量β可以表成向量组α1,α2,…,αn 的线性组合. 定理1 线性方程组⑴有解的充分必要条件为它的系数矩阵

数值分析求解非线性方程根的二分法、简单迭代法和牛顿迭代法

实验报告一:实验题目 一、 实验目的 掌握求解非线性方程根的二分法、简单迭代法和牛顿迭代法,并通过数值实验比较两种方法的收敛速度。 二、 实验内容 1、编写二分法、牛顿迭代法程序,并使用这两个程序计算02)(=-+=x e x x f 在[0, 1]区间的解,要求误差小于 4 10- ,比较两种方法收敛速度。 2、在利率问题中,若贷款额为20万元,月还款额为2160元,还期为10年,则年利率为多少?请使用牛顿迭代法求解。 3、由中子迁移理论,燃料棒的临界长度为下面方程的根,用牛顿迭 代法求这个方程的最小正根。 4、用牛顿法求方程 的根,精确至8位有效数字。比较 牛顿迭代法算单根和重根的收敛速度,并用改进的牛顿迭代法计算重根。 三、 实验程序 第1题: 02)(=-+=x e x x f 区间[0,1] 函数画图可得函数零点约为0.5。 画图函数: f un cti on Te st1() % f(x ) 示意图, f(x) = x + exp (x) - 2; f(x) = 0 r = 0:0.01:1; y = r + e xp(r) - 2 p lot(r, y); gri d on 二分法程序: 计算调用函数:[c,n um ]=bisec t(0,1,1e-4) fu ncti on [c,num ]=bisect (a,b,de lt a) %Inp ut –a,b 是取值区间范围 % -de lta 是允许误差 %O utput -c牛顿迭代法最后计算所得零点值 % -num 是迭代次数 ya = a + exp(a) - 2; yb = b + e xp(b) - 2;

牛顿迭代法求解非线性方程组的代码

牛顿迭代法求解非线性方程组 非线性方程组如下: 221122121210801080 x x x x x x x ?-++=??+-+=?? 给定初值()00.0T x =,要求求解精度达到0.00001 1.首先建立函数()F X ,方程编程如下,将F.m 保存到工作路径中: function f=F(x) f(1)=x(1)^2-10*x(1)+x(2)^2+8; f(2)=x(1)*x(2)^2+x(1)-10*x(2)+8; f=[f(1),f(2)] ; 2.建立函数()DF X ,用于求方程的jacobi 矩阵,将DF.m 保存到工作路径中: function df=DF(x) df=[2*x(1)-10,2*x(2);x(2)^2+1,2*x(1)*x(2)-10]; %jacobi 矩阵是一阶偏导数以一定方式排列成的矩阵。 3.编程牛顿迭代法解非线性方程组,将newton.m 保存在工作路径中: clear,clc; x=[0,0]'; f=F(x);

df=DF(x); fprintf('%d %.7f %.7f\n',0,x(1),x(2)); N=4; for i=1:N y=df\f'; x=x-y; f=F(x); df=DF(x); fprintf('%d %.7f %.7f\n',i,x(1),x(2)); if norm(y)<0.0000001 break; else end end ezplot('x^2-10*x+y^2+8',[-6,6,-6,6]); hold on ezplot('x*y^2+x-10*y+8',[-6,6,-6,6]); 运行结果如下: 0 0.0000000 0.0000000 1 0.8000000 0.8800000 2 0.9917872 0.9917117

线性方程组的公共解

线性方程组的公共解 问题:如何求解线性方程组的公共解? 线性方程组是高代学习的一个重点内容,它的一般形式为 ???????=+++=+++=+++bs asnxn x as x as b nxn a x a x a b nxn a x a x a ...2211... ,22...222121,11...212111 而线性方程组的求解也是这部分学习的重点和难点。其中求解线性方程组的公共解也是高等代数学习所必须掌握的一个知识点。 例1、证明:对于n 元齐次线性方程组(Ⅰ)AX=0与(Ⅱ)BX=0,有非零公共解的充要条件是r(B A )

???=-=+0 42031x x x x 又已知某齐次线性方程组(Ⅱ)的通解为 k1(0,1,1,0)’+k2(-1,2,2,1)’ 问(Ⅰ)与(Ⅱ)是否有非零公共解?若有,则求出所有公共解,若没有,则说明理由。(出自2005年中科院) 解:方法一:将(Ⅱ)的通解代入方程组(Ⅰ)得 ???=+=+0 21021k k k k 解得k1=-k2,故方程组(Ⅰ)与(Ⅱ)有非零公共解,所有非零公共解为k (1,1,1,1)’,k ≠0为任意常数 方法二:令方程组(Ⅰ)与(Ⅱ)的通解相同,即 k1(0,1,1,0)’+k2(-1,2,2,1)’=k3(-1,0,1,0)’+k4(0,1,0,1)’ 得到关于k1,k2,k3,k4的一个方程组 ???????=-=-+=-+=-0 420 422103221032k k k k k k k k k k 可求其通解为(k1,k2,k3,k4)’=k(-1,1,1,1)’ 将k1=-1,k2=k 代入(Ⅰ)的通解可得所有非零公共解为k (1,1,1,1)’,k ≠0为任意常数 方法三:方程组(Ⅱ)可以是 ? ??=+=+-041032x x x x 解(Ⅰ)与(Ⅱ)的联立方程组可得所有非零公共解为k (1,1,1,1)’,k ≠0为任意常数 韩梦雪 20132113429

(完整word版)c++求解非线性方程组的牛顿顿迭代法

牛顿迭代法c++程序设计 求解{0=x*x-2*x-y+0.5; 0=x*x+4*y*y-4; }的方程 #include #include #define N 2 // 非线性方程组中方程个数、未知量个数 #define Epsilon 0.0001 // 差向量1范数的上限 #define Max 100 //最大迭代次数 using namespace std; const int N2=2*N; int main() { void ff(float xx[N],float yy[N]); //计算向量函数的因变量向量yy[N] void ffjacobian(float xx[N],float yy[N][N]);/ /计算雅克比矩阵yy[N][N] void inv_jacobian(float yy[N][N],float inv[N][N]); //计算雅克比矩阵的逆矩阵inv void newdundiedai(float x0[N], float inv[N][N],float y0[N],float x1[N]); //由近似解向量x0 计算近似解向量x1 float x0[N]={2.0,0.25},y0[N],jacobian[N][N],invjacobian[N][N],x1[N],errornorm; int i,j,iter=0; //如果取消对x0的初始化,撤销下面两行的注释符, 就可以由键盘向x0读入初始近似解向量for( i=0;i>x0[i]; cout<<"初始近似解向量:"<

第三章线性系统状态方程的解

第三章 线性系统的运动分析 §3-1线性连续定常齐次方程求解 一、齐次方程和状态转移矩阵的定义 1、齐次方程 状态方程的齐次方程部分反映系统自由运动的状况(即没有输入作用的状况),设系统的状态方程的齐次部分为:)()(t Ax t x = 线性定常连续系统:Ax x = 2、状态转移矩阵的定义 齐次状态方程Ax x = 有两种常见解法:(1)幂级数法;(2)拉氏变换法。其解为)0()(x e t x At ?=。 其中At e 称为状态转移矩阵(或矩阵指数函数、矩阵指数),记为:At e t =)(φ。 若初始条件为)(0t x ,则状态转移矩阵记为:)(00 )(t t A e t t -=-Φ 对于线性时变系统,状态转移矩阵写为),(0t t φ,它是时刻t ,t 0的函数。但它一般不能写成指数形式。 (1)幂级数法 设Ax x = 的解是t 的向量幂级数 +++++=k k t b t b t b b t x 2210)( 式中 ,,, ,,k b b b b 210都是n 维向量,则 +++++=-1232132)(k k t kb t b t b b t x )(2210 +++++=k k t b t b t b b A 故而有: ????? ?? ????== ====003 230 2 12 01!1! 3131 2 121b A k b b A Ab b b A Ab b Ab b K K

且有0)0(b x =。 故 +++++=k k t b t b t b b t x 2210)( ++ +++=k k t b A k t b A t Ab b 02 02 00! 1! 21 )0()! 1!21(22 x t A k t A At I k k ++ ++ += 定义:∑ ∞ == ++ +++=0 2 2! 1! 1!21K k k k k At t A k t A k t A At I e 则)0()(x e t x At ?=。 (2)拉氏变换解法 将Ax x = 两端取拉氏变换,有 )()0()(s Ax x s sx =- )0()()(x s x A sI =- )0()()(1x A sI s x ?-=- 拉氏反变换,有 )0(])[()(11x A sI L t x ?-=-- 则 ])[()(11---==A sI L e t At φ 【例3.1.1】 已知系统的状态方程为x x ?? ? ???=00 10 ,初始条件为)0(x ,试求状态转移矩阵和状态方程的解。 解:(1)求状态转移矩阵 ++ ++ +==k k At t A k t A At I e t ! 1! 21)(2 2φ 此题中: ???? ??=00 10A , ?? ? ???====00 0032n A A A 所以

相关文档
最新文档