数学是知识的工具,亦是其它知识工具的泉源。 ——勒内·笛卡尔 数 学 建 模 数学是知识的工具,亦是其它知识工具的泉源。 ——勒内·笛卡尔
第9章 模型求解常用软件 宋显华 songxianhua@hrbust.edu.cn
内容 1 常用软件介绍 2 MATLAB基本操作 3 MATLAB编程及调试 4 MATLAB求解算例 数学建模
第9章 模型求解常用软件 第42讲 MATLAB求解算例
算例1 研究偶极子(Dipole)的电势(Electric potential)和电场强度(Electric field density)。 设在 处有电荷 ,在 处有电荷 。 那么在电荷所在平面上任何一点的电势和场强分别为 , 其中, 又设电荷 数学建模
算例1 clear;clf; q=2e-6;k=9e9;a=1.5;b=-1.5; x=-6:0.6:6;y=x; [X,Y]=meshgrid(x,y); rp=sqrt((X-a).^2+(Y-b).^2);rm=sqrt((X+a).^2+(Y+b).^2); V=q*k*(1./rp-1./rm); [Ex,Ey]=gradient(-V); AE=sqrt(Ex.^2+Ey.^2); Ex=Ex./AE;Ey=Ey./AE; cv=linspace(min(min(V)),max(max(V)),49); 数学建模
算例1 contourf(X,Y,V,cv,'k-') %axis('square') title('\fontname{隶书}\fontsize{22}偶极子的场'),hold on quiver(X,Y,Ex,Ey,0.7) plot(a,b,'wo',a,b,'w+') plot(-a,-b,'wo',-a,-b,'w-') xlabel('x');ylabel('y'),hold off 数学建模
算例1 数学建模
算例2 hist指令的使用示例。 randn('state',1),rand('state',31) x=randn(100,1);y=rand(100,1); subplot(1,2,1),hist(x,7) subplot(1,2,2),histfit(x,20) 数学建模
算例2 n_y1=min(y):0.1:max(y);n_y2=min(y):0.05:max(y); subplot(1,2,1),hist(y,n_y1) subplot(1,2,2),hist(y,n_y2) 数学建模
算例3 取发生信号的原始模型为 。 在[0, 4]中取值; 受到噪声0.3*(rand(n,1)-0.5) 的污染。 本例演示:如何以 为模型,通过fminsearch从受污染数据中,估计出参数 。 数学建模
算例3 (1)创建M文件xydata.m [XYDATA.M] function [x,y,STDY]=xydata(k_noise) yo=3*exp(-0.4*x)+12*exp(-3.2*x); rand('seed',234) y_noise=k_noise*(rand(size(x))-0.5); y=yo+y_noise; STDY=std(y_noise); 数学建模
算例3 创建M文件twoexps.m [TWOEXPS.M] function E=twoexps(a,x,y) %twoexps.m x=x(:); y=y(:); Y=a(1)*exp(-a(3)*x)+a(2)*exp(-a(4)*x); E=sum((y-Y).^2); 数学建模
算例3 (2)编写M脚本文件作为主文件 [EXM.M] %exm041022_1.m k_noise=0.3; [x,y,STDY]=xydata(k_noise); a0=[1 1 1 1]; options=optimset('fminsearch'); options.TolX=0.01; options.Display='off'; a=fminsearch(@twoexps,a0,options,x,y); chi_est=twoexps(a,x,y)/STDY^2; freedom=length(x)-length(a0); Q=1-chi2cdf(chi_est,freedom) 数学建模
算例3 y_est=a(1)*exp(-a(3)*x)+a(2)*exp(-a(4)*x); ych='y_e_s_t='; a1=num2str(a(1)); a2=num2str(a(2)); a3=num2str(a(3)); a4=num2str(a(4)); char_y_est=[ych,a1,'*exp(-',a3,'*x)+',a2,'*exp(-',a4,'*x)']; plot(x,y,'b+');hold on, plot(x,y_est,'r'); hold off, axis([0,4,0,16]) 数学建模
算例3 text(0.4,14,'y=3*exp(-0.4*x)+12*exp(-3.2*x)') text(0.4,12,char_y_est),text(2.5,9,['chi2=' , num2str(chi_est)]) text(2.5,7,['freedom=' , num2str(freedom)]) text(2.5,5,['Q=' , num2str(Q)]) 数学建模
算例3 (3) exm 数学建模
作业