24 views (last 30 days)

Show older comments

J on 29 Jul 2024 at 8:56

Commented: Torsten on 29 Jul 2024 at 15:23

Accepted Answer: Torsten

- W2.m

Open in MATLAB Online

I tried to solve the bvp using bvp4c procedure. But I have an error "Unable to solve the collocation equations -- a singular Jacobian encountered." How to resolve this issue?

proj()

Error using bvp4c (line 196)

Unable to solve the collocation equations -- a singular Jacobian encountered.

Error in solution>proj (line 31)

sol= bvp4c(@projfun,@projbc,solinit,options);

function sol= proj

A = 0.5;

Gr = 0.7;

Gc = 0.5;

Kp = 3.0;

beta = 0.5;

Pr = 0.3;

Df = 0.2;

Sc = 0.1;

L0 = 0.5;

Sr = 0.3;

M = 1;

Bi1 = 0.5;

Bi2 = 0.5;

K0 = 0.3;

myLegend1 = {};myLegend2 = {};

rr = [0.3 0.5 0.8];

for i =1:numel(rr)

Re = rr(i);

y0 = [1, 0, 1, 1, 0, 1, 0];

options =bvpset('stats','on','RelTol',1e-4);

x = linspace(0,10,500);

solinit = bvpinit(x,y0);

sol= bvp4c(@projfun,@projbc,solinit,options);

disp((sol.y(1,20)))

figure(1)

plot(sol.x,(sol.y(6,:)))

grid on,hold on

myLegend1{i}=['Pr = ',num2str(rr(i))];

xlabel('eta');

ylabel('(thetas-thetaf)/thetas');

i=i+1;

end

figure(1)

legend(myLegend1)

hold on

function dy= projfun(~,y)

dy= zeros(7,1);

dy(1) = y(2);

dy(2) = y(3);

dy(3) = (0.5*A*y(3) - Gr*y(4) - Gc*y(6) + (M + Kp + A)*y(2)) / (1 + (1/beta));

dy(4) = y(5);

dy(5) = (Pr*(0.5*A*y(5)) - Pr*Df*Sc*(0.5*A*y(7) + 2*A*y(6) - L0*y(6))) / (1 - (Pr*Df*Sc*Sr));

dy(6) = y(7);

dy(7) = (Sc*(0.5*A*y(7) + 2*A*dy(6)) - Sr*Pr*(0.5*A*y(5) + 2*A*y(4))) / (1 - (Sr*Df*Pr));

end

function res= projbc(ya,yb)

res= [

ya(2) - (1 + K0*ya(3));

ya(5) + Bi1*(1 - ya(4));

ya(7) + Bi2*(1 - ya(6));

yb(2);

yb(4);

yb(6);

yb(7)];

end

end

##### 1 Comment Show -1 older commentsHide -1 older comments

Show -1 older commentsHide -1 older comments

Torsten on 29 Jul 2024 at 9:11

#### Direct link to this comment

https://ms-www.mathworks.com/matlabcentral/answers/2141026-how-can-i-solve-the-error-unable-to-solve-the-collocation-equations-a-singular-jacobian-encounte#comment_3222876

How to resolve this issue?

It's not a technical issue. Usually, there is something wrong with the equations or the boundary conditions. So you should compare both of them with the mathematical description of the problem.

Sign in to comment.

Sign in to answer this question.

### Accepted Answer

Torsten on 29 Jul 2024 at 13:53

Edited: Torsten on 29 Jul 2024 at 14:53

Open in MATLAB Online

Your ODE system is linear - thus you can determine its general solution having 7 free parameters. But incorporating of your boundary conditions leads to a linear system of equations where the coefficient matrix A is rank-deficient (rank 6 instead of rank 7). Consequently, your system is not solvable.

A = 0.5;

Gr = 0.7;

Gc = 0.5;

Kp = 3.0;

beta = 0.5;

Pr = 0.3;

Df = 0.2;

Sc = 0.1;

L0 = 0.5;

Sr = 0.3;

M = 1;

Bi1 = 0.5;

Bi2 = 0.5;

K0 = 0.3;

Mat = zeros(7);

Mat(1,2) = 1;

Mat(2,3) = 1;

Mat(3,2) = (M + Kp + A) / (1 + 1/beta);

Mat(3,3) = 0.5*A/(1 + 1/beta);

Mat(3,4) = -Gr/(1 + 1/beta);

Mat(3,6) = -Gc/(1 + 1/beta);

Mat(4,5) = 1;

Mat(5,5) = Pr*0.5*A/(1 - Pr*Df*Sc*Sr);

Mat(5,6) = - Pr*Df*Sc*(2*A - L0) / (1 - Pr*Df*Sc*Sr);

Mat(5,7) = - Pr*Df*Sc*0.5*A/ (1 - Pr*Df*Sc*Sr);

Mat(6,7) = 1;

Mat(7,4) = - Sr*Pr*2*A / (1 - Sr*Df*Pr);

Mat(7,5) = - Sr*Pr*0.5*A / (1 - Sr*Df*Pr);

Mat(7,7) = Sc*(0.5*A+2*A) / (1 - Sr*Df*Pr);

syms x1(t) x2(t) x3(t) x4(t) x5(t) x6(t) x7(t)

eqns = [diff(x1);diff(x2);diff(x3);diff(x4);diff(x5);diff(x6);diff(x7)]-Mat*[x1;x2;x3;x4;x5;x6;x7]==[0;0;0;0;0;0;0];

sol = dsolve(eqns,'MaxDegree',4)

sol = struct with fields:

x2: C3*(exp((t*(2976375*((13875*3^(1/2)*342589119285075140699353^(1/2))/1883723499362492536448 - 2321348855894375/1076413428207138592256)^(1/6) - 9802324*((3364850596725*((... x1: C1 + C2*(exp((t*(2976375*((13875*3^(1/2)*342589119285075140699353^(1/2))/1883723499362492536448 - 2321348855894375/1076413428207138592256)^(1/6) - 9802324*((33648505967... x3: C6*exp(-(t*(865^(1/2) - 1))/24) + C7*exp((t*(865^(1/2) + 1))/24) - C2*(exp((t*(2976375*((13875*3^(1/2)*342589119285075140699353^(1/2))/1883723499362492536448 - 23213488... x4: C2*(exp((t*(2976375*((13875*3^(1/2)*342589119285075140699353^(1/2))/1883723499362492536448 - 2321348855894375/1076413428207138592256)^(1/6) - 9802324*((3364850596725*((... x5: C4*exp((t*(9802324*((3364850596725*((13875*3^(1/2)*342589119285075140699353^(1/2))/1883723499362492536448 - 2321348855894375/1076413428207138592256)^(1/3))/960855558009... x6: C3*(exp((t*(2976375*((13875*3^(1/2)*342589119285075140699353^(1/2))/1883723499362492536448 - 2321348855894375/1076413428207138592256)^(1/6) - 9802324*((3364850596725*((... x7: C4*exp((t*(9802324*((3364850596725*((13875*3^(1/2)*342589119285075140699353^(1/2))/1883723499362492536448 - 2321348855894375/1076413428207138592256)^(1/3))/960855558009...

eqn1 = subs(sol.x2,t,0)-(1+K0*subs(sol.x3,t,0))==0;

eqn2 = subs(sol.x5,t,0)+Bi1*(1-subs(sol.x4,t,0))==0;

eqn3 = subs(sol.x7,t,0)+Bi2*(1-subs(sol.x6,t,0))==0;

eqn4 = subs(sol.x2,t,10)==0;

eqn5 = subs(sol.x4,t,10)==0;

eqn6 = subs(sol.x6,t,10)==0;

eqn7 = subs(sol.x7,t,10)==0;

[A,b] = equationsToMatrix([eqn1,eqn2,eqn3,eqn4,eqn5,eqn6,eqn7])

A=

b=

rank(A)

ans = 6

rank([A,b])

ans = 7

##### 1 Comment Show -1 older commentsHide -1 older comments

Show -1 older commentsHide -1 older comments

Torsten on 29 Jul 2024 at 15:23

#### Direct link to this comment

https://ms-www.mathworks.com/matlabcentral/answers/2141026-how-can-i-solve-the-error-unable-to-solve-the-collocation-equations-a-singular-jacobian-encounte#comment_3223231

I suspect the reason is that you don't have any boundary condition involving y(1).

Sign in to comment.

### More Answers (0)

Sign in to answer this question.

**An Error Occurred**

Unable to complete the action because of changes made to the page. Reload the page to see its updated state.

Select a Web Site

Choose a web site to get translated content where available and see local events and offers. Based on your location, we recommend that you select: .

You can also select a web site from the following list

Americas

- América Latina (Español)
- Canada (English)
- United States (English)

Europe

- Belgium (English)
- Denmark (English)
- Deutschland (Deutsch)
- España (Español)
- Finland (English)
- France (Français)
- Ireland (English)
- Italia (Italiano)
- Luxembourg (English)

- Netherlands (English)
- Norway (English)
- Österreich (Deutsch)
- Portugal (English)
- Sweden (English)
- Switzerland
- Deutsch
- English
- Français

- United Kingdom(English)

Asia Pacific

- Australia (English)
- India (English)
- New Zealand (English)
- 中国
- 日本Japanese (日本語)
- 한국Korean (한국어)

Contact your local office