Singular Jacobian Error Message

Hi all, I'm running the following code:

function M = nonNewtonian(~) 

M = bvpinit(linspace(0,10,301),@VKinit);
sol = bvp4c(@VK,@VKbc,M);

figure;
hold all;
plot(sol.x,sol.y(2,:)); 
plot(sol.x,sol.y(4,:)); 
hold off;
xlabel('\zeta')
xlabel('\zeta')
hleg = legend('F\prime','G\prime',...
'Location','NorthEast'); %#ok<NASGU>

figure;
hold all;
plot(sol.x,sol.y(1,:)); 
plot(sol.x,sol.y(3,:)); 
plot(sol.x,(-1)*sol.y(5,:)); 
hold off;
xlabel('\zeta')
hleg = legend('F','G','-H',...
'Location','East'); %#ok<NASGU>

function yprime = VK(x,y)

n=1;

yprime = [ y(2) 
n^(-1)*((y(2)^(2)+y(4)^(2))^((n-1)/2))^(-1)*((y(1)^(2)-y(3)^(2)+(y(5)+((1-n)/(n+1))*y(1)*x)*y(2))*(1+(n-1)*(y(2)^(2)+y(4)^2)^(-1)*y(4)^(2))-(n-1)*y(2)*y(4)*(y(2)^(2)+y(4)^(2))^(-1)*(2*y(1)*y(3)+(y(5)+((1-n)/(n+1))*y(1)*x)*y(4)))
y(4)
n^(-1)*((y(2)^(2)+y(4)^(2))^((n-1)/2))^(-1)*((2*y(1)*y(3)+(y(5)+((1-n)/(n+1))*y(1)*x)*y(4))*(1+(n-1)*(y(2)^(2)+y(4)^2)^(-1)*y(2)^(2))-(n-1)*y(2)*y(4)*(y(2)^(2)+y(4)^(2))^(-1)*(y(1)^(2)-y(3)^(2)+(y(5)+((1-n)/(n+1))*y(1)*x)*y(2)))
-2*y(1)-(1-n)/(n+1)*x*y(2)];

function res = VKbc(ya,yb)

res = [ya(1);ya(3)-1;ya(5);yb(2)-(yb(5)*yb(1));yb(4)-(yb(5)*yb(3))];

function yinit = VKinit(~)

yinit = [0;0;1;0;0];

but receive the following error message:

??? Error using ==> bvp4c at 252
Unable to solve the collocation equations -- a singular Jacobian encountered.

Error in ==> nonNewtonian at 4
sol = bvp4c(@VK,@VKbc,M);

I'm using n=1 as a test case here. I know the solutions to this system for n=1 but would like to look into the solutions when n is not equal to one. Can anyone explain why I am getting this error?

Any help anyone could give would be greatly appreciated. Thanks!

Answer 1

[Star Strider]

The most likely reason you got the error was that the your boundary values are too large.

I suggest you experiment with:

N = 1;
M = bvpinit(linspace(0,N,301),@VKinit);

and increase N until you again get the error. A value of N that is <10 may not be what you want, but it may be what you have to live with, given finite machine precision. (I did not try to run your code because I cannot determine where the rows of your yprime matrix end. It always helps to format your code. Highlight your code and use the ‘{} Code’ button at the top of the text window to format it.)

Answer 2

[pxg882]

Thanks for taking the time to look at this. I didn't realise how horribly my code had been displayed on here.

Changing the boundary values is useful you were right. The problem was in the initial conditions I was using.

Thanks for the advice.

[Day Rosli]

Hello, Hi pxg882. Do you manage to solve your problem? I have the same problem like yours. Is it just by decreasing the value of N?

cc: Star Strider