How to use convolution to solve Fick's 2nd Law?

9 vues (au cours des 30 derniers jours)
Ryan Holman
Ryan Holman le 4 Juin 2021
∂C/∂t=D∇^2 C
I am tying to solve Fick's 2nd law in one dimension with a convolution kernal and convn for fast simulations. But I am not getting the correct results. The final result should appear something like the attached image, with the oxygen diffusing down concentration gradient with time. I cannot find the bug in my code. I think it is either with my boundary conditions or how I am implementing the convn() function. Can someone find the error in this code?
%This program simulates Fickian diffusion of molecular oxygen from a single
%PFOB droplet though water in radial coordinates
%Addjusting time shift seems to adjust position, not rate of reduction.
clear all; clf; close all;
tic
% % diffusivity or diffusion coefficient [m^2/s]. Oxygen in water.
D = 2.5E-9; %m2/s
xi = 1E-6; xf = 1000E-6; dx = 1E-6; %position, meters
ti = 0; tf = 1000; dt = 1; %time, seconds
cf = (2.828E6)/330; % ambient concentration value [ng/L]
ci = 2.828E6; % initial concentration value [ng/L]
X = xi:dx:xf;
rr=round(xf/dx);
R=1E-6;
dyn=round(tf/dt)
T=ci*ones(rr,dyn,'double');
for t=1:dyn
for r=1:rr
w=1/(4*D*dt*t)*(exp(-((X(r))*(X(r)))./(4*D*dt*t)));
G2(r,t)=w;
end
end
%initial condition
T(1,1)=ci;
T(:,1)=cf;
tall(G2);tall(T);
Tout=convn(G2,T(:,1),'same');
gather(G2);gather(T); gather(Tout);
for t=1:dyn
disp(t)
plot(X,Tout(:,t))
title(num2str(t))
%ylim([cf ci])
pause
end
toc

Réponses (0)

Catégories

En savoir plus sur Fluid Dynamics 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!

Translated by