| 查看: 2369 | 回复: 10 | ||
chen17307新虫 (初入文坛)
|
[求助]
求助各路大神一个偏微分方程的数值计算问题（FPK方程）
|
|
求助各位虫友大神，用MATLAB写一个求解偏微分方程的程序，方程和离散化的推导过程如图1和图2，程序贴在下面了，能跑出结果但是结果不对，与精确解的结果不相符合啊。不知道哪里错了。跪求各路大神帮忙看看，小弟不胜感激啊。。。 %隐式差分格式求解FPK方程 %二阶差分格式 clear;clc; %% 参数设置 Bata=0.1;D=0.01;delta=4;omega0=1;%参数设置 min_x=-1.5; max_x=1.5; min_y=-1.5; max_y=1.5; %坐标轴范围 Nx=300; Ny=300; %迭代步数 sx=(max_x-min_x)/Nx; sy=(max_y-min_y)/Ny; %迭代步长 x=min_x:sx:max_x; y=min_y:sy:max_y; %网格划分 for i=1:Nx+1 for j=1:Nx+1 h(i,j)=Bata*y(j)-omega0*x(i)+delta*(x(i))^3; %化简式中的h(x,y) end end %% 求精确解 %其精确解为p(x,y)=c*exp(-Bata/D*(y^2/2-omega0^2/2*x^2+delta/4*x^4)) %---精确解---% for k=1:length(x) for l=1:length(y) p(k,l)=exp(-Bata/D*(y(l)^2/2-omega0/2*x(k)^2+delta/4*x(k)^4)); end end %精确解的图像 figure(1) mesh(x,y,p) title('精确解') xlabel('x','FontSize',12) ylabel('y','FontSize',12) zlabel('p(x,y)','FontSize',12) axis([-1.5,1.5,-1.5,1.5,0,2]) %% 求数值解 %---赋初值和边值---% %---初值和边界值来自求得的精确解的结果 p1=zeros(Nx+1); p1(:,1)=p(:,1);%左边值 p1(:,end)=p(:,end);%右边值 p1(1,  =p(1,;%初值for ii=1:Nx %---线性方程组AX=B的系数矩阵A---% %A为三对角阵，主对角：Diag；上对角：Updiag %下对角：Lowdiag Diag=zeros(1,Nx-1); %Up_diag=D.*ones(1,Nx-2); %Low_diag=D.*ones(1,Nx-2); Up_diag=D/sy^2.*ones(1,Nx-2); Low_diag=D/sy^2.*ones(1,Nx-2); for jj=1:Ny-1 Diag(jj)=-2*D/sy^2-y(jj+1)/sx; end %---B---% B=zeros(Nx-1,1); for kk=1:Ny-1 % B(kk,1)=(y(kk+1)*sx-h(ii,kk+1)*sy+Bata*sy^2)*p1(ii,kk+1)+h(ii,kk+1)*sy*p1(ii,kk+2); B(kk,1)=-((y(kk+1)/sx-h(ii,kk+1)/sy+Bata)*p1(ii,kk+1)+h(ii,kk+1)/sy*p1(ii,kk+2)); end %B(1,1)= B(1,1)-D*p1(ii+1,1); % B(end,1)= B(end,1)-D*p1(ii+1,end); B(1,1)= B(1,1)-D/sy^2*p1(ii+1,1); B(end,1)= B(end,1)-D/sy^2*p1(ii+1,end); p1(ii+1,2:Nx)=zhuiganfa(Low_diag,Diag,Up_diag,B); %追赶法解线性方程组 % p1(ii+1, =( p1(ii+1,-min(p1(ii+1,))./(max(p1(ii+1,)-min(p1(ii+1,));%p1(p1<0)=0; end figure(2) mesh(x,y,p1) xlabel('x','FontSize',12) ylabel('y','FontSize',12) zlabel('p(x,y)','FontSize',12)  图1.jpg  图2.jpg |
回复此楼» 猜你喜欢
售SCI文章，我:8O5.5.1.O.54,科目齐全,可+急
已经有7人回复
-
售SCI一区T0P文章，我:8O.55.1.O.54,科目全,可伽急
已经有7人回复
-
售SCI文章，我:8O.5.5.1O.54,科目全,可十急
已经有7人回复
-
售SCI一区T0P文章，我:8.O.55.1.O.54,科目齐全,可+急
已经有6人回复
-
售SCI文章，我:8O5.5.1.O.54,科目齐全,可+急
已经有8人回复
-
售SCI一区T0P文章，我:8.O.55.1.O54,科目全,可伽急
已经有10人回复
-
售SCI一区T0P文章，我:8.O.55.1.O54,科目全,可伽急
已经有7人回复
-
售一区SCI文章T0P，我:8O.551.O54,科目全,可十急
已经有8人回复
-
售SCI一区文章，我:8.O.55.1.O.54,科目齐全,可伽急
已经有11人回复
-
售SCI文章，我:8O5.5.1.O.54,科目齐全,可+急
已经有8人回复
chen17307
新虫 (初入文坛)
- 应助: 0 (幼儿园)
- 金币: 56.2
- 散金: 22
- 帖子: 12
- 在线: 9.8小时
- 虫号: 2932786
- 注册: 2014-01-15
- 性别: GG
- 专业: 计算数学与科学工程计算
2楼2017-08-14 14:34:31
|
本帖内容被屏蔽 |
10楼2018-06-09 15:08:11
 |
3楼2017-08-18 05:43:37
 |
4楼2017-08-18 20:14:44
 |
5楼2017-08-19 04:15:36
 |
6楼2017-08-20 19:36:55
 |
7楼2017-08-21 02:41:21
|
|
8楼2017-08-21 09:54:39
sanshiyayan
新虫 (小有名气)
- 应助: 3 (幼儿园)
- 金币: 1845.5
- 红花: 4
- 帖子: 185
- 在线: 177.4小时
- 虫号: 2058588
- 注册: 2012-10-13
- 性别: GG
- 专业: 计算数学与科学工程计算
9楼2018-03-02 15:42:37
30