Why do I receive an error in calling "decic"?

1 次查看(过去 30 天)
Hey everyone,
I'am trying to solve a DAE-system with the help of ode15s. I already tried ode15i in another code, but had some problems. So I want to check if ode15s performs better.
When computing the consistent initial conditions for ode15i I receive the following error message:
Not enough input arguments.
Error in DGL_Membrane>@(t,y,yp)M(t,y)*yp-F(t,y) (line 101)
implicitDAE = @(t,y,yp) M(t,y)*yp - F(t,y)
Error in decic (line 66)
res = feval(odefun,t0,y0,yp0,varargin{:});
Error in DGL_Membrane (line 102)
[y0,yp0] = decic(implicitDAE,0,y0est,[],yp0est,[],opt);
My code is the following (following the procedure of https://de.mathworks.com/help/symbolic/solve-daes-using-mass-matrix-solvers.html)
function [ySol,tSol] = DGL_Membrane(a,k,L)
syms asym bsym csym dsym fsym hsym isym ksym jsym x(t) y(t) z(t) xx(t) e(t)
Dx = diff(x);
Dy = diff(y);
Dz = diff(z);
Dxx = diff(xx);
Eqn1 = -bsym * diff(x(t),1) - x*bsym*csym*dsym*diff(y(t),1) == asym;
Eqn2 = -fsym*diff(z(t),1) - z * fsym * hsym * dsym * diff(y(t),1) == e;
Eqn3 = -isym*Dxx(t) - xx * isym*jsym * dsym * diff(y(t),1) == ksym;
Eqn4 = csym * x(t) + hsym* z(t) + jsym*xx(t) ==0;
Eqn5 = csym * asym + hsym* e(t) + jsym*ksym == 0;
DAEs=[Eqn1,Eqn2,Eqn3,Eqn4,Eqn5];
DAEvars = [x(t);y(t);z(t);xx(t);e(t)];
[M,f] = massMatrixForm(DAEs,DAEvars);
pDAEs = symvar(DAEs);
pDAEvars = symvar(DAEvars);
extraParams = setdiff(pDAEs, pDAEvars);
M = odeFunction(M, DAEvars,bsym,csym,dsym,fsym,hsym,isym,jsym);
f = odeFunction(f, DAEvars,asym,csym,hsym,jsym,ksym);
% Declare parameters
asym = a;
bsym=1e-6*L; %P1
csym=-1; %Z1
dsym=9.64853321233100184e4/(8.3145*293.15); %F/RT
fsym=0.01e-6*L; %P2
hsym=-2;%Z2
isym=10e-6*L; %P3
jsym=1; %Z3
ksym=k;
c1f=0.2e-3; %feed concentrations mol/m^3
c2f=0.2e-3;
c3f=0.6e-3;
F = @(t,Y) f(t,Y,asym,csym,hsym,jsym,ksym);
y0est = [c1f; 0; c2f;c3f;(-csym * asym - jsym*ksym)/hsym];
yp0est = zeros(5,1);
opt = odeset('Mass', M, 'InitialSlope', yp0est,...
'RelTol', 10.0^(-6), 'AbsTol' , 10.0^(-6));
implicitDAE = @(t,y,yp) M(t,y)*yp - F(t,y)
[y0,yp0] = decic(implicitDAE,0,y0est,[],yp0est,[],opt);
opt = odeset(opt, 'InitialSlope', yp0);
time=linspace(0,L,500);
[tSol,ySol] = ode15s(F,time,y0est,opt);
plot(tSol,ySol(:,[1,3,4,5]),'LineWidth',2)
end
Does anyone have an idea what's wrong with my decic-input?
Thanks for help and best regards !
  3 个评论
Tobias Hubach
Tobias Hubach 2022-1-26
Hello Torsten,
When I refer to the decic documentation it says 7 inputs (or 6 without options)
Torsten
Torsten 2022-1-26
编辑:Torsten 2022-1-26
You are right, it was the symbolic decic function.
I think you will have to define MM the same as F
MM = @(t,y) M(t,y,asym,csym,hsym,jsym,ksym)
and use MM in the definition of implicitDAE.

请先登录,再进行评论。

回答(1 个)

Ishu
Ishu 2023-11-28
Hi Tobias Hubach,
I understand that you are trying to use the "decic" function of MATLAB to solve the differential equations.
I tried to reproduce the error at my end and it also shows me the same error. After that I checked the function handles "M" and "F" and the input arguments required by the handle "M" are a lot not just the "(t,y)". Because of this reason it gives the error of "Not enough input arguments".
Further I looked into the code and found that you are changing the defination of "M" after its initialization. Explained in the given code snippet.
[M,f] = massMatrixForm(DAEs,DAEvars); % you have initialised here
M
M = 
pDAEs = symvar(DAEs);
pDAEvars = symvar(DAEvars);
extraParams = setdiff(pDAEs, pDAEvars);
M = odeFunction(M, DAEvars,bsym,csym,dsym,fsym,hsym,isym,jsym) % here you have changed "M"
M = function_handle with value:
@(t,in2,param1,param2,param3,param4,param5,param6,param7)reshape([-param1,0.0,0.0,0.0,0.0,-param1.*param2.*param3.*in2(1,:),-param3.*param4.*param5.*in2(3,:),-param3.*param6.*param7.*in2(4,:),0.0,0.0,0.0,-param4,0.0,0.0,0.0,0.0,0.0,-param6,0.0,0.0,0.0,0.0,0.0,0.0,0.0],[5,5])
Hope it helps.

类别

Help CenterFile Exchange 中查找有关 Numerical Integration and Differential Equations 的更多信息

产品


版本

R2021a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by