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_pdepeclear 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 allfor n = linspace(1,length(t),50)plot(x,sol(n,:),'LineWidth',2)endtitle('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; endfunction u0 = DiffusionICfun(x)u0 = 0;endfunction [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;endend
Where am I going wrong?
Thank you.
Best Answer