Can we solve 3 pde simultaneously using pdepe with respect to space and time??
1 次查看(过去 30 天)
显示 更早的评论
I have a system of 3 equations with 3 variables and corresponding boundary conditions. I am using pdepe to solve these. Is it possible? Can anyone help me correcting the code I have written? Sharing the equation in an image attached. 3rd equation is : dQp/dt=Keq*Cp*(lambda-(sigma+charge)*Qp)^charge-Qp*Csin^charge
%% Main function function sma_porediffusion %% globalising variables and constants global v x F Ac Dcoe ap Cin charge lambda epsilon k_eq sigma t m L dc N Run_time CV Csin; m=0;
%% specifying system parameters ap=.0023; %particle size epsilon_e= 0.7; %porosity epsilon_p=0.8; % L =5.2; %cm dc= 0.7; %cm F=0.5;%ml/min Ac=3.142*(dc/2)^2; %cm2 v=F/Ac;%cm/min Dcoe=4.43*10^-3;%cm^2/min charge=2.14; lambda= 1833; % capacity factor in mM Cin=0.0085; %mg/ml k_eq= 0.12; % equilibrium coefficient = kads/kdes sigma= 49; % steric hinderence factor N=50; % meshing intervals Run_time =1; % min Csin= 1; %mM mobile phase salt concentration % Csout= 400;%mM end salt concentration % dq/dt=keq*C*(lambda-(sigma+charge)*q)^charge-q*(Csalt^charge) % Area= pi*(dc^2)/4; % area of cross section cm^2
%% x=linspace(0,L,N); t=linspace(0,Run_time,N); sol=pdepe(m,@sma_pdecolumnpe,@sma_pdecolumnic,@sma_pdecolumnbc,x,t);
%% Extract the first solution component as conc. conc1 = sol(:,:,1); conc2=conc1./Cin; %surf(x,t,conc2)
%% A solution profile for conc. vs length. % figure, plot(x,conc1(end,:)) % title(strcat('Solution at t = ', num2str(t))) % xlabel('Distance x') % ylabel('Concentration (C)') % % %Plot conc. vs. time % figure, plot(t,conc1) % title('conc vs time: breakthrough conc1') % xlabel('Time (min)') % ylabel('Concentration (mg/ml)')
figure, plot(t,conc1(:,end)) title('conc vs time: breakthrough') xlabel('Time (min)') ylabel('Concentration (mg/ml)')
%% Inlet conditions function u0 = sma_pdecolumnic(x) if x==0 u0=[Cin;0;0]; else u0=[0;0;0]; end end %% describing equations function [c,f,s] = sma_pdecolumnpe(x,t,u,DuDx) %Defining variables %global v Cin charge lambda epsilon k_eq sigma k_a k_d Csin; c=[1;1;1]; f=[Dcoe;0;0].*DuDx; %dq/dt=adsorption coefficient*c(ionic capacity-(steric %factor+characteristic charge)*q)^ char charge - desorption coeff*q*Csaltin^ %char charge
sq=k_eq.*u(1).*(lambda-(sigma+charge).*u(2))^charge- u(2)*Csin^charge; %disp(sq); s=[-v*DuDx(1)-(1-epsilon_e)*k_eq*ap*(u(1)-u(2));(1-epsilon_p)*sq+k_eq*ap*(u(1)-u(2));sq]; end %% boundary conditions function [pl,ql,pr,qr] = sma_pdecolumnbc(xl,ul,xr,ur,t) %global Cin; pl = [ul(1)-Cin;0;0]; ql = [0;0;0]; pr = [0;0;0]; qr = [1;1;1]; end end
Thanks for your inputs
0 个评论
回答(0 个)
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Partial Differential Equation Toolbox 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!