BVP4C:error Unable to solve the collocation equation -- a singular jacobian encountered

3 次查看（过去 30 天）
Minoli Fernando
Minoli Fernando 2023-4-22
编辑： Torsten 2023-4-23
Hello, I am facing this error below is my code
function dYdx = bvp_func_pw1(x, Y, GrT)
%defining the constants
lambda=0.4; %second grade fluid parameter
M=5; %magnetic parameter
P=5; %porous parameter
%GrT=[1:5]; %thermal grashof numbers
GrC=1; %solatal grashof numbers
Pr=1; %prandtil number
Ec=0.1; %Eckert number
Sc=1; %Schimdit number
Kr=-0.3; %Chemical reaction paramter
alpha=2; %Heat generator or absorption paramter
%defining derivatives
%Y(1)=f0, Y(2)=f1, Y(3)=f2,Y(4)=theta, Y(5)=theta1,Y(6)=phi,
%Y(7)=phi1
dYdx(1)=Y(2);
dYdx(2)=Y(3);
dYdx(3)=(1/(2*x*Y(3)))*(Y(2)^2-Y(1)*Y(3)+M*Y(2)+P*Y(2)-GrT*Y(4)-GrC*Y(6));
dYdx(4)=Y(5);
dYdx(5)=-Pr*Y(1)*Y(5)+Pr*(2*Y(2)-alpha)*Y(4)-Pr*Ec*(Y(3)^2+lambda*Y(3)*(Y(2)*Y(3)-Y(1)*dYdx(3))+M*Y(2)^2);
dYdx(6)=Y(7);
dYdx(7)=-Sc*(Y(1)*Y(7)-Kr*Y(6));
dYdx=dYdx(:);
end
function bc = boundary_conditions_pw(Ya,Yb,GrT)
%boundary conditions at x=0 and x=infinity
%boundary conditions at x=0
bc1=[Ya(1); Ya(2)-1; Ya(4)-1; Ya(6)-1];
%boundary conditions at x=infinity
bc2=[Yb(2); Yb(4); Yb(6)];
bc=[bc1;bc2];
end
  2 个评论
Torsten
Torsten 2023-4-22
Please include the part of the code in which you call bvp4c.
If your x-interval contains 0, the error can be explained by the term (1/(2*x*Y(3))) where you would get a division by zero.
Minoli Fernando
Minoli Fernando 2023-4-22
Thankyou for your respond. This is the code that i call bvp4c
%define space from 0 to 2
x=linspace(0,2,11);
GrT=0;
%initial guess for the answer
guess= bvpinit(x,[1;0;1;0;0;1;1]);
%setting options
options = bvpset('Stats','on','RelTol',10e-4,'AbsTol',10e-7);
hold on
for i = 1:5
GrT=GrT+1;
%calling bvp4c to solve the boundary value problem
sol=bvp4c(@bvp_func_pw1,@boundary_conditions_pw,guess,options,GrT);
%export solution
f = sol.y;
if GrT==1
f_1=f(1,:);
f_11=f(2,:);
X_1=sol.x;
else if GrT==2
f_2=f(1,:);
f_22=f(2,:);
X_2=sol.x;
else if GrT==3
f_3=f(1,:);
f_33=f(2,:);
X_3=sol.x;
else if GrT==4
f_4=f(1,:);
f_44=f(2,:);
X_4=sol.x;
else if GrT==5
f_5=f(1,:);
f_55=f(2,:);
X_5=sol.x;
end
end
end
end
end
end

请先登录，再进行评论。

请先登录，再回答此问题。

采纳的回答

Torsten
Torsten 2023-4-22
编辑：Torsten 2023-4-23
Ran in:
As I said before: Remove x=0 from your interval of integration to avoid a division-by-zero.
Nonetheless, bvp4c has problems with your equations. You should check for programming errors or errors in the mathematical formulation of the problem.
%define space from 0 to 2
x=linspace(1e-6,2,11);
GrT=0;
%initial guess for the answer
guess= bvpinit(x,[1;0;1;0;0;1;1]);
%setting options
options = bvpset('Stats','on','RelTol',10e-4,'AbsTol',10e-7);
hold on
for i = 1:5
GrT=GrT+1;
%calling bvp4c to solve the boundary value problem
sol=bvp4c(@(x,y)bvp_func_pw1(x,y,GrT),@boundary_conditions_pw,guess,options);
%export solution
f = sol.y;
if GrT==1
f_1=f(1,:);
f_11=f(2,:);
X_1=sol.x;
elseif GrT==2
f_2=f(1,:);
f_22=f(2,:);
X_2=sol.x;
elseif GrT==3
f_3=f(1,:);
f_33=f(2,:);
X_3=sol.x;
elseif GrT==4
f_4=f(1,:);
f_44=f(2,:);
X_4=sol.x;
elseif GrT==5
f_5=f(1,:);
f_55=f(2,:);
X_5=sol.x;
end
end
Warning: Unable to meet the tolerance without using more than 1428 mesh points.
The last mesh of 803 points and the solution are available in the output argument.
The maximum residual is 356.632, while requested accuracy is 0.001.
The solution was obtained on a mesh of 1940 points. The maximum residual is 3.566e+02. There were 71198 calls to the ODE function. There were 481 calls to the BC function.
Warning: Unable to meet the tolerance without using more than 1428 mesh points.
The last mesh of 787 points and the solution are available in the output argument.
The maximum residual is 1918.14, while requested accuracy is 0.001.
The solution was obtained on a mesh of 2152 points. The maximum residual is 1.918e+03. There were 64070 calls to the ODE function. There were 474 calls to the BC function.
Warning: Unable to meet the tolerance without using more than 1428 mesh points.
The last mesh of 735 points and the solution are available in the output argument.
The maximum residual is 370.915, while requested accuracy is 0.001.
The solution was obtained on a mesh of 1683 points. The maximum residual is 3.709e+02. There were 61149 calls to the ODE function. There were 452 calls to the BC function.
Warning: Unable to meet the tolerance without using more than 1428 mesh points.
The last mesh of 808 points and the solution are available in the output argument.
The maximum residual is 542.213, while requested accuracy is 0.001.
The solution was obtained on a mesh of 1991 points. The maximum residual is 5.422e+02. There were 56790 calls to the ODE function. There were 467 calls to the BC function.
Warning: Unable to meet the tolerance without using more than 1428 mesh points.
The last mesh of 793 points and the solution are available in the output argument.
The maximum residual is 271.252, while requested accuracy is 0.001.
The solution was obtained on a mesh of 2134 points. The maximum residual is 2.713e+02. There were 46605 calls to the ODE function. There were 424 calls to the BC function.
plot(X_1,f_1,X_2,f_2,X_3,f_3,X_4,f_4,X_5,f_5)
grid on
function dYdx = bvp_func_pw1(x, Y, GrT)
%defining the constants
lambda=0.4; %second grade fluid parameter
M=5; %magnetic parameter
P=5; %porous parameter
%GrT=[1:5]; %thermal grashof numbers
GrC=1; %solatal grashof numbers
Pr=1; %prandtil number
Ec=0.1; %Eckert number
Sc=1; %Schimdit number
Kr=-0.3; %Chemical reaction paramter
alpha=2; %Heat generator or absorption paramter
%defining derivatives
%Y(1)=f0, Y(2)=f1, Y(3)=f2,Y(4)=theta, Y(5)=theta1,Y(6)=phi,
%Y(7)=phi1
dYdx(1)=Y(2);
dYdx(2)=Y(3);
dYdx(3)=(1/(2*x*Y(3)))*(Y(2)^2-Y(1)*Y(3)+M*Y(2)+P*Y(2)-GrT*Y(4)-GrC*Y(6));
dYdx(4)=Y(5);
dYdx(5)=-Pr*Y(1)*Y(5)+Pr*(2*Y(2)-alpha)*Y(4)-Pr*Ec*(Y(3)^2+lambda*Y(3)*(Y(2)*Y(3)-Y(1)*dYdx(3))+M*Y(2)^2);
dYdx(6)=Y(7);
dYdx(7)=-Sc*(Y(1)*Y(7)-Kr*Y(6));
dYdx=dYdx(:);
end
function bc = boundary_conditions_pw(Ya,Yb)
%boundary conditions at x=0 and x=infinity
%boundary conditions at x=0
bc1=[Ya(1); Ya(2)-1; Ya(4)-1; Ya(6)-1];
%boundary conditions at x=infinity
bc2=[Yb(2); Yb(4); Yb(6)];
bc=[bc1;bc2];
end

更多回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Hypothesis Tests 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by