PDEPE: Unable to meet integration tolerances without reducing the step size below the smallest value allowed

11 views (last 30 days)
I am trying to integrate the following set of PDE
Equation 8 presents a 4th order spacial derivative, so my u vector is u = [C,rho,phi,P] where P is the second derivative of phi with respect to x. P = phi/xx is the fourth equation of the PDE system. Equations 9,10, 11 and 12 are boundary conditions valid for both boundaries, x = +1, x = -1. The initial condition for P is assumed to be P(0,t) = 0.
function [c,f,s] = pdefun_bsk(xmesh,tspan,u,dudx)
eps = 0.01;
nu = 0.5;
d_c = 10;
l_c = eps*sqrt(d_c);
c = [1; 1; 0; 0];
f = [eps/d_c*(dudx(1) + u(2)*dudx(3) + nu*u(1)/(1-nu*u(1))*dudx(1)); eps/d_c*(dudx(2) + u(1)*dudx(3) + nu*u(2)/(1-nu*u(1))*dudx(1)); dudx(3); dudx(4)];
s = [0; 0; -u(4);-u(2)/eps^2-u(4)/l_c^2];
end
function [pl,ql,pr,qr] = bcfun_bsk(xl,ul,xr,ur,tspan)
v = 0.1;
pl = [0; 0; ul(3)+v; 0];
ql = [1; 1; 0; 1];
pr = [0; 0; ur(3)-v; 0];
qr = [1; 1; 0; 1];
end
function u0 = icfun_bsk(xmesh)
v = 0.1;
u0 = [1; 0; v*xmesh; 0];
end
clear
clc
%% Data
m = 0;
xmesh = linspace(-1,1,1000);
tend = 2;
tspan = linspace(0,tend,1000);
tlength = length(tspan);
xlength = length(xmesh);
%% Solution
SOL = pdepe(m,@pdefun_bsk,@icfun_bsk,@bcfun_bsk,xmesh,tspan);
I get the following error
Warning: Failure at t=5.721307e-01. Unable to meet integration tolerances without reducing the step
size below the smallest value allowed (1.776357e-15) at time t.
> In ode15s (line 668)
In pdepe (line 289)
Warning: Time integration has failed. Solution is available at requested time points up to
t=5.705706e-01.
> In pdepe (line 303)
Is the problem due to the fact that I am trying to solve a 4-th order PDE? I tried to increase the tolerance (AbsTol e RelTol) to a few units and I do get a solution without errors, but it seems to be wrong.
  2 Comments
Alessandro Ferreri
Alessandro Ferreri on 14 Mar 2021
Hi Bill,
Thank you very much for the quick reply! You are right, the fourth variable does diverge on the right-hand side. I found this set of equations in a paper and a solution is presented. The set was solved with a different software though. Do you think the problem is numerical, or is there an issue with the equations themselves? I also tried to increase the number of steps in both time and space vectors up to 10000, without any luck.

Sign in to comment.

Answers (1)

Bill Greene
Bill Greene on 14 Mar 2021
I have occasionally seen similar problems with pdepe in the past.
I have written a PDE solver that has input similar to pdepe but
does not exhibit this spurious divergence. If you want to try it,
you can download it here, pde1dM, and simply change your call from pdepe to pde1dM.
As an aside, the number of points in tspan has no effect on the solution
accuracy so I suggest you use far less than the 1000 you have chosen.
The number of points in xmesh does affect the accuracy but 1000 is
a very fine mesh; I suggest trying say, 50, and then increase that to
say, 100 and see how it affects the solution.
  1 Comment
Alessandro Ferreri
Alessandro Ferreri on 15 Mar 2021
Edited: Alessandro Ferreri on 15 Mar 2021
Hi Bill,
thank you for the suggestion about xmesh and tspan vectors.
I tried your pde1dM function and in fact it yields a solution without any errors. It is the same solution that I get if I increase the relative tolerance in pdepe up to 30. However, this behavior differs from what I find in the original paper. I checked the way the equations are written and the parameters many times.
The right solution should be the one presented in the following paper in fig 7. EDIT: The most noticeable difference is that the third variable, phi, changes significantly between t = 0 and t = 2, whereas in my solution it barely changes.
Zhao, H., 2011. Diffuse-charge dynamics of ionic liquids in electrochemical systems. Physical Review E, 84(5), p.051504.

Sign in to comment.

Community Treasure Hunt

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

Start Hunting!