function my_pde
 
% --

m = 1;
R=2.5;% promien szalki
a = 0.318;% cm promien bibulki
r = [linspace(0,a,10) linspace(1.01*a,2.17,50) linspace(2.18,R,20)];
t = linspace(0,40000,1000);
s=0.356;
close all
sol = pdepe(m,@my_pde1pde,@my_pde1ic,@my_pde1bc,r,t);
c = sol(:,:);


figure(1);
surfl(r,t,c);
shading interp
colormap('copper')
title('c(r,t)');
xlabel('dystans r');
ylabel('czas t');
figure(2)
C=sum(c,1)/s/(t(end)-t(1));
semilogy(r,C);%/221.04
title('Srednie stezenie')
figure(3)
for i=1:length(t)
plot(r,c(i,:))
ylim([0 80])
%pause(0.1)
title(sprintf('Stezenie w czasie t= %.2f min',t(i)/60))
drawnow
%figure(4)

end
figure(4)
plot(t, sol(:,find(r==2.17)));

% --------------------------------------------------------------------------

function [c,f,s] = my_pde1pde(x,t,u,DuDx)
% tu sa moje parametry
D=7.2e-6;

tau= 2.25*3600; %s

%a tu zwracam postac funkcji c f i s
c = 1;                                 
f = DuDx*D;                   
s = 0;%-u/tau;                               

% --------------------------------------------------------------------------

function u0 = my_pde1ic(x)
a = 0.318;% cm srednica bibulki
s=0.356;
c0=221.04;%*1;  %mu g/cm3*cm
u0(find(x<a))= c0;%exp(-((x-0.5).^2)/0.01);                                 
u0(find(x>=a))=0;

% --------------------------------------------------------------------------

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

pl = 0;                                 
ql = 1;                             
pr = 0;                            
qr = 1;                                   

