29 views (last 30 days)

Hello,

I am working on solving a system of 3 PDEs using the pdepe solver. I am receiving an error message that it is "unable to meet integration tolerances" which from what I have read, could mean a number of things may be wrong with my setup for using pdepe. I have tried altering a few different things but with no success.

The .pdf above details the system of PDEs that I wish to solve along with the initial and boundary conditions.

One place that I believe may be causing the error is in the boundary conditions in my j term. The term calls for the n variable evaluate at x = 0 but I was unsure if using ul(1) was the correct way to implement that. When I removed that term altogether, I still received a warning saying "unable to meet integration tolerances" but that it failed at a value of t that was slightly greater than the value of t it failed at whilst the ul(1) term was still in place.

Any help or insight is greatly appreciated!

MATLAB code:

global d alpha phi D_0 D_I D_I3 k_edr k_eer k_rg k_ext I2_0 I_0 Keq b gamma theta beta m1 n_0 I3_0 C1 C2 C3 e0

e0 = 1.602e-19; % elementary charge (Coulombs)

d = 12.5e-4; % cm

alpha = 0.2e-4; % cm^-1

phi = 6.1e16; % cm^-2 s^-1

D_0 = 0.4; % cm^2 s^-1

D_I = 4.4e-6; % cm^2 s^-1

D_I3 = 3.6e-6; % cm^2 s^-1

k_edr = 3e-11; % cm^3 s^-1

k_eer = 2e-16; % cm^3 s^-1

k_rg = 1.3e-15; % cm^3 s^-1

k_ext = 10^12; % cm s^-1

n_0 = 1e14; % cm^-3

I2_0 = 0.05; % M

I_0 = 0.6; % M

Keq = 10^8.3;

b = 1;

gamma = 1.4; % gamma/omega^2

theta = 0.45;

beta = 0.35;

m1 = 1.55;

% pde coefficients for ease of use

C1 = gamma*(1-theta)*D_0;

C2 = gamma*(theta/1.5)*D_I;

C3 = gamma*(theta/1.5)*D_I3;

% solving for the equilibrium concentrations of I3_0 and I_0

syms I3

eqn = solve(I3/((I_0 - I3)*(I2_0 - I3)) == Keq, I3);

I3_0 = double(eqn(1));

I_0 = I_0 - I3_0;

x = linspace(0,d,101);

t = linspace(0,10,10);

m = 0;

sol = pdepe(m,@pde,@pde_ic,@pde_bc,x,t);

u1 = sol(:,:,1);

u2 = sol(:,:,2);

u3 = sol(:,:,3);

%---Functions----------------------------------------

% u(1) = n, u(2) = [I-], u(3) = [I3-]

function [c,f,s] = pde(x,t,u,dudx)

global phi alpha k_eer k_edr k_rg C1 C2 C3

S = phi.*alpha.*exp(-alpha.*x)./(k_rg.*u(2)+k_edr.*u(1));

c = [0; 0; 0];

f = [C1; C2; C3] .* dudx;

s = [phi.*alpha.*exp(-alpha.*x) - k_eer.*u(3).*u(1) - k_edr.*S.*u(1)

1.5.*k_eer.*u(3).*u(1) - 1.5.*k_rg.*S.*u(2)

-0.5.*k_eer.*u(3).*u(1) - 0.5.*k_rg.*S.*u(2)];

end

function u0 = pde_ic(x)

global n_0 I_0 I3_0

u0 = [n_0; I_0; I3_0];

end

function [pl,ql,pr,qr] = pde_bc(xl,ul,xr,ur,t)

global C1 C2 C3 e0 theta D_0 k_ext gamma

j = e0.*k_ext.*ul(1);

pl = [-j./(e0.*(1-theta).*D_0.*gamma); 0; 0];

ql = [1./C1; 1; 1];

pr = [0; 1.5.*j./(e0.*(1-theta).*D_0.*gamma); -0.5.*j./(e0.*(1-theta).*D_0.*gamma)];

qr = [1; 1./C2; 1./C3];

end

Warning messages:

When ul(1) is used:

Warning: Failure at t=2.357029e-02. Unable to meet integration tolerances without reducing the step size below the smallest value

allowed (5.551115e-17) at time t.

> In ode15s (line 653)

In pdepe (line 289)

In Fig3 (line 43)

Warning: Time integration has failed. Solution is available at requested time points up to t=0.000000e+00.

> In pdepe (line 303)

In Fig3 (line 43)

When ul(1) is removed:

Warning: Failure at t=3.909675e-02. Unable to meet integration tolerances without reducing the step size below the smallest

value allowed (1.110223e-16) at time t.

> In ode15s (line 653)

In pdepe (line 289)

In Fig3 (line 43)

Warning: Time integration has failed. Solution is available at requested time points up to t=0.000000e+00.

> In pdepe (line 303)

In Fig3 (line 43)

Bill Greene
on 25 Jun 2020

You are trying to solve a system of ODE which is something pdepe is not designed for. Specifically, the pdepe documentation says, "At least one equation must be parabolic.". This means that at least one of your c-coefficients must be non-zero. I'm surprised pdepe does not specifically test for this and issue an error message.

I suggest you look at one of the boundary value solver functions, instead (e.g. bvp4c).

Bill Greene
on 25 Jun 2020

The documentation is not very clear on this point and the error message doesn't help at all.

Opportunities for recent engineering grads.

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

Start Hunting!
## 0 Comments

Sign in to comment.