1 / 61

06 Matlab 符号计算

06 Matlab 符号计算. 一、符号计算基础 二、符号计算在微积分中的应用 三、符号表达式的化简与替换 四、求解符号方程（组） 五、积分变换. 一、符号计算基础. 1 、数值计算的误差 2 、 Matlab 符号计算功能简介 3 、创建符号量 sym / syms 4 、替换符号量 subs 5 、符号计算的运算 符 6 、符号计算的基本函数 7 、 符号计算的精度控制 8 、符号量与数值量之间的转换. 1 、数值运算的误差. 2 、 Matlab 符号计算功能简介. 符号计算以推理解析的方式进行，不受计算误差积累问题困扰；

hakan
Download Presentation

06 Matlab 符号计算

An Image/Link below is provided (as is) to download presentation Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author. Content is provided to you AS IS for your information and personal use only. Download presentation by click this link. While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server. During download, if you can't get a presentation, the file might be deleted by the publisher.

E N D

Presentation Transcript

  1. 06 Matlab 符号计算 一、符号计算基础 二、符号计算在微积分中的应用 三、符号表达式的化简与替换 四、求解符号方程（组） 五、积分变换

  2. 一、符号计算基础 1、数值计算的误差 2、Matlab符号计算功能简介 3、创建符号量 sym / syms 4、替换符号量 subs 5、符号计算的运算符 6、符号计算的基本函数 7、符号计算的精度控制 8、符号量与数值量之间的转换

  3. 1、数值运算的误差

  4. 2、Matlab符号计算功能简介 • 符号计算以推理解析的方式进行，不受计算误差积累问题困扰； • 符号计算给出解析解，当解析解不存在时，会给出数值解； • 符号计算指令的调用比较简单，与教科书上的公式相近； • 计算所需时间比数值计算要多很多； • Matlab 的符号计算功能通过Symbolic Math Toolbox（符号数学工具箱）和Extended Symbolic Math Toolbox（扩展符号数学工具箱）实现。 • Matlab 2008b以前的版本，符号计算引擎为 Maple • Matlab 2008b以后的版本（包括Maltab 2008b)，符号计算引擎为MuPAD • 常用的符号计算软件还有Maple、Mathmatica、MathCAD等。

  5. 3、创建符号量 (1)方法1：sym( ) 每次只能创建 一个 符号量。 S = sym(A) x = sym('x') x = sym('x','real') x为实数 x = sym('x','positive') x为正实数 x = sym('x','unreal') 去除x的“实”、“正”属性 S = sym(A,flag) flag为'r', 'd', 'e', 或 'f' a=sym('a') a是符号变量 b是符号常量 b=sym(1/3)

  6. b = sym('1/3+sqrt(2)') b = 1/3+sqrt(2) 用单引号括起来，是以最精确的形式存储符号常数 a = sym(1/3 + sqrt(2)) a = 7870251548315938*2^(-52) 不用单引号，则是以最接近的“有理”表示的形式存储符号常数

  7. (2) 方法2：syms 可以同时建立 多个 符号变量 syms arg1 arg2 ... syms arg1 arg2 ... real syms arg1 arg2 ... unreal syms arg1 arg2 ... positive a=sym('a'); b=sym('b'); c=sym('c'); syms a b c a b c中间一定要用空格隔开，不能用逗号 syms不能用来创建符号常量

  8. (3)创建符号表达式 • 用 sym 命令创建符号表达式： >> f = sym('a*x+b') f = a*x+b • 先定义符号变量，再按普通书写形式创建符号表达式： >> syms a x b >> f = a*x + b f = a*x+b 建议用这种方式！

  9. (4)创建符号数组 • 将数值矩阵转化成符号数组： B=[2/3,sqrt(2);5.2,log(3)]; C=sym(B) • 使用sym函数直接生成： A=sym('[1+x,sin(x); 5,exp(x)]') f=sym('[1,ab;c,d]') • 先定义符号变量，再按普通书写 形式创建符号数组: syms a b x y C= [a,b;x,y]

  10. (5) findsym( ) 找出符号表达式/数组中的符号变量 findsym(S) 按字母顺序列出 S中的所有符号变量 findsym(S,n)列出 S 中离 x最近的 n个符号变量 若表达式中有两个符号变量与 x的距离相等， 则ASCII 码大者优先。 >> f=sym('2*w-3*y+z^2+5*a') f = 2*w-3*y+z^2+5*a >> findsym(f) ans = a, w, y, z >> findsym(f,3) ans = y,w,z >> findsym(f,1) ans = y

  11. 4、替换符号量 subs( ) 用给定的数据替换符号表达式中的指定的符号变量 R = subs(S)替换 S中的所有变量 R = subs(S, new)用 new 替换 S 中的缺省变量 R = subs(S,old,new)用 new替换 S中的 old 若 new是一个由多个符号变量组成的数组， 则 old应该具有与 new相同的形状的数组。 syms a x b f = a*x+b a = 1 b = 2 c = subs(f) syms a x b f = a*x+b c = subs(f,5) syms a x b f = a*x+b c = subs(f,[a,b],[1,2]) c = 5*a+b c = x+2 缺省变量: findsym(S,1)的返回结果 c =x+2

  12. syms a x b f = a*x+b c = subs(f,[1,2;3,4]) syms a x b f = a*x+b c = subs(f,1:5) c = [ a+b, 2*a+b] [ 3*a+b, 4*a+b] c = [ a+b, 2*a+b, 3*a+b, 4*a+b, 5*a+b] syms a x b f = a*x+b c = subs(f,{a,b},{[1,2],[3,4]}) syms a x b f = a*x+b c = subs(f,{a,b},{[1,2;3,4],5}) c = [ x+3, 2*x+4] c = [ x+5, 2*x+5] [ 3*x+5, 4*x+5]

  13. 练习：指出下面各条语句的输出结果 >>f=sym('2*u'); >>subs(f,'u',2) >>f2=subs(f,'u','u+2') >>a=3; >>subs(f2,'u',a+2) >>subs(f2,'u','a+2') >>syms x y >>f3=subs(f,'u',x+y) >>subs(f3,[x,y],[1,2]) >>subs(f3,[x,y],[x+y,x+y])

  14. 5、符号计算的运算符 Matlab 符号计算中采用的运算符，在名称和使用上，都与数值计算中的运算符完全相同。 • 矩阵运算：+ - * \ / ^ • 数组运算：+ - .* .\ ./ .^ • 转置运算：' .' X=sym('[a,b;c,d]') Y=sym('[a,1;b,0]') Z1=X*Y Z2=X.'.*Y Z2 = [ a^2, c] [ b^2, 0] Z1 = [ a^2+b^2, a] [ c*a+d*b, c]

  15. 6、符号计算的基本函数 Matlab 符号计算中，一些基本的函数，在名称和使用上，与数值计算中的函数基本相同。 三角函数，双曲函数及它们的反函数，无论在数值计算还是在符号计算中，使用方法相同。但是没有 atan2函数。 指数、对数函数：与数值计算一样 复数函数：符号计算中，没有angle函数（求相角） 矩阵代数函数：与数值计算一样 与数组相关的函数，如size，length，括号[ ]、( )的作用等等，都与数值计算一样

  16. 练习：用det( )函数计算下列符号行列式：

  17. 7、符号计算的精度控制 digits(n) 设置符号数值计算以n位相对精度进行，默认为32位 xs = vpa(x) 在digits指定的精度下，给出x的数值型符号结果xs xs = vpa(x,n) 在n位相对精度下，给出x的数值型符号结果xs a = sym(pi) digits(20) b = vpa(a*a) digits(32) c = vpa(a*a) d = vpa(a*a,40) a = pi b = 9.8696044010893586191 c = 9.8696044010893586188344909998761 d = 9.869604401089358618834490999876151135313

  18. 8、符号量与数值量之间的转换 (1) double( ) 把符号常数转化为16位相对精度的浮点数值对象。 a = sym('1/123') + sym(sqrt(2)) + sym('pi') b = double(a) a = 1/123+2^(1/2)+pi b = 4.5639 (2) s = poly2sym(p) 把多项式转化相应的符号表达式 p = [1,2,3] s = poly2sym(p,x) s = x^2+2*x+3

  19. (3) p = sym2poly(s) 把多项式的符号表达式转换为数值多项式 syms x s = x + 2*x^3 - 6 p = sym2poly(s1) p = 2 0 1 -6

  20. 二、符号计算在微积分中的应用 • 1、极限 • 2、微分 • 3、积分 • 4、级数求和 • 5、函数的泰勒级数展开

  21. 1、极限 limit(F) limit(F,a) limit(F,v,a) limit(F,v,a,’right’) limit(F,v,a,’left’) 求默认变量→0时F的极限 求默认变量→a时F的极限 求v →a时F的极限 求v →a+时F的右极限 求v →a-时F的左极限 F为符号表达式，没有指定变量时，由 findsym确定默认变量 syms x a h F = (sin(x+h) - sin(x))/h f = (1 + a/x)^x lim_1 = limit(F,h,0) lim_2 = limit(f,inf)

  22. syms x limit(x/abs(x), x, 0, 'left') ans = -1 syms x limit(x/abs(x), x, 0, 'right') ans = 1 syms x limit(x/abs(x), x, 0) ans = NaN

  23. 2、导数（偏导数） diff(S) diff(S,n) diff(S,v) diff(S,v,n) 求 S对默认变量的一阶导数 求 S对默认变量的 n阶导数 求 S对变量 v的一阶导数 求 S对变量 v的 n阶导数，也可写成 diff(f,n,v) syms x y a b f = sqrt(exp(x)+x*sin(x)) df = diff(f) z = exp(a*x)*cos(b*y) yy = diff(z,y,2) xy = diff(diff(z,x),y)

  24. 3、积分 int(S) int(S,v) int(S,a,b) int(S,v,a,b) 求 S对默认变量的不定积分 求 S对变量 v的不定积分 求 S对默认变量的定积分 求 S对变量 v的定积分 a：积分下限，b：积分上限。 当 a或 b取 inf（或 -inf）时，计算的就是广义积分。 syms x y A1 = int(sin(x*y),y) A2 = int(sin(x)+2,0,pi/6) A3 = int(1/x^2,1,inf) A4 = int(1/(x^2-2*x+3),-inf,inf)

  25. 练习：求下面的积分，给出50位精度的数值 先对 z 变量的积分 然后对 y 变量的积分 最后对 x 变量的积分

  26. 4、级数求和 symsum(s,v) 求s的前v项的和，变量为 v = 0,1,2,…,v-1 功能同上，对默认变量求和 symsum(s) symsum(s,v,a,b) 变量v = a,a+1,…,b，求这些项的和 symsum(s,a,b) 默认变量 = a,a+1,…,b，求这些项的和 syms k a = symsum(k,1,100) b = symsum(1/k^2,1,inf)

  27. syms q n s1 = symsum(q^n,n,0,n-1) s2 = symsum(q^n,n) s1 = q^n/(q-1)-1/(q-1) s2 = q^n/(q-1) √√√ 按照symsum函数的说明，这里s2应该与s1相同，但结果却并不一致，故建议使用symsum函数时使用下面的完整的形式： symsum(s,v,a,b) ？？？

  28. 5、函数的泰勒级数展开 taylor(f) taylor(f,n) taylor(f,a) taylor(f,n,a) taylor(f,n,v) taylor(f,n,v,a) 求 f在默认变量 = 0处的 5阶泰勒展开式 求 f在默认变量 = 0处的 n-1 阶泰勒展开式 求 f在默认变量 = a处的 5 阶泰勒展开式 求 f在默认变量 = a处n-1 阶泰勒展开式 求 f在变量 v = 0处的 n-1 阶泰勒展开式 求 f在变量 v = a 处的 n-1阶泰勒展开式 当 a 为正整数时，请指定 n 的值。例如： a = 4时， taylor(f,4) : 返回的是f 在 0 处的 3阶泰勒展开式 taylor(f,7,4):返回的是f 在 4 处的 6阶泰勒展开式

  29. syms x f = sin(x)/x t1 = taylor(f,8) t2 = taylor(exp(x), 8) t3 = taylor(sin(x), 8) t4 = taylor(cos(x), 8)

  30. 三、符号表达式的化简 1、factor 2、expand 3、collect 4、simplify 5、simple 6、horner 7、numden

  31. 1、factor factor(f)对 f进行因式分解，也可用于正整数的分解 ans = (x-y)*(x^2+x*y+y^2) syms x y factor(x^3 - y^3) a = 1234567890 b1 = factor(a) b2 = factor(sym(a)) 输入数值量，返回一个行数组 输入符号量，返回符号表达式 b1 = 2 3 3 5 3607 3803 b2 = (2)*(3)^2*(5)*(3803)*(3607)

  32. syms a b c d r = [a,b,c,d] A = [r.^0;r;r.^2;r.^4] B = det(A) C = factor(B)

  33. 2、expand expand(S) 把 S展开为多项式形式，也可展开复合的三角函数、指数函数 syms x expand((x-2)*(x-4)) ans = x^2-6*x+8 syms x y; expand(cos(x+y)) ans = cos(x)*cos(y)-sin(x)*sin(y) syms a b expand(exp((a + b)^2)) ans = exp(a^2)*exp(a*b)^2*exp(b^2)

  34. 3、collect collect(S) 按默认变量对S合并同类项 collect(S,v)按指定变量v对 S合并同类项 syms x y; f = (x+y)*(x^2+y^2+1) f1 = collect(f) f2 = collect(f,y) f1 = x^3+y*x^2+(y^2+1)*x+y*(y^2+1) f2 = y^3+x*y^2+(x^2+1)*y+x*(x^2+1)

  35. 4、simplify simplify(S):用内置的化简规则对S进行化简 syms x; f = sin(x)^2 + cos(x)^2 ; f1 = simplify(f) f1 = 1 syms a b c; g = exp(c*log(sqrt(a+b))); g1 = simplify(g) g1 = (a+b)^(1/2*c)

  36. 5、simple simple(S)对S尝试多种不同的算法（包括simplify）进行简化，并显示相应的结果，最后返回其中最简短的形式并赋值给默认变量ans r = simple(S)对S尝试多种不同的算法（包括simplify）进行简化，但不显示相应的结果，只是返回其中最简短的形式并赋值给变量 r [r, How] = simple(S)对S尝试多种不同的算法（包括simplify）进行简化，但不显示相应的结果，只是返回其中最简短的形式并赋值给变量 r，同时返回字符串How ，记录得到 r所使用的化简方法。

  37. syms x; f = (1/x^3+6/x^2+12/x+8)^(1/3); f1 = simple(f) f2 = simple(f1) f3 = simple(f2) f1 = (2*x+1)/x f2 = 2+1/x f3 = 2+1/x 多次使用 simple ,直到结果与前一次相同，以达到最简单的表示形式。

  38. 6、horner horner(P)将多项式 P改写为horner形式（这种形式可以快速计算多项式的值） syms x; f = x^4+2*x^3+4*x^2+x+1; g = horner(f) g = 1+(1+(4+(x+2)*x)*x)*x

  39. 7、numden [N,D]=numden(f) 返回值： N 为 f 的分子，D 为 f 的分母 syms x y; f = x/y + y/x; [N,D]=numden(f) N = x^2+y^2 D = y*x

  40. 四、求解符号方程 • 1、线性方程组的符号解 • 2、符号方程(组)求解 • 3、常微分方程（组）的求解

  41. 1、线性方程组的符号解 clc; clear all; A = sym('[1,1,1;3,-1,6;0,1,3]') b = sym('[a;7;4]') s = A\b s = -14/15+3/5*a -3/5+3/5*a 23/15-1/5*a

  42. 2、符号方程（组）求解 s = solve(eq,var) 求方程指定变量var的解 s = solve(eq) 求方程默认自变量的解 eq是用字符串表示的方程，或符号表达式。 (1)用 字符串创建方程 eq1 = 'a*x^2+b*x+c=d' sx = solve(eq1) sa = solve(eq1,'a')

  43. eq1 = 'a*x^2+b*x+c-d' sx = solve(eq1) sa = solve(eq1,'a') 若方程 eq 中不含等号，则表示解方程 eq = 0 sx = -1/2*(b-(b^2-4*a*c+4*a*d)^(1/2))/a -1/2*(b+(b^2-4*a*c+4*a*d)^(1/2))/a sx是2行1列的符号数组 sa = -(b*x+c-d)/x^2 以a为方程的变量求出的解

  44. (2)用符号表达式创建方程 syms a b c d x eq1 = a*x^2+b*x+c-d; sx = solve(eq1) sa = solve(eq1,a) 符号表达式不能写成： eq = a*x^2+b*x+c=d

  45. 解方程 syms x; f = x^3-3*x+1; s = solve(f,x) s = solve('x^3-3*x+1','x') s = solve('x^3-3*x+1=0','x') 上面的三种方法都可以求得方程的解。

  46. solve也可以解方程组 s = solve( f1 , f2 , ... , fn , v1, v2 , ... , vn) 求解由 f1 , f2 , ... , fn确定的方程组关于 v1, v2 , ... , vn的解 若不提供v1, v2, ... , vN ，Matlab会自动确定默认变量。 返回值 s是一个架构数组，如果要显示(引用)每个解，使用 s.v1，s.v2

  47. eq1 = 'x + y - z = 1' eq2 = 'x + z = 2' eq3 = 'x^3 + y = 2' s = solve(eq1,eq2,eq3) x = s.x y = s.y z = s.z s = x: [3x1 sym] y: [3x1 sym] z: [3x1 sym] x = 1 1/2*5^(1/2)-1/2 -1/2*5^(1/2)-1/2 y = 1 -5^(1/2)+4 5^(1/2)+4 z = 1 -1/2*5^(1/2)+5/2 1/2*5^(1/2)+5/2 x(1),y(1),z(1)是方程组的第 1 组解 x(2),y(2),z(2)是方程组的第 2 组解 x(3),y(3),z(3)是方程组的第 3 组解

  48. clc； clear all; syms u v w y z eq1 = u*y^2 + v*z + w; eq2 = y + z + w; s = solve(eq1,eq2,y,z)

  49. 求范德瓦尔斯气体的三个临界参量 clc; clear all; syms a b R T V P = R*T/(V-b) - a/V^2 PV = diff(P,V) PV2 = diff(PV,V) S = solve(PV,PV2,V,T) Vc = S.V Tc = S.T Pc = subs(P,{V,T},{S.V,S.T})

  50. solve 在得不到解析解时，会给出数值解 s = solve('(x+2)^x=2') s = .69829942170241042826920133106081 s = solve('x= cos(x)') s = .73908513321516064165531208767387

More Related