Matlab 自适应 FDM Burgers

来源:互联网 发布:数控机床编程入门书籍 编辑:程序博客网 时间:2024/06/05 15:11
function sec1_2_burgersFDM(jmax)%% example driver for Burgers' equation in section 1.2.% it calls movfdm().%%% Copyright (C) 2010 Weizhang Huang and Robert D. Russell% all rights reserved.%% This program is provided "as is", without warranty of any kind.% Permission is hereby granted, free of charge, to use this program% for personal, research, and education purposes. Distribution or use% of this program for any commercial purpose is permissible% only by direct arrangement with the copyright owner.%cpu0=clock;   % job = 1 for solution   % job = 2 for mesh trajectories   job = 1;   %jmax = 41;   npde = 1;   nn = npde*jmax;   x=zeros(jmax,1);   u=zeros(npde,jmax);% define initial solution   x=linspace(0,1,jmax)';   u(1,:)=(sin(2*pi*x)+0.5*sin(pi*x))';   xt=zeros(jmax,1);   ut=zeros(npde,jmax);% call the moving mesh function% monitor = 0 (fixed mesh), 1, 2, 3, 4, or 5% mmpde = 4, 5, 6, 7 (mmpde = 7 ==> modified MMPDE5)% alpha_type = 1, 2, or 3: (integral def, integral def with flooring, or alpha = constant)   monitor = 3;   mmpde = 7;   alpha_type = 2;   alpha = 1.0;   reltol=1e-6; abstol=1e-4;   if job==1 % for solution      tspan = [0 0.2 0.4 0.6 0.8 1.0];      [t,X,U]=movfdm(npde,jmax,tspan,x,xt,u,ut,@PDE_F,@PDE_G,@BC_L,@BC_R,...                   [],monitor,reltol,abstol,[],mmpde,alpha_type,alpha);      fprintf('\n cpu time used = %e \n', etime(clock,cpu0));      % output the solution      N=size(t,1);      figure      labels={};      marks='+osdv^<>';      lines1='--:-';      lines2='-  .';      colors='mkbgrc';      for n=1:N         u=U(:,:,n);         x=X(:,n);         hold on         labels{n} = ['t = ' num2str(t(n)) ];         style=['-' marks(1+rem(n,8))];         plot(x,u(1,:)',style);%,'LineWidth',2);      end      hold off      legend(labels{:});      xlabel('x');      ylabel('u');      axis([0 1 -1 2]);      box on;   else % for trajectories      tspan = [0 1.0];      [t,X,U]=movfdm(npde,jmax,tspan,x,xt,u,ut,@PDE_F,@PDE_G,@BC_L,@BC_R,...                   [],monitor,reltol,abstol,[],mmpde,alpha_type,alpha);      fprintf('\n cpu time used = %e \n', etime(clock,cpu0));      plot(X',t);      axis([0 1 0 1]);      xlabel('x')      ylabel('t')   endfprintf('cpu time used = %e \n', etime(clock,cpu0));% -----------------------------------------------------function f=PDE_F(t,x,u,ux,ut)  f(1,:) = ut(1,:);function f=PDE_G(t,x,u,ux,ut)  f(1,:) = 1e-4*ux(1,:)-u(1,:).*u(1,:)*0.5;function f=BC_L(t,x,u,ux,ut)  f(1) = ut(1)-0.0;function f=BC_R(t,x,u,ux,ut)  f(1) = ut(1)-0.0;
原创粉丝点击