How to write a functional coefficient f in Parabolic function
5 views (last 30 days)
Show older comments
I'm solvig a system of 3 nonlinear parabolic PDEs (u1, u2, and u3) on rectangular domain using the parabolic function in matlab. I have some difficulties writing the nonlinear interaction terms. How can I write a function of this form u3/(u2+u1) to be my f1? Any hep is appreciated
xmin=0;xmax=20;ymin=0;ymax=1;
gd = [3;4;xmin;xmax;xmax;xmin;ymax;ymax;ymin;ymin];
sf = 'R1';
ns = [82; 49];
[dl, bt] = decsg(gd, sf, ns);
fig1=figure(1);clf
axis([xmin xmax ymin ymax]);
pdegplot(dl, 'edgeLabels', 'on');
axis equal
title 'Domain With Edge Labels Displayed'
pg = pdeGeometryFromEdges(dl);
es = pg.Edges([1,2,3,4]);
bc1 = pdeBoundaryConditions(es(1),'q',[0;1;1],'g',[1;20;155]);
bc2 = pdeBoundaryConditions(es(2),'q',[0;0;1],'g',[1;1;48]);
bc3 = pdeBoundaryConditions(es(3),'q',[0;0;1],'g',[1;1;24]);
bc4 = pdeBoundaryConditions(es(4),'q',[0;0;0],'g',[1;1;1]);
[p,e,t] = initmesh(dl,'hmax',0.2); % hmax gives max length of an edge
[p,e,t] = refinemesh(dl,p,e,t);
fig2=figure(2);clf
pdemesh(p,e,t);
axis equal
title 'Domain With Mesh Displayed'
np=size(p,2);
u0s=zeros(np,numberOfPDE);
u0s(:,1) = ones(np,1);
ix = find(and(and(and(p(1,:)>5,p(1,:)<8),p(2,:)>0.2),p(2,:)<=1))
u0s(ix,1) = zeros(size(ix));
u0s(:,2) = 1*ones(size(np));
u0s(:,3) = 1*ones(size(np));
u0=zeros(np*numberOfPDE,1);
u0(1:np)=u0s(:,1);
u0(np+1:np+np)=u0s(:,1);
u0(np+1+np:np+np+np)=u0s(:,1);
tmin=0;tmax=0.1;nframes = 20;
tlist = linspace(tmin,tmax,nframes);
numberOfPDE =3;
pb = pde(numberOfPDE);
c1 = 1.6;
c2 = 10.01;
c3 = 0.1;
%--------------
a1 = 0;
a2 = 0;
a3 = 0;
%--------------
d1 = 0.6;
d2 = 0.6;
d3 = 0.1;
%--------------
f1 = ??;
f2 = 0.6;
f3 = 0.1;
u1 = parabolic(u0,tlist,pb,p,e,t,c1,a1,f1,d1);
u2 = parabolic(u0,tlist,pb,p,e,t,c2,a2,f2,d2);
u3 = parabolic(u0,tlist,pb,p,e,t,c3,a3,f3,d3);
fig3=figure(3);clf
colormap(cool);
for j = 1:nframes,
pdesurf(p,t,u1(1:np,j));
axis([xmin xmax ymin ymax 0 1.1]);
shading interp;
Mv(j) = getframe;
end
0 Comments
Answers (0)
See Also
Categories
Find more on Eigenvalue Problems in Help Center and File Exchange
Products
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!