数 学 软 件 选 讲
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,定义如下的函数:



2
20
00
2
xx
xx
x
① 使用 /; 定义:
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 图形的宽、高比
AxesOrigin 设置原点坐标形为间隔时,绘出这组图,以变化到从当

157515
)2008.9(
2
s e c
)( t a n 0
2
0
22
vg
v
gx
xy
程序:
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的形式应为:
}},{,},,{},,{{ 1100 nn yxyxyx
例:在同一坐标系下绘制下列两组散点图
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)的函数的密度图。
例:绘制函数 f=sinx·siny的等高线图和密度图
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)}的参数图形。
)200,8.9(60
)0,0(
22
0
2
2
0
2
0


vgy
yxx
v
g
g
v
y
所得的图形。轴旋转绕例:画出抛物线
2
0
2
0
2
2
,
2
],0[],
3
2
,0[s i n
c os
v
g
b
g
v
a
b
a
r
bray
rz
rx



其中
,其中参数方程为:解:旋转所得的抛物面

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的极限某些函数在一点处的极限随逼近方向不同而不同,可用 Direction选择方向:
Limit[expr,x-> xo,Direction -> 1] 左极限
Limit[expr,x-> xo,Direction-> -1] 右极限例:求 1/x 的左右极限
x
x
x
s inlim
0
xxx
x


32l i m例:
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,积分
Integrate[ f,x] 求 f 的不定积分
Integrate[ f,{x,xmin,xmax}] 求 f 的定积分
Integrate[ f,{x,xmin,xmax},{y,ymin,ymax}]
求 f 的多重积分例:
ax 2 1? )s in ( s in x
21 2 1 ax?
2
0 )s in ( s in 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只能给出多项式方程(组)的解,
因此它们只适用于幂次不高、规模不大的多项式方程(组)。


2
1
yx
byax
NSolve[ eqns,vars] 求多项式方程的数值解
NSolve[{eqn1,…eqnn},{var1,…varn}]
求多项式方程组的数值解对于数值解,可以直接用 NSolve求解例:求解以下方程(组)
x2+ax=2 x3+34x+1=0
x5-1331x+11= 0


1
33
xyyx
xyyx
2,求解微分方程(组)
DSolve[ eqns,y[x],x] 求解 y[x]的微分方程
DSolve[ eqns,y,x]
以纯函数的形式给出 y的解
DSolve[{eqn1,eqn2,…},{y1,y2,…},x]
求解微分方程组例:求解以下微分方程(组)
y’ = y y’’- k y =1

xy
yx
'
' 的值时函数并求 yx
y
yy
8.2
1)0(
2'
第五章 数值处理
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 x600+5x+3=0

1
c o ss in
yx
yx
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的近似值。
上在区间 ]1,0[
4
0)(
0
2
22


xy
y
dx
dy
例:求解以下微分方程(组)并画出函数 y的图形上在区间 ]10,0[
1
2
00
2



tt
yx
yx
dt
dy
xy
dt
dx
NDSolve[{eqn1,eqn2,…},y,{x,xmin,xmax},
{t,tmin,tmax}] 求由函数 y构成的偏微分方程的数值解
NDSolve[{eqn1,eqn2,…},{y1,y2,… },
{x,xmin,xmax},{t,tmin,tmax}]求由函数 yi构成的偏微分方程组的数值解例:求下面微分方程的数值解并绘图。
]5,0[],5,5[,
0
0
55
0
2
2
2
2
2


tx
t
y
yy
ey
x
y
t
y
t
xx
x
t
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完全相同。
内的极小值在例:求
]2,0[
s i n100300c o s100 222


t
tt
内的极小值在例:求 ]1,0[,53 224 yxyxyyxx
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,方阵的行列式,det
8,方阵的逆,inv
9,矩阵的迹,trace
10,矩阵的秩,rank
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,特殊二维图形
polar 画极坐标形式函数 r = f (θ )的 极坐标图用法如下:
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,…eq n) 求方程组的符号解例:
solve('x^2+3x-6')
solve('-x^2*y+3*x-6','x+y^2-1')
2.方程(组)的数值解
fzero (fun,x0) 用数值方法求方程根
fsolve(fun,x0) 用数值方法求方程根例:求下列方程的根
2
1
x
21
x
21
ex2x
exx2


解:先建立方程函数文件,并保存为 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 = fourier( f ) 对符号单值函数 f 中的缺省变量 x(由命令 findsym确定)计算 Fourier变换形式例:
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 = laplace(F) 输出参量 L = L(s)为有缺省符号自变量 t的标量符号对象 F的 Laplace变换例:
syms x s t v
f1= sqrt(t);
L1 = laplace(f)
F = ilaplace(L) 逆 Laplace变换
3,Z变换
F = ztrans(f ) 对缺省自变量为 n的单值函数 f
计算 z-变换八、求解微分方程(组)
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方程(一类特殊的椭圆型方程),能求解的方程形如:


0u
1u
G
2
}1yx)y,x{(G 22
② hyperbolic 仅能求解如下形式的双曲型方程:


0
t
u
0u
0
z
u
y
u
x
u
t
u
0t
0t
2
2
2
2
2
2
2
2
}1z,y,x0)z,y,x{(G
③ parabolic 仅能求解如下形式的抛物型方程:






0u
0
z
u
y
u
x
u
t
u
G
2
2
2
2
2
2
}1z,y,x0)z,y,x{(G
九、极值问题(优化工具箱)
1,无条件极值问题
fminu ( fun,x0,options)
2,条件极值问题
constr ( fun,x0,options)
3,有界条件问题
constr ( fun,x0,options,VLB,VUB)