Error using pdepe - Spatial discretization has failed. Discretization supports only parabolic and elliptic equations, with flux term involving spatial derivative.

3 次查看（过去 30 天）
Jennifer Yang
Jennifer Yang 2018-7-25
回答： Bill Greene 2018-7-28
Hello,
I'm using pdepe to solve for a reaction-diffusion equation. When I set the initial concentration to anything between 0-250, it runs fine, but the moments, I get this error saying
"Error using pdepe Spatial discretization has failed. Discretization supports only parabolic and elliptic equations, with flux term involving spatial derivative."
Am I setting my boundary conditions incorrectly? I want the left hand side (at x = 0, C_L = C_L0) and on the right hand side I want the flux to equal to 0 (at x = x_f, dC_L/dx = 0).
function time_dependent_pdepe
clear all; close all; clc;
D_ij = 30;
L0 = 1000; %c0 [nM]
x_f =40; %Length of domain [um]
maxt = 2; %Max simulation time [s]
m = 0;
x = linspace(0,x_f,1000); %xmesh
t = linspace(0,maxt,50) %tspan
sol = pdepe(m,@DiffusionPDEfun,@DiffusionICfun,@DiffusionBCfun,x,t,[]);
u = sol;
% Plotting
hold all
for n = linspace(1,length(t),50)
plot(x,sol(n,:),'LineWidth',2)
end
title('Time Dependent (Diffusion Only)')
xlabel('Distance (\mum)')
ylabel('Concentration (nM)')
axis([0 x_f 0 L0])
function [c,f,s] = DiffusionPDEfun(x,t,u,dudx)
D = D_ij;
%Rate constants
k_1 = 0.00193;
k_2 = 0.00255;
k_3 = 4.09;
d_1 = 0.007;
d_2 = 3.95*10^-5;
d_3 = 2.26;
%Equilibrium constants
K_1 = 3.63;
K_2 = 0.0155;
K_3 = 0.553;
K_4 = 9.01;
K_i = 0.139;
%Total Recptors
N_T = 1.7;
N_r = -(K_3.*K_4.*K_i.*(K_2 + u - ((8.*K_2.*u.^2.*N_T + 8.*K_2.*K_3.*u.*N_T + K_2.^2.*K_3.*K_4.*K_i + K_3.*K_4.*K_i.*u.^2 + 8.*K_2.^2.*K_3.*K_i.*N_T + 2.*K_2.*K_3.*K_4.*K_i.*u + 8.*K_2.*K_3.*K_i.*u.*N_T)./(K_3.*K_4.*K_i)).^(1./2)))./(4.*u.^2 + 4.*K_3.*u + 4.*K_2.*K_3.*K_i + 4.*K_3.*K_i.*u);
N_R = (u./K_1).*N_r;
N_r2 = (1./K_4).*(N_r).^2;
N_rR = (u./(2.*K_2.*K_4)).*(N_r).^2;
N_pR = (u./(2.*K_2.*K_4.*K_i)).*(N_r).^2;
N_R2 = ((u.^2)./(K_2.*K_3.*K_4.*K_i)).*(N_r).^2;
R_L = ((2.*d_1.*(N_T-N_r-N_r2-N_rR-N_pR-N_R2))-(2.*k_1.*u.*N_r))+...
((2.*d_2.*(N_T-N_r-N_R-N_r2-N_pR-N_R2)))-((k_2.*u.*(N_T-N_r-N_R-N_rR-N_pR-N_R2)))+...
(d_3.*(N_T-N_r-N_R-N_r2-N_rR-N_pR))-(2.*k_3.*u.*(N_T-N_r-N_R-N_r2-N_rR-N_R2));
% PDE
c = 1;
f = D_ij.*dudx;
s = R_L;
end
function u0 = DiffusionICfun(x)
u0 = 0;
end
function [pl,ql,pr,qr] = DiffusionBCfun(xl,ul,xr,ur,t)
c0 = L0;
% BCs: No flux boundary at the right boundary and constant concentration on
% the left boundary
% At x = 0, c0 = L0
pl = ul-c0;
ql = 0;
% At x = L, dc/dx = 0
pr = 0;
qr = 1;
end
end
Where am I going wrong?
Thank you.

请先登录，再回答此问题。

采纳的回答

Bill Greene
Bill Greene 2018-7-28
The fundamental cause of your problem is that your initial conditions are inconsistent with your boundary condition at x=0. Mathematically this is referred to as an ill-posed problem and the PDE has no solution.
Sometimes, if the inconsistency is not too large, the numerical algorithm in pdepe can adjust the initial conditions to make them consistent with the boundary conditions. But, in your case, for large values of L0, it is not able to do this.
One simple way to fix this problem is to modify your initial condition function slightly:
function u0 = DiffusionICfun(x)
if(x==0)
u0 = L0;
else
u0 = 0;
end
end

更多回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Eigenvalue Problems 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by