Error using bvp4c: "Warning: Unable to meet the tolerance without using more than 5000 mesh points. The last mesh of 4921 points and the solution are available in the output argument."
15 vues (au cours des 30 derniers jours)
Afficher commentaires plus anciens
Hello,
I'm using bvp4c to solve for a reaction-diffusion equation. When I set the initial concentration to 100, I get this error saying
"Warning: Unable to meet the tolerance without using more than 5000 mesh points.
The last mesh of 4921 points and the solution are available in the output argument.
The maximum residual is 0.512336, while requested accuracy is 0.001.
> In bvp4c (line 323)
In steady_state_bvp4c_FINAL (line 38)
Warning: Imaginary parts of complex X and/or Y arguments ignored
> In steady_state_bvp4c_FINAL (line 44)
Warning: Imaginary parts of complex X and/or Y arguments ignored
> In steady_state_bvp4c_FINAL (line 50) "
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 steady_state_bvp4c_FINAL
clear all; clc;
%Diffusion coefficient
D_ij = 30; %3.0*10^-7 cm^2/s -> 30 um^2/s
%Initial concentration at x=0
L0 = 100;
%Total length
x_f = 3500;
%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 Receptors
N_T =1.7;
L = L0./2;
solinit = bvpinit(linspace(0,x_f,11),[L 0]);
sol = bvp4c(@(x,y)odefcn(x,y,D_ij,k_1,k_2,k_3,d_1,d_2,d_3,K_1,K_2,K_3,K_4,K_i,N_T),@twobc,solinit);
figure(1)
plot(sol.x,(sol.y(1,:)),'LineWidth',1)
xlabel('Distance (\mum)')
ylabel('Concentration (nM)')
axis([0 x_f 0 L0])
figure(2)
plot(sol.x,(sol.y(1,:))./L0,'LineWidth',1)
xlabel('Distance (\mum)')
ylabel('Normalized Concentration (C_L/C_L_0)')
axis([0 x_f 0 1])
function dy = odefcn(x,y,D_ij,k_1,k_2,k_3,d_1,d_2,d_3,K_1,K_2,K_3,K_4,K_i,N_T)
C_L = y(1);
dC_Ldx = y(2);
N_r = (-1 - (C_L./K_1) + sqrt((-1 - (C_L./K_1)).^2 -4.*(-(2./K_4) - ((2.*C_L)./(K_2.*K_4)) - ((2.*C_L)./(K_2.*K_4.*K_i)) - ((2.*C_L.^2)/(K_2.*K_3.*K_4.*K_i))).*N_T))/(4.*((1./K_4) + (C_L./(K_2.*K_4)) + (C_L./(K_2.*K_4.*K_i)) + (C_L.^2./(K_2.*K_3.*K_4.*K_i))));
N_R = (C_L./K_1).*N_r;
N_r2 = (1./K_4).*(N_r).^2;
N_rR = (C_L./(2.*K_2.*K_4)).*(N_r).^2;
N_pR = (C_L./(2.*K_2.*K_4.*K_i)).*(N_r).^2;
N_R2 = ((C_L.^2)./(K_2.*K_3.*K_4.*K_i)).*(N_r).^2;
r = N_T-(2.*N_R2)-(4.*N_pR)-N_R-(4.*N_rR)-(2.*N_r2);
R = N_T-(2.*N_R2)-(4.*N_pR)-N_r-(4.*N_rR)-(2.*N_r2);
r_2 = (N_T./2)-(N_R./2)-N_R2-(2.*N_pR)-(N_r./2)-(2.*N_rR);
rR = (N_T./4)-(N_R./4)-(N_R2./2)-N_pR-(N_r./4)-(N_r2./2);
pR = (N_T./4)-(N_R./4)-(N_R2./2)-(N_r./4)-N_rR-(N_r2./2);
R_2 = (N_T./2)-(N_R./2)-(2.*N_pR)-(N_r./2)-(2.*N_rR)-N_r2;
R_L = (2.*d_1.*(R))-(2.*k_1.*C_L.*r)+...
((2.*d_2.*(rR)))-((k_2.*C_L.*(r_2)))+...
(d_3.*(R_2))-(2.*k_3.*C_L.*(pR));
dy = zeros(2,1);
dy(1) = dC_Ldx;
dy(2) = (-R_L/D_ij);
end
function res = twobc(ya,yb)
res = [ ya(1)-(L0); yb(2) ];
end
end
0 commentaires
Réponses (1)
Torsten
le 16 Nov 2018
If the code works for values of L0 smaller than 100, try approaching 100 by subsequently calling bvp4c with values L0start, L01, L02,...,L0=100.
Take a look at
https://de.mathworks.com/matlabcentral/answers/428271-how-to-write-matlab-bvp4c-code-for-moving-boundary-surface-problem
to see how this can be implemented.
Best wishes
Torsten.
0 commentaires
Voir également
Catégories
En savoir plus sur Boundary Value Problems dans Help Center et File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!