Discover MakerZone

MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi

Learn more

Discover what MATLAB® can do for your career.

Opportunities for recent engineering grads.

Apply Today

solving coupled PDEs using method of lines

Asked by tejas grewal on 18 Jul 2013

I have 6 coupled pde (mass balance eqn for 6 component coupled through rate of reaction Ri)of type:

dFi/dt = dFi/dz + Ri with initial condition: Fi(0,z) = 0

     boundary cond. Fi(t,0) = Fi0

i tried using method of lines

initial cond.

F0 = [FinE repmat(0,Z-1); % FinE inlet flow rate of ethanol

      FinW repmat(0,Z-1);  % Finw inlet flow rate of water
      zeros(1,Z)                    % products inlet flow rate
      zeros(1,Z)
      zeros(1,Z) 
      zeros(1,Z)] 
tspan = linspace(0,0.2,N)
[t,F] = ode45('rates',tspan,F0)
******************************************
function dFdt = rates(t,F)
i = 1
k1 = k01*exp(-Ea1/RT)
k1'= ...

dFdt(i,:) = (-((F(i,:)-FinE)./(dx(1,1)*Area))+eta*(-k1*F(i+1,:)*F(i,:)+k3*F(i+3,:))

dFdt(i+1,:) = (-((F(i+1,:)-Finw)./(dx(1,1)*Area))+eta*(-k1*F(i+1,:)*F(i,:))

dFdt(i+2,:) = (-((F(i+2,:)-0)./(dx(1,1)*Area))+eta*(k1*F(i+1,:)*F(i,:)+k2*F(i+4,:))

dFdt(i+3,:) = ...

dFdt(i+4,:) = ....

dFdt(i+5,:) = ....

i = 7:6:6*Z-5

dFdt(i,:) = (-((F(i,:)-F(i-6,:))./(dx(1,1)*Area))+eta*(-k1*F(i+1,:)*F(i,:)+k3*F(i+3,:))

dFdt(i+1,:) = (-((F(i+1,:)-F(i-5,:))./(dx(1,1)*Area))+eta*(-k1*F(i+1,:)*F(i,:))

dFdt(i+2,:) = (-((F(i+2,:)-F(i-4,:))./(dx(1,1)*Area))+eta*(k1*F(i+1,:)*F(i,:)+k2*F(i+4,:))

dFdt(i+3,:) = ...

dFdt(i+4,:) = ....

dFdt(i+5,:) = ....

this code is giving NAN values for the above initial conditions but running with initial condition Fi(0,z) = Fi0 why is this so.. also if any one can suggest the alternative method to solve for Fi vs t as well as vs z

0 Comments

tejas  grewal

Products

No products are associated with this question.

1 Answer

Answer by Deepak Ramaswamy on 19 Jul 2013
Edited by Deepak Ramaswamy on 19 Jul 2013

This is with regards to the second part of your question: This problem should be solvable via pdepe. I could not infer the form of Ri and the right BC is not specified. But the form for pdepe (check out pdepe doc for what below terms such as m, c, f, s etc. mean) would look something like:

m = 0

c

+-                  -+ 
|  1, 0, 0, 0, 0, 0  | 
|                    | 
|  0, 1, 0, 0, 0, 0  | 
|                    | 
|  0, 0, 1, 0, 0, 0  | 
|                    | 
|  0, 0, 0, 1, 0, 0  | 
|                    | 
|  0, 0, 0, 0, 1, 0  | 
|                    | 
|  0, 0, 0, 0, 0, 1  | 
+-                  -+

f

+-          -+ 
|  F1(z, t)  | 
|            | 
|  F2(z, t)  | 
|            | 
|  F3(z, t)  | 
|            | 
|  F4(z, t)  | 
|            | 
|  F5(z, t)  | 
|            | 
|  F6(z, t)  | 
+-          -+
*s*
+-                                                                -+ 
|  R1(F1(z, t), F2(z, t), F3(z, t), F4(z, t), F5(z, t), F6(z, t))  | 
|                                                                  | 
|  R2(F1(z, t), F2(z, t), F3(z, t), F4(z, t), F5(z, t), F6(z, t))  | 
|                                                                  | 
|  R3(F1(z, t), F2(z, t), F3(z, t), F4(z, t), F5(z, t), F6(z, t))  | 
|                                                                  | 
|  R4(F1(z, t), F2(z, t), F3(z, t), F4(z, t), F5(z, t), F6(z, t))  | 
|                                                                  | 
|  R5(F1(z, t), F2(z, t), F3(z, t), F4(z, t), F5(z, t), F6(z, t))  | 
|                                                                  | 
|  R6(F1(z, t), F2(z, t), F3(z, t), F4(z, t), F5(z, t), F6(z, t))  | 
+-                                                                -+

pLeft (p for left BC)

+-               -+ 
|  F1(z, t) - F0  | 
|                 | 
|  F2(z, t) - F0  | 
|                 | 
|  F3(z, t) - F0  | 
|                 | 
|  F4(z, t) - F0  | 
|                 | 
|  F5(z, t) - F0  | 
|                 | 
|  F6(z, t) - F0  | 
+-               -+
qLeft (q for left BC)
+-   -+ 
|  0  | 
|     | 
|  0  | 
|     | 
|  0  | 
|     | 
|  0  | 
|     | 
|  0  | 
|     | 
|  0  | 
+-   -+

0 Comments

Deepak Ramaswamy

Contact us