Matlab程序设计打印版
第五讲MATLAB程序设计ppt课件
语句组m
otherwise
语句组n
end
(exswitch.m)
第五讲 MATLAB程序设计
18
(3)try语句 语句格式为: try
语句组1 catch
语句组2 end
try语句先试探性执行语句组1,如果语句组1在执行 过程中出现错误,则将错误信息赋给保留的lasterr 变量,并转去执行语句组2。
第五讲 MATLAB程序设计
14
2、选择结构
(1) 条件分支语句——if语句 在MATLAB中,if语句有3种格式。 1) 单分支if语句: if 条件 语句组
end
第五讲 MATLAB程序设计
15
2) 双分支if语句: if 条件
语句组1 else
语句组2 end
第五讲 MATLAB程序设计
16
第五讲 MATLAB程序设计
24
三、程序调试
1 错误分类
一般来说,应用程序的错误有两类:
一类是语法错误,例如函数名的拼写错、表达式 书写错等。
另一类是运行时的错误。指程序的运行结果有错 误,这类错误也称为逻辑错误。
第五讲 MATLAB程序设计
25
2、查找逻辑错误的方法:
◆ 删去语句行末的分号,使显示其运行中间结果 ◆ 利用keyboard 命令实现,return继续程序执行 ◆ 注释掉M 函数文件的函数定义行,使函数文件转
第五讲 MATLAB程序设计
19
例: 矩阵乘法运算要求两矩阵的维数相容,否则会 出 错。先求两矩阵的乘积,若出错,则自动转去 求两矩阵的点乘。(extry.m)
第五讲 MATLAB程序设计
20
3、 循环结构
(1)硬循环语句——for语句
MATLAB程序设计
MATLAB程序设计MATLAB提供了一个完善的程序设计语言环境,使用户能够方便地编制复杂的程序,完成各种计算。
本节先介绍关系运算、逻辑运算,再介绍M-文件(即程序文件)的结构及MATLAB的程序控制流语句。
一、关系运算和逻辑运算1.关系运算(1)关系运算符:< ;< = ;> ;> = ;= = ;~ =(2)关系表达式:用关系运算符将两个同类型的量(表达式)连接起来的式子。
【注】①关系运算本质上是标量运算,关系表达式的值是逻辑值(0-假1-真);②当作用于两个同样大小矩阵时,则分别对两个矩阵的对应元素运算,结果是一个0-1矩阵。
例1.对向量进行关系运算。
>> A=1:5,B=5:-1:1 % 输入向量A = 1 2 3 4 5B = 5 4 3 2 1>> C=(A>=4) % 对向量进行关系运算C = 0 0 0 1 1>> D=(A==B) % 对向量进行关系运算D = 0 0 1 0 02.逻辑运算(1)逻辑运算符:& (and,与)、| (or,或)、~ (not,非)(2)逻辑表达式:用逻辑将两个逻辑量连接起来的式子。
【注】①逻辑运算本质上是标量运算,它将任何非零元素视为1(真);②当作用于两个同样大小矩阵时,则分别对两个矩阵的对应元素运算,结果是一个0-1矩阵。
(真值表见P27)例2.对向量进行逻辑运算。
>> a=1:9,b=9-aa = 1 2 3 4 5 6 7 8 9b = 8 7 6 5 4 3 2 1 0>> c=~(a>4) % 非运算c = 1 1 1 1 0 0 0 0 0>> d=(a>=3)&(b<6) % 与运算d = 0 0 0 1 1 1 1 1 13.逻辑函数any(x) 向量x 中有非零元返回1,否则返回0。
(向量函数) all(x) 向量x 中所有元素非零返回1,否则返回0。
《MATLAB程序设计》复习资料,DOC
Matlab习题及复习要点1.Matlab的英文组成;程序和变量的命名规则;在MATLAB语言中是区分字母大小的,也就是说,大写字母和小写字母代表的东西是不同的。
234510.读懂逻辑表达式,会写出逻辑表达式的结果(0或1)11.掌握集中循环结构,尤其if..elseif…else…end和swich…case结构的语法,要准确。
12.会编写分段函数的程序;x和y满足如下关系:编写函数y=f(x),用于计算上述分段函数。
13.绘图时采用的几个命令的掌握:holdon、plot、plot3 14*.用语句[x,y]=meshgrid(a:b)构建网格数据,例子如下:第一讲概论1.简述matlab基本特点。
(0.5分)交互式操作界面;高效的数值计算功能;演算式语言;可视化输出;代码、数据文件的集成管理环境;支持用户界面开发,自定义创建工具(GUIDE);丰富的外部接口——支持C/C++、Java、Excel/Word、Ansys,COM、DDE(动态数据交换)和ActiveX……。
删除工作空间的变量a:cleara;清空工作空间:clear或clearall;删除命令行:esc;查询函数sin的帮助文档:helpsin;1.分别用直接输入法和存储变量法求1+cos(pi)*(2+2i)。
直接输入法:>>1+cos(pi)*(2+2i); 存储变量法:>>a=cos(pi);>>b=2+2i;>>c=1+a+b;2.a=int8(100),b=int8(50)a+b=127;a-b=50;第三讲数组1.生成一个3*3随机矩阵,将其对角线元素的值加1。
(写出代码)rand(3)+eye(3)1.生成一个元素值在1和10之间的3*3随机矩阵,将其重新排列,使得:(1)每列按降序排列;(2)每行按降序排列。
(3)C<=D=[0,0;1,1].(10)已知A为如下4*4矩阵:则运行B=A([1:2],[1:2])后,B为2行2列矩阵,其值为__[12;56]_______。
Matlab程序设计1
x=linspace(0,2*pi,60);
y=sin(x);z=cos(x);
t=sin(x)./(cos(x)+eps); ct=cos(x)./(sin(x)+eps);
subplot(2,2,1);
%选择2×2个区中的1号区
stairs(x,y);title('sin(x)-1');axis ([0,2*pi,-1,1]);
当前窗口分成几个矩形部分,不同部分按 行方向以数字进行标号。调用格式如下 :
subplot(m,n,p)%将一个窗口分成m×n 个小窗口,在第p个小窗口中创建坐标轴。
例如: 在一个图形窗口中以子图形式同 时绘制正弦、余弦、正切、余切曲线。 程序如下:
x=linspace(0,2*pi,60);
y=sin(x);z=cos(x);
t=sin(x)./(cos(x)+eps); ct=cos(x)./(sin(x)+eps);
subplot(2,2,1); plot(x,y);title('sin(x)');axis([0,2*pi,-1,1]); subplot(2,2,2); plot(x,z);title('cos(x)');axis([0,2*pi,-1,1]); subplot(2,2,3); plot(x,t);title('tangent(x)');axis([0,2*pi,40,40]); subplot(2,2,4); plot(x,ct);title('cotangent(x)');axis([0,2*pi ,-40,40]);
二、 绘制图形的辅助操作
1. 图形注释:通过选择图形窗口主菜单
(打印)实验四 MATLAB 高级图形绘制
实验四MATLAB 高级图形绘制一、实验目的及要求:1.熟悉各种绘图函数的使用;2.掌握图形的修饰方法和标注方法;3.了解MATLAB 中图形窗口的操作。
二、实验内容:1.用图形表示连续调制波形Y=sin(t)sin(9t)及其包络线。
程序代码如下:包络线:2.x=[-2π,2π],y1=sinx、y2=cosx、y3=sin2x、y4=cos 2x①用MATLAB语言分四个区域分别绘制的曲线,并且对图形标题及横纵坐标轴进行标注。
程序:结果:②另建一个窗口,不分区,用不同颜色、线型绘出四条曲线,并标注图例注解。
程序:结果:③绘制三维曲线:⎪⎩⎪⎨⎧=≤≤==)cos()sin()200()cos()sin(t t t z t t y t x π程序:结果:3.绘制极坐标曲线ρ=asin(b+nθ),并分析参数a、b、n对曲线形状的影响。
(1)a=1;b=1;n=1(2)a=10;b=1;n=1(3)a=10;b=10;n=1 (4)a=10;b=10;n=10参数a、b、n对曲线形状的影响:由上面绘制的图形可知:a决定图形的大小,当a为整数时,图形半径大小就是a;b决定图形的旋转角度,图形的形状及大小不变;n决定图形的扇叶数,当n 为奇数时,扇叶数为n,当n为偶数时,扇叶数为2n。
三、结论本次实验用到了曲线绘图、三位曲线绘图的知识,与老师上课的内容一致,让我学的matlab绘图的知识得到了巩固,我还学会了如何使用title、subplot、plot、axis等函数。
在做实验的过程复习了hold on指令是覆盖函数继续绘图的意思。
Matlab命令汇总打印
Matlab命令汇总一、常用对象操作:除了一般windows窗口的常用功能键外。
1、!dir 可以查看当前工作目录的文件。
!dir& 可以在dos状态下查看。
2、who 可以查看当前工作空间变量名,whos 可以查看变量名细节。
3、功能键:功能键快捷键说明方向上键Ctrl+P 返回前一行输入方向下键Ctrl+N 返回下一行输入方向左键Ctrl+B 光标向后移一个字符方向右键Ctrl+F 光标向前移一个字符Ctrl+方向右键Ctrl+R 光标向右移一个字符Ctrl+方向左键Ctrl+L 光标向左移一个字符home Ctrl+A 光标移到行首End Ctrl+E 光标移到行尾Esc Ctrl+U 清除一行Del Ctrl+D 清除光标所在的字符Backspace Ctrl+H 删除光标前一个字符Ctrl+K 删除到行尾Ctrl+C 中断正在执行的命令4、clc可以命令窗口显示的内容,但并不清除工作空间。
二、函数及运算1、运算符:+:加,-:减,*:乘,/:除,\:左除^:幂,‘:复数的共轭转置,():制定运算顺序。
2、常用函数表:sin( ) 正弦(变量为弧度)Cot( ) 余切(变量为弧度)sind( ) 正弦(变量为度数)Cotd( ) 余切(变量为度数)asin( ) 反正弦(返回弧度)acot( ) 反余切(返回弧度)Asind( ) 反正弦(返回度数)acotd( ) 反余切(返回度数)cos( ) 余弦(变量为弧度)exp( ) 指数cosd( ) 余弦(变量为度数)log( ) 对数acos( ) 余正弦(返回弧度)log10( ) 以10为底对数acosd( ) 余正弦(返回度数)sqrt( ) 开方tan( ) 正切(变量为弧度)realsqrt( ) 返回非负根tand( ) 正切(变量为度数)abs( ) 取绝对值atan( ) 反正切(返回弧度)angle( ) 返回复数的相位角atand( ) 反正切(返回度数)mod(x,y) 返回x/y的余数sum( ) 向量元素求和3、其余函数可以用help elfun和help specfun命令获得。
数值计算与金融仿真matlab程序总结(双面打印)
变量合法命名:只含字母、数字、下划线,并以字母开头 who 显示定义变量 whos 显示定义变量及其信息 lookfor “?”查询与关键词有关的所有函数Integer 整数 real 实数 complex 复数 inf 无限大 NaN 非数字 format long 显示多位小数 short 4位 bank 2位 round()四舍五入 fix()去小数部分 prod() 连乘 asin()transpose()或矩阵 矩阵转置 eye(4)4行4列单位矩阵 q13(2,4) 矩阵第二行第四列数据 q13(2:3,2:3) 矩阵q13 2到3行,2到3列的数据 .*对应元素相乘 同理./ .^a :b :c 以a 开始,以b 为间隔,最大数小于等于c 的数列sum(A)对列求和(得行) sum(A ,2)对行 sum(sum(A))对列求和再对行求 prod(A)对列连乘 min(A)输出每一列的最小值 [m,j]=max(B),m=最大值取值,j=最大值位置 判断:a==b a 等于b a~=b a 不等于b 正确输出1 . plot(x,y, '(颜色)(形状)') linspace(a,b,n)把起点为a 终点为b 的直线等分为n 份 plot(x,y,'.',x,cos(3*pi*x),'g*')两条线一图 legend('Sin curve', 'Cos curve')加图标grid 加格子 hold on 保留原曲线,可用于画多线在一图subplot(121), plot(x,y) 将作图区域分为1行2列,作第1个区域的图 N = 100; h = 1/N; x = 0:h:1; y = cos(3*pi*x); plot(x,y) x = linspace (0,1,101); y = sin(3*pi*x); plot(x,y)rand(3)3行方阵随机矩阵 randn r6 =random('Poisson',6,1,6) 文档读写IBM = xlsread('IBM.xls');IBM = textread('IBM.txt');IBM = dlmread('IBM.txt','',1,1); IBM = importdata('IBM.xls');m 文件which mfile(文件名) 查看路径edit mfile 编辑脚本mfile 与函数mfile 的区别,前者是全局变量,后者是局部变量function [var1 var2 …] = functionname(arg1,arg2, …)多个图形x = -1:.05:1; for n = 1:8 subplot(4,2,n); plot(x,sin(n*pi*x)); end 求和s = 0 for i = 1:100 s = s + i; end sum 能以1^2 + 2^2 + … + n^2表示并小于100的?S = 1; n = 1; while S+ (n+1)^2 < 100 n = n+1; S = S + n^2; end [n, S] If...Then LoopS = 1; n = 1; for i=1:100;if S+(n+1)^2 < 100 S = S+ (n+1)^2; n = n+1; end end [n, S] Chapter 2 现金流分析pvvar 求变动现金流量的现值pvfix fvvar fvfix PV = pvvar(CashFlow, Rate, CFDates)> CashFlow = [100 300 450] CFDates = ['01/12/1987'; '02/14/1988'; '03/03/1988'] PV = pvvar(CashFlow, 0.09, CFDates) 日期列数据 >> pvvar([-15000 3000 4500 5000 6800], 0.08) ans = 603.1667 没有现金流日期时,默认相隔一年Irr 固定周期的内部收益率Xirr 变动周期的内部收益率 Return = irr(CashFlow) Return = xirr(CashFlow, Dates)>CF = [-10000 3000 4500 5000]; Return = irr(CF) >>Return =0.1105 >Dates=['01/12/00'; '02/14/01'; '09/03/01'; '12/31/02']; Return = xirr(CF, Dates) >>Return =0.1177 cfdur and cfconv 久期和凸性>CashFlow=[5 5 5 5 5 105];> [Dur, ModDur] = cfdur(CashFlow, 0.05)收益Dur =5.3295ModDur =5.0757 > Conv = cfconv(CF, 0.05) Conv =95.5410 债券价格和收益率 [P, I 应付利息] = bndprice(Yield, CouponRate 票面利率, ‘Settle ’交收日, ‘Maturity ’到期) p=price+accruedint Yield = bndyield(Price, CouponRate, Settle, Maturity) Duration = bnddurp(Price, CouponRate, Settle, Maturity) Duration = bnddury(Yield, CouponRate, Settle, Maturity) Convexity = bndconvp(Price, CouponRate, Settle, Maturity) p-y Price-yield curve >> yields = 0.01: 0.01: 0.20;>> [P I] = bndprice(yields, 0.1, '08/10/07', '12/31/20');>> plot (yields, P+I);>> grid on;>> xlabel('Yield');>> ylabel('Price') ;>> Title('Price-Yield Curve'); 利率免疫 duration=holding period >> settle = '28-Aug-2007';>> maturity = ['15-Jun-2012' ; '31-Oct-2017' ; '01-Mar-2027'];>> couponRate = [0.07 ; 0.06 ; 0.08];>> yield = [0.06 ; 0.07 ; 0.075];>> duration = bnddury(yield, couponRate, settle, maturity); >> convexity = bndconvy(yield, couponRate, settle, maturity);% COMPUTE PORTFOLIO WEIGHTS >> A = [duration'; convexity'; 1 1 1];>> b = [ 10; 160; 1];>> weights = A\b 有效前沿 [PortRisk, PortReturn, PortWts] = frontcon(ExpReturn 历史收益率, ExpCovariance 历史协方差, NumPorts 返回结果的个数, PortReturn 目标收益, AssetBounds 资产界限, Groups, GroupBounds 组合界限) r = [0.2 0.1];两种资产历史收益率> s = [0.2 -0.1; -0.1, 0.4];历史协方差 >> [Risk, Return, Wts] = frontcon(r, s, 5);>> [Risk, Return, Wts]ans =0.2958 0.1625 0.6250 0.3750 0.3075 0.1719 0.7187 0.2813 0.3400 0.1812 0.8125 0.18750.3883 0.1906 0.9062 0.09380.4472 0.2000 1.0000 0 最后两列为权重含约束条件的有效前沿[PortRisk, PortReturn, PortWts] = portopt (Return, Cov, [], PortReturn, ConSet) 约束条件ConSet = portcons ('Default', Num, 'AssetLims', Min, Max,'GroupLims', Group, GroupMin, GroupMax);>> r = [0.1 0.2 0.15]; >> s = [0.0100,-0.0061,0.0042; -0.0061, 0.0400,-0.0252; 0.0042,-0.0252,0.0225];>> ObjReturn = 0.17; >> [Risk, Return, Wts] = portopt(r, s, [], ObjReturn)Return =0.1700 Wts =0.0278 0.4278 0.5445>> Return = [0.1 0.2 0.15]; >> Cov = 0.01*[0.5,-1,0.4; -1,4, -0.2; 0.4,-0.2,2.3]; >> PortReturn = [0.15 0.16];目标收益率在0.15~0.16之间>> Num = 3; >> Min = [0.20 NaN NaN];投资组合中每种资产的下限>> Max = [0.5 0.5 0.5]; ……上限>> Group = [1 1 0; 0 1 1];分为两组,第一组由第1、2种资产构成;第二组……>> GroupMin = [0.2, 0.2]; 第一组的投资比例下限为0.2;二…… >> GroupMax = [0.8, 0.8]; 第一组的投资比例上限为0.8;二……>> ConSet = portcons('Default', Num, 'AssetLims', Min, Max,'GroupLims', Group, GroupMin, GroupMax);>> [PortRisk, PortReturn, PortWts] = portopt(Return, Cov, [], PortReturn, ConSet) PortRisk = 0.0724; 0.0910 PortReturn = 0.1500; 0.1600 PortWts = 0.4000 0.4000 0.20000.2606 0.4606 0.2789可借贷无风险资产时的有效前沿[RiskyRisk, RiskyReturn, RiskyWts, RiskyFraction 所有风险资产占总资产的比例, OverallRisk, OverallReturn] = portalloc(PortRisk, PortReturn, PortWts, RisklessRate, BorrowRate, RiskAversion 风险厌恶程度) >> Ret = [0.1 0.2 0.15];>> Cov = [ 0.005 -0.010 0.004; -0.010 0.040 -0.002; 0.004 -0.002 0.023]; >> [Risk, Return, Wts] = frontcon(Ret, Cov, 20);>> RiskfreeR = 0.08; BorrowR = 0.12; >> RiskAversion = 3;>> portalloc (Risk, Return, Wts, RiskfreeR, BorrowR, RiskAversion); (return the graphic) >>[RiskyRisk, RiskyReturn, RiskyWts, RiskyFraction, OverallRisk, OverallReturn] = portalloc (Risk, Return, Wts, RiskfreeR, BorrowR, RiskAversion) 最优化问题[x, fval] = fmincon(fun,目标函数, x0,从x0开始试值, A, b, Aeq, beq, lb, ub, nonlcon 非线性约束,options) options = optimset('LargeScale','off')注意:x1要加括号x(1);目标函数要建立m 文件;非线性约束要建立m 文件 程序 1: function [f] = Myobj(x) f=exp(x(1))*(4*x(1)^2+2*x(2)^2+4*x(1)*x(2)+x(2)+1);2: 写约束条件function [c, ceq] = Mycon(x)% Nonlinear inequality constraints c=[1.5+x(1)*x(2)-x(1)-x(2); -x(1)*x(2)-10];% Nonlinear equality constraints ceq = []; Step 3:调用最优化路径:A = [1 1; 2 -1]; b = [20; 10];Aeq = [];beq = []; lb = [-10 -10];ub = [10, 10];x0 = [-1,1]; % Make a starting guess at the solution options = optimset('LargeScale','off'); [x, fval] = fmincon(@Myobj, x0, A, b, Aeq, beq, lb, ub, @Mycon, options) B-S 模型 [Call, Put] = blsprice(Price, Strike, Rate, Time, V olatility, Yield) [CallDelta, PutDelta] = blsdelta(Price, Strike, Rate, Time, V olatility, Yield) [CallEl, PutEl] = blslambda,[CallRho, PutRho]= blsrho, V olatility = blsimpv(Price, Strike, Rate, Time, Value, Limit, Yield, Tolerance, Class) Chapter 3 二叉树模型 European and American Call European Put American Put 欧式二叉树 function [] = BinomialEuro() T,n,K,r,sigma,S0;deltat=T/n; u=exp(sigma*sqrt(deltat)); d=exp(-sigma*sqrt(deltat)); p=(exp(r*deltat)-d)/(u-d); q=1-p; for i=1:n+1 s(i)=u^(n+1-i)*d^(i-1)*S0; 股价 w(i)=nchoosek(n,i-1)*p^(n+1-i)*q^(i-1); x(i)=max(s(i)-K,0); y(i)=max(-s(i)+K,0); cv(i)=w(i)*x(i); pv(i)=w(i)*y(i); end C=sum(cv(:))*exp(-r*T) P=sum(pv(:))*exp(-r*T)美式二叉树function [c]=BinomialAmer() deltat=T/n; u=exp(sigma*sqrt(deltat)); d=exp(-sigma*sqrt(deltat)); p=(exp(r*deltat)-d)/(u-d); q=1-p; s(1,1)=S0; for t=2:n+1 (r 扣除红利) for i=1:ts(t,i)=u^(t-i)*d^(i-1)*S0; end; end for i=1:n+1; x(n+1,i)=max(K-s(n+1,i),0); y(n+1,i)=max(s(n+1,i)-K,0); end for t=n:(-1):1for i=1:tx(t,i)=(p*x(t+1,i)+q*x(t+1,i+1))/exp(r*deltat); x(t,i)=max(x(t,i),K-s(t,i));y(t,i)=(p*y(t+1,i)+q*y(t+1,i+1))/exp(r*deltat);y(t,i)=max(y(t,i),-K+s(t,i)); end; end c=y(1,1); p=x(1,1) Chapter 4 蒙特卡罗模拟 计算面积: y = cos x [-pi/2, pi/2] function [area] = MC1(n)X=rand(n,1)*pi-pi/2;Y=rand(n,1) ;Z=cos(X); counter=0; for i = 1 : n if Y(i) <= Z(i)counter=counter+1; end end area=pi*counter/n;估计pi function [Mypi] = Mypi() n = 10000;count = 0; for i = 1:n x = rand; y = rand; if x^2+y^2 <= 1count = count + 1; end end Mypi = 4*count/n; 计算积分 count = 0;for i = 1:n ;x = 2*rand+1; count = count + cos(x)*sin(x); end Myinte = (b-a)*count/n;欧式看涨期权定价Chapter 5 随机数生成器线性同余产生器 x(i+1)=ax(i)mod m; u(i+1)=x(i+1)/m function[]=line(a,m,n) x(1)=2 %余数mod(被除数,除数) for i=1:n x(i+1)=mod(a*x(i),m); u(i+1)=x(i+1)/m; end x u 混合线性同余产生器 x(i+1)=(ax(i)+b)mod m; u(i+1)=x(i+1)/m function [u] = lcg(seed, n) a = 1229;b = 1;m = 32768;x = zeros(1,n); x(1)=mod(a*seed+b, m); for i = 2:n x(i)= mod(a*x(i-1)+b, m); end u = x./m; % make the randomunmber between 0 and 1 plot(u(1:n-1), u(2:n), '.'); % plot the pair of point (u(n), u(n+1)) title('Autocorrelaton of LCG'); 计算循环的周期 function [p]=lcg(a,b,m,x0,n) x(1)=mod((a*x0+b),m)221122121212121212: min (4241)..: 1.50 100 20 210-10,10x x obj e x x x x x s t x x x x x x x x x x x x ++++--+≤--≤+≤-≤≤≤(),r q t t te d p u e d eu dσσ-∆∆-∆-⇒===-0,0,,0max{,0)T rTT i T i i f e p S K -==-∑0,0,,0max{,0)T rTT i T i i f e p K S -==-∑,01,1,1max{max{,0}, ()}t i i r tt i t i t i f K S u d e pf qf --∆+++=-+1()[()]()[()]()bbNi i aab a f x dx E f x dx b a E f x f x N =-≈=-≈1()[()]()[()]()b b Ni i a a b a f x dxE f x dx b a E f x f x N =-≈=-≈∑⎰⎰()()()21 1,, ()(0) ()ˆ ir T T Z i rT i i nnfor i n generate Z set S T S e set c e S T K set c c c nσ-++-===-=+u(1)=x(1)/mfor i=2:nx(i)=mod((a*x(i-1)+b),m)u(i)=x(i)/mif u(i)==u(1)p=iendend洗牌程序u = lcg(seed, n+20);Seq = ceil(20*lcg(10, n))+1;w = zeros(1,n);for i = 1:nh = Seq(i);w(i) = u(h);u(h) = u(i+20);endplot(w(1:n-1), w(2:n), '.'); % plot the pair of point (u(n), u(n+1)) title('Autocorrelaton of LCG');掷骰子程序function [outcome] = dice(n)U = rand(1, n);outcome = zeros(1,n);for i = 1:nif (0 <= U(i) & U(i) <= 1/6)outcome(i) = 1;elseif (1/6 < U(i) & U(i) <= 2/6)outcome(i) = 2;elseif (2/6 < U(i) & U(i) <= 3/6)outcome(i) = 3;elseif (3/6 < U(i) & U(i) <= 4/6)outcome(i) = 4;elseif (4/6 < U(i) & U(i) <= 5/6)outcome(i) = 5;elseoutcome(i) = 6;endend生成符合possion过程的随机数function [] = Poisson(seed, lambda, n)m = 1000;P = ones(1,m);for i = 1:mP(i) = exp(-lambda)*lambda^(i-1)/factorial(i-1);endpoisson = ones(1,n);rand('state',seed)for i = 1:ns1 = 0;s2 = P(1);U = rand;for j = 1:mif (s1 <= U & U <= s2)poisson(i) = j-1;breakelses1 = s2;s2 = s2 + P(j+1);endendendmean(poisson)var(poisson)生成符合连续概率分布的随机数>> rand(‘state',1);>> -2*log(rand)ans =1.34Box-Mueller方法生成符合正态分布的随机数function [w] = Boxmuller(seed, n)u = lcg(seed,n); G1 = zeros(1, n); G2 = zeros(1, n); w = zeros(1, 2*n);for i = 1: n,G1(i) = sqrt(-2 * log(u(i))) * cos(2 * pi * u(n-i+1));G2(i) = sqrt(-2 * log(u(i))) * sin(2 * pi * u(n-i+1)); endw = [G1 G2]; plot(G1, G2, '.'); title('Plot of G1, G2'); xlabel('G1');ylabel('G2'); [H,P,KSSTAT,CV] = kstest(w,[],0.05,0)figurenormplot(w);Chapter6 减小方差的方法控制变量法function []= ContrVariate()S0=90; K=90; r=0.05; sigma=0.3; T=0.25; n=10000;for i=1:n% simulate stock price using Continous risk neutral dynamics whendeltat=0.1deltaw = random('Normal',0,sqrt(T),1,1);ST(i) =S0* exp((r-1/2*sigma^2)*T+sigma*deltaw);Y1(i)=max(ST(i)-K,0)*exp(-r*T);endcovM=cov(ST,Y1); corrM=corrcoef(ST, Y1);b=-covM(1,2)/covM(1,1);for i=1:nY2(i)=Y1(i)+b*(ST(i)-S0*exp(r*T)); % payoff of European Callusing control variateendEurC1=mean(Y1(:)) EurC2=mean(Y2(:))VarEurC1=var(Y1); VarEurC2=var(Y2);ratio=VarEurC2/VarEurC1 % the ration of the variance of the optimallycontrolled estimator to the uncontrolled estimator is对偶变量法>> rand('seed',0) ;>> U=rand(1,100);>> X=exp(U);>> I=mean(X)I = 1.7461 >> S2=var(X) S2 =0.2499Chapter7 布朗运动function BrownianMotion()n=100000;deltat=1/n;mu=0.1; sigma=0.2; X=zeros(4,n);randnumber=normrnd(0,1,4,n);for i=1:nX(1,i+1)=X(1,i)+sqrt(deltat)*randnumber(1,i);X(2,i+1)=X(2,i)+mu*deltat+sigma*sqrt(deltat)*randnumber(2,i);X(3,i+1)=X(3,i)+mu*deltat+sigma*X(3,i)*sqrt(deltat)*randnumber(3,i);X(4,i+1)=X(4,i)+mu*X(3,i)*deltat+sigma*X(4,i)*sqrt(deltat)*randnumber(4,i); end跳跃过程function [c] = jumpdiffusion()n=10; r=0.1; sigma=0.4612; S0=90.44; lambda=1*4; deltat = 1/252;T=10/252;w = 1000;for s = 1:wX = random('Poisson',T*lambda);Y=1.1+rand(1,X)*(1.3-1.1);M=prod(Y);deltaw = random('Normal',0,sqrt(T),1,1);ST(s) = S0*exp((r-sigma^2/2)*T+sigma*deltaw)*M; endChapter8 美式期权的定价function AmericanPut()T=3; n=3; S0=1;K=1.1; r=0.06;sigma=0.3;deltat=T/n;lambda=0; w = 10000; State=[zeros(w,n),ones(w,1)];for s = 1:wSt(s,:) = zeros(1,n+1);St(s,1) = S0;for k = 1:nX = random('Poisson',deltat*lambda);Y=1.1+rand(1,X)*(1.3-1.1);M=prod(Y);deltaw = random('Normal',0,sqrt(deltat),1,1);St(s,k+1) = St(s,k)*exp((r-sigma^2/2)*deltat+sigma*deltaw); endVt(s,:) = zeros(1,n+1);V(s,n+1)=max(K-St(s,n+1),0);endfor k = n:-1:2xdata=[]; ydata=[]; for s = 1:wif St(s,k)<K xdata=[xdata,St(s,k)];ydata=[ydata,V(s,k+1)*exp(-r*deltat)]; end endx0 = [0, 0 , 0]; x = lsqcurvefit(@Regfun,x0,xdata,ydata);for s = 1:w if St(s,k)<KC(s,k)=x(1)+x(2)*St(s,k)+x(3)*St(s,k)^2;if K-St(s,k)>C(s,k);State(s,k)=1; V(s,k)=max(K-St(s,k),0);for s = 1:8for k = 1:n+1if State(s,k)==1;V0(s)=max((K-St(s,k)),0)*exp(-r*(k-1));break end end endmeanV0=mean(V0)function F = Regfun(x,xdata)F = x(1) + x(2)*xdata +x(3)*xdata.^2;Chapter9 有限差分的方法向前差分function imfdamput(Smax, dS, T, dT,K, R, SIG);M = ceil(Smax/dS);ds = Smax / M;N = ceil(T/dT);dt = T / N;J = 1:M-1;a = .5*R*dt*J - .5*SIG^2*dt*J.^2;b = 1 + SIG^2*dt*J.^2 + R*dt;c = -.5*R*dt*J - .5*SIG^2*dt*J.^2;A = diag(b) + diag(a(2:M-1), -1) +diag(c(1:M-2), 1);put = zeros(N+1, M+1);put(N+1, :) = max(K - [0:ds:Smax],0);%赋初值put(:, 1) = K; %赋初值put(:, M+1) = 0;%赋初值for i = N:-1:1y = put(i+1, 2:M)';y(1) = y(1) - a(1)*K;put(i, 2:M) = [A \ y]';put(i, :) = max(K - [0:ds:Smax],put(i,:));endx=0:dT:T;y=0:dS:Smax;mesh(y,x,put)向后差分程序function exfdamput(Smax, dS, T, dT,K, R, SIG);M = ceil(Smax/dS);ds = Smax / M;N = ceil(T/dT);dt = T / N;J = 1:M-1;a = (-.5*R*dt*J + .5*SIG^2*dt*J.^2)/ (1+R*dt);b = (1 - SIG^2*dt*J.^2) / (1+R*dt);c = (.5*R*dt*J + .5*SIG^2*dt*J.^2) /(1 + R*dt);A = diag(b) + diag(a(2:M-1), -1) +diag(c(1:M-2), 1);put = zeros(N+1, M+1);put(N+1, :) = max(K - [0:ds:Smax],0);put(:, 1) = K;put(:, M+1) = 0;for i = N:-1:1y = zeros(1, M-1);y(1) = a(1)*put(i+1, 1);y(M-1) = c(M-1)*put(i+1,M+1);put(i, 2:M) = put(i+1, 2:M) * A'+ y;put(i, :) = max(K - [0:ds:Smax],put(i,:));endx=0:dT:T;y=0:dS:Smax;mesh(y,x,put) *(,)()Cov X YbVar X=-()([])i i iY b Y b X E X=--。
matlab的打印文件
1、程序如下:function[y,y1,y2]=mwave(f1,m1,f2,m2)if(f1<2)|(f1>20) error('f1 is over!'), return,endif(f2<2)|(f2>20) error('f2 is over!'), return,endif(m1<0.5)|(m1>2) error('m1 is over!'), return,endif(m2<0.5)|(m2>2) error('m2 is over!'), return,endt=0:2*pi/(500-1):2*pi;y1=m1*sin(2*pi*f1*t);y2=m1*sin(2*pi*f2*t);y=y1+y2;figuresubplot(311);plot(t,y1);title('y1 waveform');subplot(312);plot(t,y2);title('y2 waveform');先以mwave.m形式保存起来,在command window 中输入mwave(5,1,10,2)回车键结束,结果如右图所示。
2、程序如下:a=1;s=0;for i=1:64s=s+aa=2*aendn=s/1.4/10^83、程序如下:for x=0:19for y=0:33for z=0:100if (x+y+z==100)&(5*x+3*y+z/3==100) d=[x,y,z]endendendend实验2:a=3的时候图形,见上图所示。
图二:a=7的时候图形,见上图所示。
实验6:1.a=rand(1,10)subplot(221); plot(a,'r');title('连线图');subplot(222); stem(a,'y');title('脉冲图');subplot(223); stairs(a,'b');title('阶梯图');subplot(224); bar(a,'g');title('条形图');2.x=-5:0.1:5;y1=2.^x;y2=(1/2).^x;plot(x,y1,'r');text(2,10,'y1=2^x'); hold;plot(x,y2,'b');text(-3,10,'y2=(1/2)^x');3.实验7:(1)>> A=[2,-3,7,-1]A =2 -3 7 -1>> AA=poly(A)AA =1 -5 -19 29 42(2)>>P=[1,-5,-30,150,273,-1365,-820,4100,576,-1880]; >>R = roots(P)R =4.9857-4.01074.0880-2.9157-2.18572.5561 + 0.2260i2.5561 - 0.2260i-0.77810.7043(3)>>P=[1,0,-2,1];R=roots(P)n=1;for x=-2:0.01:2y(n)=x^4-2*x^2+1;n=n+1;endx=-2:0.01:2;plot(x,y)结果:R =-1.61801.00000.6180(4)>>f1=[1,3,5,7]f1 =1 3 5 7>>f2=[8,-6,4,-2]f2 =8 -6 4 -2>>F=conv(f1,f2)F =8 18 26 36 -28 18 -14 >>f3=F-[0 0 0 f1]f3 =8 18 26 35 -31 13 -21 >>[q,r]=deconv(f3,f2)q =1.0000 3.0000 5.0000 6.8750r =0 0 0 0 -3.7500 -4.5000 -7.2500(5)>>Y='x^5+tan(4*x^2)+3'Y =x^5+tan(4*x^2)+3>> diff(Y)%求导ans =8*x*(tan(4*x^2)^2 + 1) + 5*x^4(6)>> f1=sym('x^3+3*x^2+5*x+7');>> f2=sym('8*x^3-6*x^2+4*x-2');>> f=f1*f2f =(x^3 + 3*x^2 + 5*x + 7)*(8*x^3 - 6*x^2 + 4*x - 2) >>collect(f)ans =8*x^6 + 18*x^5 + 26*x^4 + 36*x^3 - 28*x^2 + 18*x - 14>> (f-f1)/f2ans =-(5*x - (x^3 + 3*x^2 + 5*x + 7)*(8*x^3 - 6*x^2 + 4*x - 2) + 3*x^2 + x^3 + 7)/(8*x^3 - 6*x^2 + 4*x - 2) >>collect(ans)ans =(8*x^6 + 18*x^5 + 26*x^4 + 35*x^3 - 31*x^2 +13*x - 21)/(8*x^3 - 6*x^2 + 4*x - 2)。
第4章 MATLAB程序设计 [MATLAB大学教程][肖汉光,邹雪,宋涛]
case {2,3,4}
%价格大于等于200但小于500
rate=3/100;
case num2cell(5:9) %价格大于等于500但小于1000
rate=5/100;
case num2cell(10:24) %价格大于等于1000但小于2500
rate=8/100;
case num2cell(25:49) %价格大于等于2500但小于5000
7 16 27 40 55 72 ans = Error using ==> mtimes Inner matrix dimensions must agree.
命令文件可以直接运行,在MATLAB命令窗口输入命令 文件的名字,就会顺序执行命令文件中的命令,而函数 文件不能直接执行,而要以函数调用的方式来调用它。
4.1.2 M文件的建立与打开
M文件是一个文本文件,它可以用任何编辑程序来建立和编辑, 而一般常用且最为方便的是使用MATLAB提供的文本编辑器。 1.建立新的M文件
例4.6 矩阵乘法运算要求两矩阵的维数相容,否则会出错。 先求两矩阵的乘积,若出错,则自动转去求两矩阵的点乘。
程序如下:
A=[1,2,3;4,5,6]; B=[7,8,9;10,11,12];
try
C=A*B;
catch
C=A.*B;
end
C lasterr
%显示出错原因
4.2 结构化程序设计
>> fexch C=
1000≤price<2500 8%折扣
2500≤price<5000 10%折扣
5000≤price 14%折扣
输入所售商品的价格,求其实际销售价格。
程序如下:
matlab教程(全)09Matlab程序设计
2020/11/8
Application of Matlab Language
14
5.5 Matlab矩阵分析与处理
5.5.1 特殊矩阵 常见的特殊矩阵有零矩阵、幺矩阵、单位矩阵等,这类特殊矩阵在应用
中具有通用性。 1、通用的特殊矩阵 常用的产生通用殊矩阵的函数有: zeros:产生全0矩阵(零矩阵)。 ones: 产生全1矩阵(幺矩阵)。 eye: 产生单位矩阵。 rand:产生0~1间均匀分布的随机矩阵。 randn:产生均值为0,方差为1的标准正态分布随机矩阵。
5.4.1 程序调试概述 一般说来,应用程序的错误有两类,一类是语法错误,另一类是运行时
的错误。语法错误,给出相应的错误信息,并标出错误在程序中的行 号。例如:输入下列程序: A = 87;
B = 9.3;
C = A+*B; 系统将给出错误信息:
??? Error: File: Untitled1.m Line: 3 Column: 7
2020/11/8
Application of Matlab Language
4
说明:
将以上函数文件以文件名fcircle.m保存,然后在命令窗口调用。
[s,p] = fcircle(10) 输出结果是: s=
314.1593 p=
62.8319 采用help命令或lookfor命令可以显示出注释说明部分的内容。 help fcircle 屏幕显示
进行存取和修改。
全局变量用global命令定义,格式为:
global 变量名
例5.13 全局变量应用示例。
先建立函数文件wadd.m,该函数将输入的参数加权相加:
function f = wadd(x,y)
《MATLAB程序设计与应用》刘卫国高等教育出版社课后答案 已解锁-去水印 适合打印
SY404 clear all for n=1:4 if n==1 f1=1; elseif n==2 f2=0; elseif n==3 f3=1; else a=f3-2*f2+f1; b=a-2*f3+f2; c=b-2*a+f3; d=c-2*b+a; H=[1,0,1,a,b,c,d]; for m=8:4:99 a=d-2*c+b; b=a-2*d+c; c=b-2*a+d; d=c-2*b+a; H=[H,a,b,c,d]; end f100=d-2*c+b; end end max=max(H);
[f1,f2]=f(n); a=f1; b=f2; elseif n==30; [f1,f2]=f(n); c=f1; d=f2; else [f1,f2]=f(n); e=f1; f=f2; end end y1=e/(a+c); y2=f/(b+d); disp(['(1) y=',num2str(y1)]) disp(['(2) y=',num2str(y2)]) f function [f1,f2]=f(n) f1=n+10*log(n^2+5); x=0; for a=1:n b=a*(a+1); x=x+b; end f2=x; fushu function [e,l,s,c]=fushu(x) e=exp(x); l=log(x); s=sin(x); c=cos(x); disp(['复数e的指数是:',num2str(e)]) disp(['复数e的对数是:',num2str(l)]) disp(['复数e的正弦是:',num2str(s)]) disp(['复数e的余弦是:',num2str(c)])
实验二MATLAB程序设计含实验报告
实验二 MATLAB 程序设计一、 实验目的1.掌握利用if 语句实现选择结构的方法。
2.掌握利用switch 语句实现多分支选择结构的方法。
3.掌握利用for 语句实现循环结构的方法。
4.掌握利用while 语句实现循环结构的方法。
5.掌握MATLAB 函数的编写及调试方法。
二、 实验的设备及条件计算机一台(带有MATLAB7.0以上的软件环境)。
M 文件的编写:启动MATLAB 后,点击File|New|M-File ,启动MATLAB 的程序编辑及调试器(Editor/Debugger ),编辑以下程序,点击File|Save 保存程序,注意文件名最好用英文字符。
点击Debug|Run 运行程序,在命令窗口查看运行结果,程序如有错误则改正三、 实验内容1.编写求解方程02=++c bx ax 的根的函数(这个方程不一定为一元二次方程,因c b a 、、的不同取值而定),这里应根据c b a 、、的不同取值分别处理,有输入参数提示,当0~,0,0===c b a 时应提示“为恒不等式!”。
并输入几组典型值加以检验。
(提示:提示输入使用input 函数)2.输入一个百分制成绩,要求输出成绩等级A+、A 、B 、C 、D 、E 。
其中100分为A+,90分~99分为A ,80分~89分为B ,70分~79分为C ,60分~69分为D ,60分以下为E 。
要求:(1)用switch 语句实现。
(2)输入百分制成绩后要判断该成绩的合理性,对不合理的成绩应输出出错信息。
(提示:注意单元矩阵的用法)3.数论中一个有趣的题目:任意一个正整数,若为偶数,则用2除之,若为奇数,则与3相乘再加上1。
重复此过程,最终得到的结果为1。
如:2?13?10?5?16?8?4?2?16?3?10?5?16?8?4?2?1运行下面的程序,按程序提示输入n=1,2,3,5,7等数来验证这一结论。
请为关键的Matlab 语句填写上相关注释,说明其含义或功能。
MATLAB语言程序设计基础
3.2 matlab语言基本运算及输入输出
3.2.5 输入与输出语句
input A=input(提示字符串)要求输
入矩阵
A=input(提示字符串,‘s’) 要求字符串eg:
n=input('how much')
n=input('ho第w23页m/共4u1页ch','s')
3.2 matlab语言基本运算及输入输出
关系运算和逻辑运
算
表3-6 关系运算和逻辑运算函数
函数 any all find
exist isnan
意义 逻辑条件任何一个
逻辑条件全部 寻找逻辑值的向量元素下 标
检查某变量是否存在 检查非数值量
函数 finite isempty isstr
strcmp
随机数元素矩阵 设三维绘图基底坐
第5页/共41页
单位矩阵
3.1.3构造多维数组
cat( ) a=cat(n,a1,a2,….) n:多维函数的维数 n=1:
cat(a1,a2,a3..)=[a1;a2;a3…] n=2:
cat(a1,a2,a3..)=[a1,a2,a3…] n=3: 图3-1示
重新定义维数
end 表示某一维末尾元素下标
2、复数矩阵: b=[1 2;3 4]+i*[5 6;7 8] b=[1+5i 2+6i;3+7i 4+8i]
第3页/共41页
3、空矩阵[]
0×0阶
与clear不同之处:clear删除变量
[]删除矩阵中的元素
A(:,[2,3])=[] 第2,3列元素删除
函数 abs angle sqrt real imag conj round fix
matlab程序设计基础-PPT
在 MATLAB 中,除了可以在命令窗口中输入命 令逐句执行外,也可以和其他形式的 C、FORTRAN 等高级语言一样采用编程的方式,这就是 M 文件编 程。
MATLAB 程序设计原则 ➢ 百分号“%”后面的内容是程序的注解,要善于运用
注解使程序更具可读性; ➢ 养成在主程序开头用 clear 指令清除变量的习惯,以消
else leap=1;}
else leap=0;
if(leap) printf(“%d is”,year);
else printf(“%d is not”,year);
printf(“a leap year.\n”);}
MATLAB程序: year=input('year='); if rem(year,4)==0
end
if 表达式 A 语句组1
elseif 表达式B 语句组2
else 语句组3
end
注意:除直接应用上述三种形式外,第3种结构可扩展,if 还可以嵌套。
举例:
例5.1 输入数n,判断其奇偶性。
程序式书写法:* n=input(‘n=’); if rem(n,2)==0 A=‘even’ else A=‘odd’ end
MATLAB程序:
t=1;pi=0; n=1;s=1; while abs(t)>1e-6
pi=pi+t; n=n+2; s=-s; t=s/n; end pi=4*pi
6、MATLAB 中的函数及调用
MATLAB 函数
匿名函数可匿以名每接函一受数个多实M个文例输件:入第输一出行参定数义。的创文建件匿就名是函M文数的格式: fhandle=>@>(件amr主gyflih函sdt)数1e=,x@p一r(x个)(xM+x文.^件2)只能包含一个主函 其中:“ex>p>r数”m通,y常fh通d是常1一(2将)个M简文单件的名M和AMTL文A件B变主量函表数达名式设,实现函数 的功能;“aanrs为g=li一6st”致是。参数列表;“@”是MATLAB中创建函数句柄
《MATLAB程序设计》实验指导书
三、
实验仪器和设备
1、 计算机一台。 2、 MATLAB7.0 以上集成环境。
四、
预习要求
2
《MATLAB 程序设计》实验指导书
1、 复习 MATLAB 的启动与Байду номын сангаас出,熟悉 MATLAB 运行环境。 2、 复习 MATLAB 中矩阵的生成以及矩阵运算的基本原理。
五、
实验内容及步骤
实验内容:
1、 求下列表达式的值 1) z1
5
《MATLAB 程序设计》实验指导书
2) 函数调用 函数文件编制好后,就可调用函数进行计算了。函数调用的一般格式为 [输出实参表]=函数名(输入实参表) 注意:函数调用时各实参出现的顺序、个数,应与函数定义时形参的顺序、个 数一致,否则会出错。函数调用时,实参先传递给形参,然后再执行函数功能。 5、 选择结构 1) if 语句 a) 单分支 if 语句 if 条件 语句块 end b) 双分支 if 语句 if 条件 语句块 1 else 语句块 2 end c) 多分支 if 语句 if 条件 1 语句块 1 elseif 条件 2 语句块 2 …… elseif 条件 n 语句块 n else 语句块 n+1 end 2) switch 语句 switch 表达式 case 结果表 1 语句块 1 case 结果表 2 语句块 2 ……
二、
实验原理
1、 M 文件 用 MATLAB 语言编写的程序,称为 M 文件,它们的扩展名均为.m。M 文件根据 调用方式的不同分为两类,命令文件(Script file)和函数文件(Function file) 。 2、 建立新的 M 文件 启动 MATLAB 文本编辑器有 3 种方法: 1) 单击工具栏上的“New M-File”命令按钮。 2) 从 MATLAB 主窗口的“File”菜单中选择“New”菜单项,再选择“M-file”命令。 3) 在 MATLAB 命令窗口输入命令“edit”。 3、 打开已有的 M 文件 1) 在当前目录窗口选中要打开的 M 文件,双击鼠标左键。 2) 单击 MATLAB 主窗口工具栏上的“Open File”命令按钮,再从弹出的对话框中 选择所需打开的 M 文件。 3) 从 MATLAB 主窗口的“File”菜单中选择“Open”命令,在“Open”对话框中选中 所需打开的文件。 4) 在 MATLAB 命令窗口输入命令“edit 文件名” 。 4、 函数文件 1) 函数文件的基本结构 函数文件由 function 语句引导,其基本结构为: function 输出形参表=函数名(输入形参表) 注释说明部分 函数体语句 注意:函数名的命名规则与变量名相同。当输出形参多于一个时,应用方括号 括起来。
matlab打印函数
matlab打印函数MATLAB印函数是MATLAB编程中常用的工具,它可以把编辑的程序代码打印出来,并轻松分享给别人。
MATLAB的打印函数包括print 和fprintf,它们都有自己独特的功能,因此通常会被编程人员用来满足多种不同的应用需求。
本文将详细分析MATLAB中print函数和fprintf函数,具体涵盖以下内容:一、Print函数1.基本概念Print函数是MATLAB中非常常用的一种打印函数,它可以把编辑的程序代码打印出来,并使用标准文件输出。
它包括三种参数:第一个参数是文件名,第二个参数是格式,第三个参数是换行参数。
它的使用方法是,在程序中使用语句“print(文件名,格式,换行参数)”进行操作,例如:a=2;b=3;print(test.txt%f %fa,b);以上代码将以“test.txt”的文件名打印出a和b的数值,a和b以格式“%f %f”输出,最后一个参数“a”代表每输出一行,即使用换行符“a”分隔。
2.应用实例(1)实际应用中,print函数可以用来将程序输出的数值存储在txt文件中。
例如,在MATLAB中生成一定的随机数序列后,可以将这些数序列保存在txt文件中,从而方便下次使用。
x = rand(10,1);print(‘rand_array.txt’,’%f‘,x);以上代码将以“rand_array.txt”的文件名,存储生成的随机数序列,格式为“%f”,最后一个参数“x”代表每输出一行,即使用换行符“x”分隔。
(2)时,print函数还可以被用来在不同的程序中传递变量,例如:a=2;b=3;print(data.txt%f %fa,b);data = load(data.txtc=data(1);d=data(2);以上程序将以“data.txt”的文件名,存储变量a和b的值,然后以load函数读取存储的变量值,并将变量a和b的值存入新变量c和d中。
二、fprintf函数1.基本概念fprintf函数是MATLAB中一种非常常用的打印函数,它可以用来把格式化的文本输出到文件或者命令行窗口等。
MATLAB程序设计范文精简版
MATLAB程序设计范文精简版MATLAB程序设计MATLAB程序设计简介基本语法MATLAB的基本语法与其他编程语言类似,包括变量定义、运算符、控制流程等。
以下是一些常用的基本语法:变量定义MATLAB中的变量不需要预先声明类型,直接使用即可。
变量名是大小写敏感的,并且不能使用保留字作为变量名。
matlabx = 5;y = 'Hello MATLAB!';z = [1 2 3 4 5];运算符与其他编程语言一样,MATLAB支持各种数学运算符和逻辑运算符,可以进行加减乘除、比较和逻辑操作等。
matlabMATLAB程序设计范文精简版a = 5 + 3;b = 7 2;c = (a > b) && (b < 10);控制流程MATLAB提供了各种控制流程语句,如条件语句、循环语句等,可以根据条件执行不同的操作。
matlabif x > 0disp('x is positive');elseif x < 0disp('x is negative');elsedisp('x is zero');endfor i = 1:5disp(i);endwhile x < 10x = x + 1;end函数定义和调用函数是MATLAB程序设计的重要组成部分,可以封装一些常用的操作和算法,并在需要时调用。
以下是函数的定义和调用示例:matlabfunction result = add(a, b)result = a + b;endx = 3;y = 4;z = add(x, y);disp(z);数据处理和可视化MATLAB提供了丰富的数据处理和可视化工具,可以帮助用户对数据进行分析和展示。
以下是一些常用的数据处理和可视化操作示例:加载和保存数据matlabdata = load('data.txt');save('result.mat', 'data');统计分析MATLAB提供了丰富的统计函数,可以进行各种统计分析操作,如求平均值、标准差、相关系数等。
matlab程序画几何图形并进行网格划分(可打印修改)
编写程序,对下图几何区域按有限元法进行网格划分,确定单元中结点数目及位置,并对单元和结点分别进行编号。
二、程序思路基于MATLAB软件平台,首先画出几何图形,然后对几何区域按照有限元法思想进行三角形网格剖分,最后按顺序依次输出单元及结点编号和坐标。
三、编写源程序:clear;%第一步:画出几何区域%为按逆时针顺序绘制几何区域x_up=[3.5,0];y_up=[2,2]; %上固体边界数组x_l=[0,0];y_l=[2,0]; %入流左边界数组x_low=[0,2.5];y_low=[0,0]; %下固体边界数组n=9; %圆弧等分变量nc=linspace(pi,pi/2,n);r=1; %圆弧半径变量rx_c=r*cos(c)+3.5;y_c=r*sin(c); %三角函数表示圆柱左上部分x_r=[3.5,3.5];y_r=[1,2]; %出流右边界数组x=[x_low,x_c ,x_r,x_up,x_l]; %整体边界X坐标数组y=[y_low,y_c ,y_r,y_up,y_l]; %整体边界Y坐标数组plot(x,y,'m'); %绘制几何区域xlabel('x轴');ylabel('y轴');title('几何区域');hold onl=k*(pi/2)*r;x_up(3)=3.5-l;y_up(3)=2;cut_point=[x_up(3), y_up(3)];x_up(3)=x_up(2);y_up(3)=y_up(2);x_up(2)= cut_point(1); y_up(2)= cut_point(2); %将上边界点坐标按顺序重排x_cutline=[ x_up(2),x_low(2)];y_cutline=[ y_up(2),y_low(2)]; %分割线数组plot(x_cutline, y_cutline); %绘制分割线n1=6; %上边界左半部分n1等分变量x_upcp1=linspace(x_up(2),x_up(3),n1+1); %上边界左半部分n1等分y_upcp1=linspace(y_up(2),y_up(3),n1+1);x_lowcp1=linspace(x_low(2),x_low(1),n1+1); %下边界n1等分y_lowcp1=linspace(y_low(2),y_low(1),n1+1);x_upcp2=linspace(x_up(1),x_up(2),n+1); %上边界右半部分n等分y_upcp2=linspace(y_up(1),y_up(2),n+1);x_upper=[ x_upcp2,x_upcp1]; %上固体边界坐标数组y_upper=[ y_upcp2,y_upcp1];c1=linspace(pi*3/2,pi,n+1);x_c1=cos(c1)+3.5;y_c1=abs(sin(c1)) ;x_lower=[x_c1,x_lowcp1];y_lower=[y_c1,y_lowcp1];for i=1:n1+n+2 %画竖线plot([x_upper(i) x_lower(i)],[ y_upper(i) y_lower(i)],'m');endm=7;for i=n+n1+2:-1:1plot([x_upper(i) x_lower(i)],[ y_upper(i) y_lower(i)],'m');%绘制竖线x1(i,:)=linspace(x_upper(i),x_lower(i),m+1); %竖线m等分y1(i,:)=linspace(y_upper(i),y_lower(i),m+1);endfor i=1:mfor j=n+n1+2:-1:2 plot([x1(j,i) x1(j-1,i)],[y1(j,i) y1(j-1,i)],'m');%plot([x1(j,i) x1(j-1,i+1)],[y1(j,i) y1(j-1,i+1)],'m');endendfor j=1:n+n1for i=1:mp1=i+(j-1)*(m+1);p2=i+1+j*(m+1);p3=i+j*(m+1);t=2*i+2*(j-1)*m;p4=i+(j-1)*(m+1);p5=i+1+(j-1)*(m+1);p6=i+1+j*(m+1);e(t-1,:)=[p1,p2,p3]; %整体编号与局部编号关系e(t,:)=[p4,p5,p6];endend四、程序各步运行结果:具体程序如上所写第一步:画出几何区域。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
实验四 Matlab程序设计[实验目的]1.掌握字符串数组的创建和构造方法及常用字符串函数的使用。
2.熟练掌握MATLAB 控制流的使用方法。
3.熟悉M 脚本文件、函数文件的编写方法和技巧。
[实验原理]与数值数组相比,串数组在MATLAB 中的重要性较小,但不可缺少。
如果没有串数组及相应的操作,那么数据可视化、图形用户界面的制作将会遇到困难。
字符串与数值数组是两种不同的数据类,它们的创建方式也不同。
字符串的创建方式是:将待建的字符放在“单引号对”中。
注意,“单引号对”必须是在英文状态下输入,其作用是MATLAB 识别送来内容“身份”所必需的,如A=’This is an example!’;就创建了一个字符串A。
注意创建带单引号的字符串时,每个单引号符用“连续2 个单引号符”标识。
字符串的标识同数值数组同,而且也可以使用size 指令观察串数组的大小。
串数组的ASCII 码可以通过指令abs 和double 来获取,而用char 指令可以把ASCII 码变为串数组,另外,MATLAB 可以很好的支持中文字符串数组。
对于复杂串数组的创建,一是可以直接创建,但是要保证同一串数组的各行字符数相等,即保证各行等长,不推荐,太繁琐。
二是可以利用串操作函数创建多行数组,比如char, str2mat, strvcat 等,具体操作自己通过帮助体会。
另外还可以通过转化函数产生数码字符长,比如A_str=int2str(A) 就是把整数数组A 转换成串数组,如果是非整数将被四舍五入后再转换,类似的函数还有num2str(把非整数数组转换为串数组,常用于图形中数据点的标识)、mat2str(把数值数组转换成输入形态的串数组,常与eval 指令配用)。
假如想灵活运用MATLAB 去解决实际问题,想充分调动MATLAB——科学技术资源,想理解MATLAB 版本升级所依仗的基础,那么掌握M 脚本文件合函数的编写规则将十分有用。
用户通过本次实验,感受抽象概念的内涵、各指令间的协调,从感知上领悟MATLAB 编程的优越和要领。
编写M 脚本文件的步骤:●点击MATLAB 指令窗工具条上的New File图标,就可打开如上图所示的MATLAB 文件编辑调试器MATLAB Editor/Debugger。
其窗口名为untitled ,用户即可在空白窗口中编写程序。
●点击编辑调试器工具条图标,在弹出的Windows 标准风格的“保存为”对话框中,选择保存文件夹,键入新编文件名(如newfile.m),点动【保存】键,就完成了文件保存。
●运行可有两种方法,一种是直接点击编辑调试工具条图标,即可直接运行;或者使newfile.m 所在目录成为当前目录,或让该目录处在MATLAB 的搜索路径上,然后在命令窗口键入指令newfile+回车,便可得到运行结果。
●调试程序方法有多种,常见的是设置断点的方法,将光标移到程序欲执行到的位置,点击编辑调试工具条图标,保存后运行,程序将停止在该语句位置并弹出编辑器界面等待用户下一步运行的指令,只有再次点击按钮,才继续向下执行。
相应的按下按钮,表示清除所有断点。
如果不设置断点,也可以在程序中加入pause 指令,使得程序在此处暂停,只有用户按任意键程序才依次向下执行。
则在pause 指令的前面位置我们可以通过交互的方式得到我们想要的信息,以检测程序的正确性。
编写MATLAB 脚本文件或函数文件时要区分开与C 语言格式的不同。
MATALB 使用变量前不需要声明数据类型,对于所有的数值型数据MATLAB 均以Double 型存储。
另外编程时尽量使用MATLAB 向量(数组)编程方式,可大大提高编程效率,尽量避免过多使用for 循环等语句。
MATLAB 提供了五种控制流的结构:for 循环结构,while 循环结构,if-else-end 分支结构,以及switch-case 结果、try-catch 结构。
这些控制指令用法与其他语言十分类似,这里只给出简要说明。
For 循环:while 循环结构if-else-end 结构单分支(常用)双分支(常用)多分支(常被swith-case 取代)上面几条控制语句中,for 循环结构中x 称为循环变量,组命令(commands)被称为循环体,循环体被重复执行的次数是确定的,该次数由for 指令后面的数组array 的列数决定。
换言之,循环变量依次取数组的各列,对于每个变量值,循环体被执行一次。
while 循环是首先检测expression 的值,如其值为逻辑真(非0),则执行组命令,当组命令执行完毕,继续检测表达式的值,仍为真,循环执行组命令,一旦表达式值为假,就结束循环。
一般情况下,表达式的值是标量值,但MATLAB 允许其为一个数组,此时只有该数组所有元素均为真时,MATLAB 才会执行循环体。
若表达式为空数组,则不执行循环体。
if 指令判决和break 指令的配合使用,可以强制中止for 循环或while 循环。
switch-case 结构 try-catch 结构switch 指令后面的表达式应为一个标量或者为一个字符串。
对于标量形式的表达式, 比较这样进行:表达式==检测值i 。
对于字符串,MATLAB 将调用函数strcmp 来实现比较: strcmp(表达式,检测值i).try-catch 结构,只有当MATLAB 在执行组命令1 时出现错误后,组命令2 才会被执行。
当执行组命令2 时又出错,MATLAB 将中止该结构。
随指令数的增加或随控制流复杂度的增加,以及重复计算要求的提出,采用M 脚本文件 进行编程较为适宜。
这种脚本文件的构成比较简单,它是一串按照用户意图排列而成的 MATLAB 指令集合。
脚本文件运行后,所产生的所有变量都驻留在MATLAB 基本工作空间中, 只要用户不使用clear 指令加以清除,且MATLAB 指令窗口不关闭,这些变量将一直保存在 基本工作空间中。
与脚本文件不同的,函数文件犹如一个“黑箱”。
从外界只能看到传给它的输入量和送 出的计算结果,而内部运作是藏而不见的,特点是:● 从形式上看,与脚本文件不同,函数文件的第一行总是以“function ”引导的“函数声 明行”。
该行还罗列出函数与外界的联系的全部“标称”输入输出宗量。
但对“输入输出宗量”的标称数目并没有限制,即可以完全没有输入输出宗量,也可以是任意数目。
形如unction sa=circle(r,s)。
这里r 、s 称为输入宗量,sa 称为输出宗量,函数名为circle ,同时注意保存的函数文件名应与这里的函数名一致,即存为circle.m 文件。
● MATLAB 允许使用比“标称”数目较少的输入输出宗量实现对函数的调用,但前提是函数中应该有相应的处理程序。
● 从运行上看,与脚本文件不同,每当函数文件运行时,MATLAB 就会专门为它开辟一个临时的工作空间,称之为函数工作空间。
所有中间变量都存放在函数工作空间中。
当执行完文件最后一条指令或遇到return 时,就结束该函数文件的运行,同时该临时函数空间及其所有的中间变量就立即被清除。
●假如在函数文件中,发生对某脚本文件的调用,那么该脚本文件运行产生的所有变量都存放在该函数空间中,而不是存放在基本空间。
[实验内容]一.串数组的创建和寻访1.先请实际操作下例,以体会数值量与字符串的区别clear %清除所有内存变量a=12345.6789 %给变量a 赋数值标量 class(a) %对变量a 的类别进行判断 a_s=size(a) %数值数组a 的“大小” b='S' %给变量b 赋字符标量(即单个字符) class(b) %对变量b 的类别进行判断 b_s=size(b) %符号数组b 的“大小 whos %观察变量a,b 在内存中所占字节2.已知串数组a=”This is an example.”,试将其到序输出。
3.接上题,试执行ascii_a=double(a),观察其ASCII 码,并将ASCII 码变回字符串。
4.设A=”这是一个算例”,重复上面的2-3。
5.尝试用直接输入法在命令窗口创建字符串s ,第一行时“This string array ”,第二行是“has multiple rows.”。
6.利用串操作函数char 、str2mat 、strvcat 分别写出使以下这段文字成为字符串的程序,注意保持这段文字的格式。
在英式用法中,引号通常是单引号,如‘Fire!’。
In GB usage quotation marks are usually single:’Fire!’.二.脚本文件实现)3cos(14.0t e y t --=,0 ≤t ≤3π,并在图上标出图名和极大值点坐标。
如下图所示。
可能用到的函数:num2str, char, text, hold on, 具体应用自己查找help 文档。
三.编程实现分别用for 或while 循环语句计算:6326322212+⋅⋅⋅+++==∑=i i K的程序,并给出运行结果。
此外,实现一种避免使用循环的的计算程序。
四.函数文件1.详读并运行下面的circle.m 函数文件。
体会M 函数文件的编写结构及方法。
%后面的内容称为注释行,不被执行,起注释说明作用。
2.编写一个简单的函数文件,它具有如下性质:该函数被调用时,如果不指定输入变量,则自动输出“用户,你忘记给定输入变量了!”;当输入大于1 的整数时,则输出“你是一个合法用户!”;当输入的是一个非正整数时,函数文件会给出一个错误提示“你是非法用户!”【提示:可能用到disp,error 等指令,使用方法自己查询帮助】实验五 Matlab符号计算[实验目的]1.掌握定义符号对象的方法;2.掌握符号表达式的运算法则以及符号矩阵运算。
3.掌握求符号函数极限及导数的方法。
4.掌握求符号函数定积分和不定积分的方法。
[实验原理]1.利用符号法求表达式的值已知x=6,y=5,利用符号表达式求z=>> syms x y>> x=6;>> y=5;>> z=(x+1)./(sqrt(x+3)-sqrt(y)) z =9.163118960624632. 分解因式x y-44>> clear>> syms x y>> factor(x.^4-y.^4) ans =(x-y)*(x+y)*(x^2+y^2)3.化简表达式sin cos cos sin ββββ-1212>> syms bata1 bata2>> simplify(sin(bata1).*cos(bata2)-cos(bata1).*sin(bata2))4.求sin tan()()limsinx xxx e ex→+--3121。