getting error while solving PDE using ode45
显示 更早的评论
m=50;
delx=1/(m-1);
h_int=0.0005;
final_time=1;
N_samples=(final_time/h_int)+1;
t0=0;
tf=final_time;
p_mat=zeros(N_samples,m);
p_mat(1,:)=0;
for n=2:N_samples-1
[t(n),p_mat(n,:)]=ode45(@(t(n),p_mat(n,:)) GH(t(n),p_mat(n,:)), [t0 tf], p_mat(1,:));
end
plot(t, p_mat);
xlabel('time'),ylabel('pressure');
%%
function F_X= GH(t,p)
m=50;
delx=1/(m-1);
F_X = zeros(1,m);
% for i=1 and m
k=0.98*10^-12;
mu=0.001;
phi=0.4;
K1=k/((1-phi)*mu);
pg=2*10^6;
ph=30*10^6;
p(1)=pg;
p(m)=ph;
for i=2:m-1
F_X(i)=K1*((p(i+1)-2*p(i)+p(i-1))/(delx)^2);% K1 is k/(1-phi)*mu
end
end
2 个评论
Jan
2021-12-23
This is the code. Please add a copy of the complete error message.
MANISH KUMAR GUPTA
2021-12-23
采纳的回答
m=50;
delx=1/(m-1);
h_int=0.0005e7;
final_time=1.0e7;
N_samples=(final_time/h_int)+1;
t0=0;
tf=final_time;
tspan=linspace(t0,tf,N_samples);
k=0.98*10^-12;
mu=0.001;
phi=0.4;
K1=k/((1-phi)*mu);
pg=2*10^6;
ph=30*10^6;
p0=zeros(m,1);
p0(1)=pg;
p0(m)=ph;
[t p_mat]=ode45(@(t,p)GH(t,p,m,delx,K1),tspan,p0);
plot(t,p_mat)
xlabel('time'),ylabel('pressure');
function F_X= GH(t,p,m,delx,K1)
F_X = zeros(m,1);
for i=2:m-1
F_X(i)=K1*(p(i+1)-2*p(i)+p(i-1))/delx^2;% K1 is k/(1-phi)*mu
end
end
更多回答（0 个）
类别
在 帮助中心 和 File Exchange 中查找有关 Ordinary Differential Equations 的更多信息
产品
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
