Risultati per

Hello: your help is greatly appreciated, I am trying to solve the following pde numerically with ode45 (if possible). There are three equations in total. I keep getting an error. know parameters: v, D(j), KL(j), cs(j) The IC and BC are provided in the m.file. pde:

    % code
% Initial, final values of independent variable
tspan = [0 7];
P=1, F2=2, F3=3;
c_initial_1= zeros(N+1,1);
c_initial_2= zeros(N+1,1);
c_initial_3= zeros(N+1,1)
if t == 0
c_initial_2(1,2) = cs(2)
c_initial_3(1,3) = cs(3)
c_initial_1(1,1) = cinj
end
  [t, c, kF2, kF3, kP, n2, m2, n3, m3, q, r] = ode45(@ode, tspan, c_initial_1,c_initial_2, c_initial_3);
function [dcdt] = ode(t, c, k2, k3, kP, n2, m2, n3, m3, q, r)
global N dx dxs
dcdt = zeros(N,N)
for i = 1:N
    for j = P:F3 % j== 1:3 1=persulfate, 2 Fraction 2, and 3 Fraction 3
        if j == 1
            epsilon = 1
        else 
            epsilon = 0
        end 
        dcdt(i, j) = (-v/(2*dx))*(c(i+1,j)-c(i-1,j))+(D(j)/dxs)*(c(i+1,j)...
            -2*c(i,j)-c(i-1,j))+(1-epsilon)*kL(j)*(cs(j)-c(i,j))- ...
            (1-epsilon)*k(j)*(c(i,j)^n(j)*(c(i,1)^m(j))- ...
            epsilon*kP*(c(i,j+epsilon)+c(i,j+2*epsilon))^q *c(i,1)^r
        dcdt(N+2,j) = dcdt(N-1,j)
        dcdt(0,j) = dcdt(2,j)
        if t == 0
            dcdt(1,P) = cinj 
        end
    end
end
end
end

dc(i,j)/dt=-v*(dc(i,j)/dx)+D(j)*d/dx(dc(i,j)/dx)+(1-epsilon)*kL(j)*(cs(j)-c(i,j))+(1-epsilon)*k(j)*c(i,j)^n(j) * c(i,1)^m(j) + epsilon*kp*[c(i,j+epsilon)+c(i,j+2*epsilon)]^r * c(i,1)^q