Transcript 数学软件选讲

数学软件选讲
• Mathematica
• Matlab
• SAS
第 一 篇
Mathematica
 基础知识
 作为一门新的编程语言
 图形处理（二维、三维及其参数方程的形式）
 极限、微分与积分
 求解方程（组）、微分方程（组）
 在线形代数方面的应用
 数值处理
 文件及其它高级操作
第一章 基础知识
一、Mathematica3.0 界面及运行介绍
二、基本数值运算
1. 整数运算：加、减、乘、除、幂、阶乘
2. 数学常量：E、Pi、I、Degree、Infinity
3. 函数及数学函数
4. 浮点数及复数运算：N函数
三、变量及表达式
1. 变量的定义及清除
◆ 变量的特点
(1) 变量的默认作用域是全局的
(2) 全局变量不需事先定义或声明
(3) 尽量避免使用下划线定义变量
2. 多项式及其操作
(1) 定义、替换符操作
(2) 常用操作：
Expand、Factor、Together、Part
Simplify、Collect、Coefficient、
Exponent
四、序列及其操作
1. 序列的定义
2. 序列的生成：Table函数
3. 序列的操作
(1) 添加删除：Append、Prepend、Insert、
Delete、DeleteCases
(2) 取元素：Part、Take、Drop、Select
(3) 检测：Length、Count、Position
五、表达式“头”的概念：
Head及Apply函数
六、自定义函数
1. 一元函数
例： Clear[f,x]
f[x_]:= x^2+4x-2
2. 多元函数
例： f[x_,y_]:=
x^2+y^2-3
3. 迭代函数
例：f[n_]:=
f[n-1]+f[n-2];
f[0]= 1; f[1]=1;
第二章 编程语言
1· 条件语句
◆ 逻辑判断符
== >= <= >
===
=!=
<
!=
◆ 逻辑运算符
!
||
&&
◆ /；运算符
x = a /;test
仅当test为True时才执行赋值语句
◆ If 语句
语法：If [test, then, else]
若test为 True，则执行then，若test为
False，则执行else.
◆ Which 语句
语法：Which [test1, value1, test2,…]
依次计算testi，给出对应第一个test
为True 的value
◆ Switch[expr,form1,value1,form2,…]
比较expr与formi，给出与第一个form
值匹配的value
例1. 定义如下的函数：
0

 x
 2
 x
x  0
0 x2
x2
① 使用 /; 定义：
f [x_]:= 0 /;x<=0
f [x_]:= x /; x>0&&x<=2
f [x_]:= x^2 /; x>2
② 使用 If 定义：
f [x_]:= If [ x<=0, 0, If [x>2, x^2, x ] ]
③ 使用Which定义：
f [x_]:= Which [ x<=0, 0, x>2, x^2, True, x ]
2· 输出语句Print
3· 循环语句
◆ Do 语句
语法：Do[expr, {i, imin, imax, di}]
计算expr，i=imin,…,imax，步长为di
◆ While 语句
语法：While[test, body]
当test为True时，计算body
◆ For 语句
语法：For[start, test, incr, body]
以start为起始值，重复计算body和
incr，直到test为False时为止
◆ 循环控制语句Break和Continue
Break[]
退出最里面的循环
Continue[]
转入当前循环的下一步
第三章 图形处理
1. 基本二维图形
① Plot[ f, { x, xmin, xmax}]，用于
绘制形如y =f (x)的函数的图形。
当将多个图形绘制在同一坐标系上时，
形如： Plot[{ f1,…, fn},{x, xmin,
xmax}]
注意：有时需要使用Evaluate函数。
例：在同一坐标系下绘出
sinx, sin2x, sin3x, sin4x,
sin5x
的图形。
常用的选项：
PlotStyle－>Hue[a]
设置线条颜色
PlotRange－>{a,b} 控制显示范围
DisplayFunction
控制图形显示
AspectRatio
比
图形的宽、高
例：有如下的抛物线簇：
gx 2 sec 2 
y  (tan  ) x 
2v02
( g  9.8，v0  200)
当从15 变化到75，以15 为间隔时，绘出这组图形
程序：
Clear[a,y,x]
v=200;g=9.8;
y[a_,x_]:=Tan[a]*x-g*x^2*Sec[a]^2/(2v^2)
Plot[Evaluate[Table[y[i,x],{i,Pi/12,5Pi/12,
Pi/12}]],{x,0,4000}]
② ListPlot [List]，用于绘制散点
图。
{{注意，List的形式应为：
x0 , y0 },{x1 , y1},,{xn , yn }}
例：在同一坐标系下绘制下列两组散点图
p1={{0,0},{0,45},{5.3,89.6},{22.6,131.2}};
p2={{0,0},{2.68,44.8},{12.57,88.28},{27,130.3}};
程序：
g1=ListPlot[p1,PlotJoined->True,
DisplayFunction -> Identity];
g2=ListPlot[p2,PlotJoined -> True,
DisplayFunction -> Identity];
Show[g1,g2,DisplayFunction -> $DisplayFunction];
③ ParametricPlot [{ fx ,
fy},{t,tmin,tmax}]
用于绘制形如{x = fx(t) , y =
fy(t)}的参数方程图形。
例：绘制以点（3,4）为圆心，半径为2的圆。
ParametricPlot[{3+2Cos[t],4+2Sin[t]},
{t,0,2Pi}]
可增加如下选项：
AspectRatio->1,
AxesOrigin->{0,0}
2. 其它二维图形
① ContourPlot[ f, {x,xmin,xmax},
{y,ymin,
ymax}]，用于绘制形如z =f (x, y)的函
数的等高线图。
② DensityPlot[ f, {x,xmin,xmax},
{y,ymin,
ymax}]，用于绘制形如z =f (x, y)的函
数的密度图。
3. 三维图形
① Plot3D[ f,{x,xmin,xmax},{y,ymin,ymax}]
绘制形如Z = f (x, y)的三维图形。
例：绘制以下的函数图形：
Z = 10sin(x+siny)
命令：Plot3D[10
Sin[x+Sin[y]],{x,-10,10},
{y,-10,10}]
可增加选项：
PlotPoints->40
② ParametricPlot3D [{ fx , fy ,
fz},
{t,tmin,tmax} ,{u,umin,umax}]
用于绘制形如{x = fx(t) , y = fy(t) ,
z = fz(t)}的参数图形。
2
v0
g 2
例：画出抛物线y 
 2 x ( x  0, y  0)绕
2 g 2v0
y轴旋转60 所得的图形。
( g  9.8, v0  200)
解：旋转所得的抛物面参数方程为：
 x  r cos
2

 z  r sin  ，其中  [0,  ], r  [0,
3
 y  a  br 2

2
0
v
g
其中a 
,b  2
2g
2v0
a
]
b
4. 利用函数包绘制特殊图形
载入图形函数包的方法：
<<类名`包名`
例：<<Graphics`Graphics`
PolarPlot[r,{t,tmin,tmax}]
绘制极坐标图形
LogPlot[f,{x,xmin,xmax}]
画对数线性图
BarChart[list]
画出list的条形图
PieChart[list]
画出list的百分图
例：<<Graphics`ImplicitPlot`
ImplicitPlot[eqn,{x,xmin,xmax}]
绘制形如
f (x,y)=0的隐函数图形
例：绘制以点（3,4）为圆心，半径为2的圆。
ImplicitPlot[(x-3)^2+(y-4)^2==2,{x,0,5}]
第二章
幂级数、极限、微分与积分
1. 幂级数展开
Series[expr,{x, xo ,n}]
次的幂级数展开
求在点 x=xo 处至多n
例：求ex 在点 x=0处 x4 级幂级数展开
注：使用Normal函数可以去掉级数中的极小
项，从而转变成一般表达式。
2. 极限
Limit[expr,x-> xo] 求 x 逼近 xo时expr的极限
sin x
例：lim x
x  0
lim
x 2  3x  x
x   
某些函数在一点处的极限随逼近方向不同而
不同，可用Direction选择方向：
Limit[expr, x-> xo, Direction -> 1] 左极限
Limit[expr, x-> xo,Direction-> -1] 右极限
例：求1/x 的左右极限
3. 微分
D[ f ,{x,n}]
求f 的n阶偏微分
Dt[ f ]
求f 的全微分
例：D[x^n,{x,3}]
Dt[x^2+y^2]
例：y = xarctgx，求其100阶导数及其在0
点的值
4. 积分
求f 的不定积分
Integrate[ f , x]
Integrate[ f ,{x,xmin,xmax}]
求 f 的定积分
Integrate[ f ,{x,xmin,xmax}, {y,ymin,ymax}]
求 f 的多重积分
1
x2  a
例：

2
1
1
x2  a
 sin(sin
2
 sin(sin
0
x)
x)
第三章 线性代数
1. 构造矩阵和向量
Table[ f ,{i,m} ,{j,n}]
构造m×n矩阵，f 是
i, j的函数，给出[i, j]项值
Array[ f ,{m, n}]
构造m×n矩阵，[i, j]
项的值是 f [i, j]
DiagonalMatrix[ List] 生成对角线元素为
List的对角矩阵
IdentityMatrix[n]
构造n阶单位阵
2. 截取矩阵块
M[[i]]
取矩阵M的第 i 行
Map[#[[i]]&, M]
取矩阵M的第 i 列
M[[i, j ]]
取矩阵M的i, j 位置的元素
M[[{i1,…,ir}, {j1,…,js}]]
矩阵M的r×s子
矩阵，元素行标为ik，列标为jk
M[[Range{i0,i1}, Range{j0,j1}]]
矩阵M的从
i0到i1行， j0到j1列元素组成的子矩阵
3. 矩阵及向量的运算
M.N
对M、N做矩阵乘法（向量内积）
M*N
将M、N的对应位置元素相乘
Outer[Times,M,N] 求M、N的外积
Dimensions[ M ]
给出矩阵M的维数
Transpose[ M ]
转置
Inverse[ M ]
求逆
Det[ M ]
方阵M的行列式值
MatrixPower[M,n]
n阶矩阵幂
MatrixExp[M]
矩阵指数
Eigenvalues[ M ]
M的特征值
Eigenvectors[M]
M的特征向量
第四章
求解方程（组）、微分方程（组）
1. 求解多项式方程（组）
Solve[ eqns ,vars]
求解多项式方程
Solve[{eqn1,…eqnn}, {var1,…varn}]
求解多项式方程组
注：Solve只能给出多项式方程（组）的解，
因此它们只适用于幂次不高、规模不大的多
项式方程（组）。
NSolve[ eqns ,vars]
求多项式方程的数值解
NSolve[{eqn1,…eqnn}, {var1,…varn}]
求多项式方程组的数值解
对于数值解，可以直接用NSolve求解
例：求解以下方程（组）
x2+ax=2
x3+34x+1=0
x5-1331x+11= 0
ax  by  1

 x y 2
 x 3  y 3  xy

 x  y  xy  1
2. 求解微分方程（组）
DSolve[ eqns ,y[x], x]
求解y[x]的微分方程
DSolve[ eqns ,y, x]
以纯函数的形式给出y的解
DSolve[{eqn1, eqn2,…}, {y1, y2, …}, x]
求解微分方程组
例：求解以下微分方程（组）
y’ = y
 x'   y

 y'  x
y’’－ k y =1
 y'  2 y
并求x  2.8时函数y的值

 y (0)  1
第五章 数值处理
1. 数值积分
NIntegrate[expr , {x,xmin,xmax}]
注意，NIntegrate直接计算数值积分，不先给
出符号结果，而Integrate[…]//N会尽可能的
先求精确解的形式。
2. 数值根求解
FindRoot[lhs=rhs , {x, x0}]
以x0为初始点求方程的数值解
FindRoot[lhs=rhs , {x, {x0 ,x1}}] 给出两个
初值求数值根（方程的符号导数无法求出
时，必须使用此形式）
FindRoot[{eqn1, eqn2,…}, {x, x0},{y, y0 }, …]
对联立方程 eqni 求数值解
例：求解下列方程（组）
cosx =x
sin x  cos y

 x  y 1
x600+5x+3=0
3. 微分方程数值解
NDSolve[{eqn1, eqn2,…}, y,{x,xmin,xmax}]
求函数y的数值解，x的范围为{xmin,xmax}
NDSolve[{eqn1, eqn2,…},{y1, y2, …},
{x,xmin,xmax}]
求函数yi的数值解
注：以上两种形式用于求解常微分方程（组）
NDSolve以InterpolatingFunction目标生成函
数yi的解。 InterpolatingFunction目标提供独
立变量x在xmin到xmax范围内yi的近似值。
例：求解以下微分方程（组）并画出函数y的图形
 dy 2
2
( )  y  0
在区间[0,1]上
 dx
 y 2 x 0  4
 dx
2
 dt   y  x
 dy
在区间[0,10]上
  2x  y
 dt

 x t 0  y t 0  1
NDSolve[{eqn1, eqn2,…}, y,{x,xmin,xmax},
{t,tmin,tmax}] 求由函数y构成的偏微分方
程的数值解
NDSolve[{eqn1, eqn2,…},{y1, y2, …},
{x,xmin,xmax} ,{t,tmin,tmax}]求由函数yi构
成的偏微分方程组的数值解
例：求下面微分方程的数值解并绘图。
 2 y 2 y
 2  2
x
 t


 x2
 y t 0  e


 y x  5  y x 5

 y
 t 0  0
 t
, x  [5,5], t  [0,5]
4. 极大极小值
ConstrainedMax[ f, {inequalities}, {x, y, …}]
ConstrainedMax[ f, {inequalities}, {x, y, …}]
求由目标函数 f 和不等式约束inequalities
构成的线形规划
例：ConstrainedMax[x+y,{x<1,y<2},{x,
y}]
LinearProgramming[ c, m, b]
求使cx在
约束mx≥b和x≥0下取最小值的矢量x
FindMinimum[ f, {x, x0}]
求函数的局部极小值
以x0为初始点，
注：FindMinimum的用法与FindRoot完全相同。
例：求100 cos t  300  100 sin t
2
2
2
在
t  [0,2 ]内的极小值
例：求x 4  3x 2 y  5 y 2  x  y在x, y  [0,1]
内的极小值
5. 曲线拟合
Fit [ data, funs, vars]
用变量为vars的函数funs拟合一组数据data
第 二 篇
Matlab
第一章 矩阵及其基本运算
一、矩阵的表示
1．实数值矩阵生成
2．复数矩阵生成
3. 符号矩阵的生成
用sym函数或syms函数
4. 大矩阵的生成
.m文件及函数的定义
5. 特殊矩阵的生成
全零阵、全1阵、单位阵：zeros,eye,ones
随机矩阵：
均匀分布：
rand
标准正态分布：
randn
线性等分向量：
linspace
Hilbert矩阵：
hilb
魔方矩阵：
magic
二、矩阵操作
1．取矩阵中的元素
2．增加及删除矩阵中的元素
3．矩阵的旋转与变形
三、矩阵运算
1. 加减法运算
2. 乘法运算
① 矩阵乘法
② 数组乘法（数乘）
③ 向量内积、外积、叉乘
*④
矩阵的卷积与解卷、张量积
3. 集合运算
并：union
返回a、b的并集，即c = a∪b
交：intersect
返回向量a、b的公共部分，即c= a∩b
差：setdiff
返回属于a但不属于b的不同元素的集合，C = a-b
交集的非： setxor
检测集合中的元素： ismember
4. 除法运算
左除：\
右除：/
x=A\b是方程 Ax =b的解
x=b/A是方程 xA=b的解。
5. 矩阵乘方
6. 矩阵函数
expm
logm
sqrtm
7. 方阵的行列式：
8. 方阵的逆：
inv
9. 矩阵的迹：
trace
10. 矩阵的秩：
rank
det
11. 矩阵和向量的范数
norm
欧几里德范数
norm(x,inf )
无穷范数
12. 其它运算
四、矩阵分解
1．LU分解：
[L,U]=lu(X)
U为上三角阵，L为下三角阵或其变换形式，满足LU=X
2．QR分解：
[Q,R]=qr(A)
求得正交矩阵Q和上三角阵R，Q和R满足A=QR
3．特征值分解
[V,D]=eig(A)
计算A的特征值对角阵D和特征向量V，使AV=VD成立
五、其它
二次型、秩与线性相关性、稀疏矩阵
第二章 Matlab语言基础
一、M文件
1．脚本文件：在Matlab的工作空间内
对数据进行操作。
2．函数文件：可接受输入参数并返回
输出参数，其内的变量不占用Matlab工作空
间，第一行包含function
注： M文件的调用以文件名为准。
%为Matlab的注释符，其后的语句
不执行（只对当前行有效）。
二、Matlab语言
1．逻辑判断符
>= <= >
isequal函数
<
2．逻辑运算符
&
|
~
3．条件语句
① if-else语句
② switch-case语句
==
~=
4．循环语句
① for语句
② while语句
三、编程技巧
1.调试程序
2.输入输出参数
nargin、nargout
第三章 Matlab图形处理
一、二维图形
1. 基本二维图形
Plot
用法如下：
a.
Plot (X)
b.
Plot (X,Y)
c.
Plot (X1,Y1,X2,Y2,…)
d.
Plot (X1,Y1,LineSpec1,X2,Y2, X3,Y3, …)
其中参数LineSpec定义线条的属性。Matlab
中可以对线条定义如下的特性：
a.
线型： -(实线) -- (划线) :(点线) -. (点划线)
b.
线条宽度： LineWidth
c.
颜色
d.
标记类型
e.
标记大小：Markersize
fPlot
在指定的范围limits内画出一元函数y=f (x)的图形
用法：fplot('function',limits)
注意：函数function必须是一个M文件函数
或者是一个包含变量 x，且能用函数eval计算的
字符串。
例：在同一坐标系下绘制tgx和的sinx图形
fplot(‘[tan(x),sin(x)]’,[-1,1,0,2*pi])
注意坐标系调整函数axis的作用和用法
2. 图形标注
title
为图形添加标题
xlabel
为x轴加标注
ylabel
为y轴加标注
text
在指定位置上添加文本字符串
gtext
用鼠标在图形上放置文本
legend
为图形添加图例
3. 特殊二维图形
画极坐标形式函数r = f (θ)的极坐标图
polar
用法如下：
polar(theta,rho,LineSpec)
例：
t = 0:.01:2*pi;
polar(t,sin(3*t).*cos(2*t),'--r')
4. 其它二维图形
pie
用x中的数据画一饼形图
semilogx x轴对数图形
loglog
bar
双对数图形
用二维垂直条形显示向量或矩阵中的值
barh 用二维水平条形显示向量或矩阵中的值
hist
二维条形直方图，可以显示出数据的分
配情形
二、三维图形
1. 曲面与网格图形命令
mesh
生成由X，Y和Z指定的网线面
在使用该命令前应先用meshgrid函数生成可用
于计算函数值的矩阵网格。
通常用法如下：
[X,Y ]=meshgrid(a)
Z= f (X,Y)
mesh(X,Y,Z)
2. 三维图形的其它形式
contour
曲面的等高线图
pie3
三维饼图
surf
在矩形区域内显示三维带阴影曲面图
quiver
矢量图或速度图
surfnorm
计算与显示三维曲面的法线
第四章 Matlab应用
一、多项式运算
二、极限
limit (F, x, a, ‘right’ )
x趋向于a时F的极限
三、导数
diff (S, v, n)
四、积分
1. 符号积分
a. 不定积分
int (S, v)
b. 定积分
int (S, v, a,b)
2. 数值积分
a. 一元函数
quad ( fun,a,b)
自适应Simpson法
trapz ( X, Y )
梯形法
b. 二元函数
dblquad ( fun,xmin,xmax,ymin,ymax)
在矩形区域[xmin,xmax,ymin,ymax]上计算二
元函数z=f (x,y)的二重积分
quad2ggen ( fun,xlower,xupper,ylower,yupper)
在任意区域[xlower,xupper,ylower,yupper]上计
算二元函数z=f (x,y)的二重积分
五. 插值
a. interp1( X,Y,xi,method)
一维数据插值
b. interp2( X,Y,Z,xi, yi,method) 二维数据插值
例：已知1900年到2010年每隔十年的数据如下：
75.995
91.972
105.711
123.203
131.669
150.697
179.323
203.212
226.505 249.633
256.344
267.893
用插值法求1995年的数据。
六、方程（组）求解
1. 方程（组）的符号解
solve (eq)
求方程的符号解
solve (eq1,eq2,…eqn) 求方程组的符号解
例：
solve('x^2+3x-6')
solve('-x^2*y+3*x-6','x+y^2-1')
2.方程（组）的数值解
fzero (fun,x0)
用数值方法求方程根
fsolve(fun,x0)
用数值方法求方程根
例：求下列方程的根
2 x1  x 2  ex1
 x1  2 x 2  e x 2
解：先建立方程函数文件，并保存为myfun.m
function F = myfun(x)
F = [2*x(1) - x(2) - exp(-x(1));
-x(1) + 2*x(2) - exp(-x(2))];
然后调用优化程序
x0 = [-5; -5];
% 初始点
[x,fval] = fsolve(@myfun,x0,options)
七、积分变换
1. Fourier积分变换
对符号单值函数 f 中的缺
省变量 x（由命令findsym确定）计算Fourier变换形式
F = fourier( f )
例：
syms x w u v
f = sin(x)*exp(-x^2)
F = fourier(f)
注：用eval函数计算得出的表达式
f = ifourier(F) 逆Fourier积分变换
Y = fft(X)
快速Fourier变换
2. Laplace变换
输出参量L = L(s)为有缺
省符号自变量t的标量符号对象F的Laplace变换
L = laplace(F)
例：
syms x s t v
f1= sqrt(t);
L1 = laplace(f)
F = ilaplace(L)
逆Laplace变换
3. Z变换
F = ztrans(f )
计算z-变换
对缺省自变量为n的单值函数f
八、求解微分方程（组）
1. 常微分方程（组）符号解
dsolve(eq1,eq2,… )
缺省独立变量为t
例：
dsolve(‘Dy=1+y^2’,’y(0)=1’)
dsolve('D3u=u','u(0)=1','Du(0)=-1',
'D2u(0)=pi')
2. 常微分方程（组）数值解
ode45、ode23、ode113、ode15s、ode23s、de23t、
ode23tb
3. 偏微分方程数值解
① assempde
单的Poission方程（一类特
殊的椭圆型方程），能求解的方程形如：
，
，
 2u  1

 u G  0
② hyperbolic
型方程：
G  {(x, y) x2  y2  1}
仅能求解如下形式的双曲
 2u  2u  2u  2u 
 t 2   x 2  y2  z2   0



u t 0  0


u
t 0  0

t

G  {(x, y, z) 0  x, y, z  1}
③ parabolic
型方程：
，
仅能求解如下形式的抛物
u   2u  2u  2u 
  2  2  20
y
z 
 t  x
 u G  0
G  {(x, y, z) 0  x, y, z  1}
九、极值问题（优化工具箱）
1. 无条件极值问题
fminu ( fun, x0 ,options)
2. 条件极值问题
constr ( fun, x0 ,options)
3. 有界条件问题
constr ( fun, x0 ,options, VLB, VUB)