How to set this time derivative boundary condition using MOL?

4 次查看（过去 30 天）

显示更早的评论

Hashim 2023-7-17

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1997113-how-to-set-this-time-derivative-boundary-condition-using-mol

编辑： Torsten 2023-7-17

在 MATLAB Online 中打开

LHS = Left hand side aka

BC = Boundary condition

Hi,

I have a system of coupled raection diffusion PDEs as such.

For my

species I have a LHS at

such as

I am trying to solve this using the method of lines by using ODE15s solver for which I have written the function as such. But I do not know how to incorporate the

term into the script and how to assign it to be my LHS boundary condition.

k1=0.073317274764052;
k2=1.052512778689219e-04;
tau=0.466500000000000;
tau_2=0.699700000000000; 
dchi=0.100000000000000;
N=11; 
tspan_2=linspace(tau, tau_2, N);
[t, CmBC] = ode15s(@(t, Cm2) BCodefun(t, Cm2, k1, k2, dchi), tspan_2, 0.9, []);
function dCdt2BC = BCodefun(t, CmBC, k1, k2, dchi)
         dCdt2BC  = (k1/k2).*(((CmBC).*(CmBC^2-CmBC))/dchi); 
end

This gives me a 1*11 matrix (code till above should function). You can see below how I am using a multi point approximation to set a no flux BC at my left hand. I am not sure how to assign my ODE output as my left hand side BC in my MOL implementation.

% What I have been using to set a no flux boundary condition on the LHS
Cm2(1)  = (-Cm2(3)+4*Cm2(2))./3;
% WHat I think I should do to set the output from my ODE solution to set my
% LHS BC
Cm2(1)  = CmBC;
Function definitions in a script must appear at the end of the file.
Move all statements after the "BCodefun" function definition to before the first local function definition.

9 个评论
显示 7更早的评论隐藏 7更早的评论

Hashim 2023-7-17

在 MATLAB Online 中打开

There is no

in this boundary condition implementation only

. Also, just to clarify I am using N to discretize both the spatial as well as the temporal domain.

This is what I am trying to do Cm2 being

[t, Cm2(1)] = ode15s(@(t, Cm2) BCodefun(t, Cm2, k1, k2, dchi), tspan_2, 1, []);
function dCdt2BC = BCodefun(t, Cm2, k1, k2, dchi)
            dCdt2BC  = (k1/k2).*((Cm2-Cm2)/dchi)*(Cm2^2-Cm2); 
end

If I go your route I get the following error.

Index exceeds the number of array elements. Index must not exceed 1.
Error in PcWise_NDOCP_Pulse2v0/fpde2/BCodefun (line 133)
            dCdt2BC  = (k1/k2).*((Cm2(2)-Cm2(1))/dchi)*(Cm2(1)^2-Cm2(1));
Error in PcWise_NDOCP_Pulse2v0>@(t,Cm2)BCodefun(t,Cm2,k1,k2,dchi) (line 123)
    [t, Cm2(1)] = ode15s(@(t, Cm2) BCodefun(t, Cm2, k1, k2, dchi), tspan_2, 1, []);
Error in odearguments (line 92)
f0 = ode(t0,y0,args{:});   % ODE15I sets args{1} to yp0.
Error in ode15s (line 153)
    odearguments(odeIsFuncHandle, odeTreatAsMFile, solver_name, ode, tspan, y0, options, varargin);
Error in PcWise_NDOCP_Pulse2v0/fpde2 (line 123)
    [t, Cm2(1)] = ode15s(@(t, Cm2) BCodefun(t, Cm2, k1, k2, dchi), tspan_2, 1, []);
Error in PcWise_NDOCP_Pulse2v0>@(t,C2)fpde2(t,C2,alpha,kappa,eta,gamma,mu,k1,k2,dchi,N,tspan_2) (line 104)
[t, C2]      = ode15s(@(t, C2) fpde2(t, C2, alpha, kappa, eta, gamma, mu, k1, k2, dchi, N, tspan_2), tspan_2, IC2, []);
Error in odearguments (line 92)
f0 = ode(t0,y0,args{:});   % ODE15I sets args{1} to yp0.
Error in ode15s (line 153)
    odearguments(odeIsFuncHandle, odeTreatAsMFile, solver_name, ode, tspan, y0, options, varargin);
Error in PcWise_NDOCP_Pulse2v0 (line 104)
[t, C2]      = ode15s(@(t, C2) fpde2(t, C2, alpha, kappa, eta, gamma, mu, k1, k2, dchi, N, tspan_2), tspan_2, IC2, []);