How to pass solution of pdefun1 as IC-2 for pdefun2
조회 수: 2(최근 30일)
표시 이전 댓글
my pde is

Given below is my effort so far...
function pdepe_philip_v3b
% 1-b Extend the above calculation to model a potential step to any chosen
% potential assuming that the ratio of the surface concentrations of Os2+
% and Os3+ obeys the Nernst equation both before and after the potential
% step. Now we want to further add a potential step which will make it a
% double potential step.
global D n F R
n = 1; % No. of Electrons Involved
F = 96485.3329; % sA/mol
R = 8.3145; % kgcm^2/s^2.mol.K
D = 5e-06; % cm^2/s
Area = 1; % cm^2
Os_bulk = 1e-05; % mol/cm^3
N = 10;
m = 0; % Cartesian Co-ordinates
% logarithmized xmesh for semi infinite boundary condition
x = logspace(log10(0.000006), log10(0.06), N); x(1) = 0;
% Time Span-1
t1 = linspace(0, 10, N); % s
E = 0.8;
sol1 = pdepe(m, @pdev3bpde, @(x)pdev3bic(x,Os_bulk), @(xl, ul, xr, ur, t)...
pdev3bbc(xl, ul, xr, ur, t, E, Os_bulk), x, t1);
c1 = sol1(:, :, 1);
I_num = ( n*F*Area) * c1(1,:) .* sqrt(D./(t1*pi));
% I Flux For Time Span-1
figure(1);
plot(t1, I_num, 'r--', 'LineWidth',1);
xlabel('t (s)'); ylabel('I (A/cm^2s)');
hold on
% Time Span-2
t2 = linspace(10, 20, N); % s
c_new = (c1(N,:))';
E = -1.6;
sol2 = pdepe(m, @pdev3bpdeNEW, @(x)pdev3bicNEW(x,c_new),@(xl, ul, xr, ur, t)...
pdev3bbcNEW(xl, ul, xr, ur, t, E, c_new), x, t2);
end
%% pdepe Function
%%
function [a, f, s] = pdev3bpde(x, t, u, DuDx)
global D
a = 1;
f = D*DuDx;
s = 0;
end
%% pdepe Function
%%
function [a, f, s] = pdev3bpdeNEW(x, t2, u, DuDx)
global D
a = 1;
f = D*DuDx;
s = 0;
end
%% Initial Condition-1
%%
function u0 = pdev3bic(x, Os_bulk)
u0 = Os_bulk;
end
%% Initial Condition-2
%%
function u0 = pdev3bicNEW(x, c_new)
u0 = c_new.*size(x);
end
%% Boundry Condition-1
%%
function [pl, ql, pr, qr] = pdev3bbc(xl, ul, xr, ur, t, E, Os_bulk)
global n F R
E0 = 0.2;
alpha = exp((E-E0).*((n*F)/(R*298.15)));
pl = ul - (Os_bulk./(1+alpha));
ql = 0;
pr = ur - Os_bulk;
qr = 0;
end
%% Boundry Condition-2
%%
function [pl, ql, pr, qr] = pdev3bbcNEW(xl, ul, xr, ur, t, E, Os_bulk)
global n F R
E0 = 0.2;
alpha = exp((E-E0).*((n*F)/(R*298.15)));
pl = ul - (Os_bulk./(1+alpha));
ql = 0;
pr = ur - Os_bulk;
qr = 0;
end
%%
%%
채택된 답변
추가 답변(0개)
참고 항목
범주
Find more on Test Model Components in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!