Краевая задача

[6AYMAHEU,]

Помогите пожалуйста, проект завис. Не могу решить систему диффуров.

Сама прога в MatLab
S = [0,0,0;0,0,0;0,0,0];
options = bvpset('SingularTerm',S);
guess = [1;1;1];
solinit = bvpinit(linspace(0,1,20),guess);
sol = bvp4c(@emdenode,@emdenbc,solinit,options);
x = linspace(0,1);
plot(sol.x,sol.y(1,:));
Функции m-файлы
function dydx = emdenode(x,y)
dydx = [ -1000^2*x+ tan(y(3))/cos(y(2))
(1000^2*x*sin(y(2)) - 2*1000*10*cos(y(3)))/y(1)/cos(y(3))^2/cos(y(2)) - tan(y(2))/(x + 0.00001)
tan(y(3))/y(1)*1000^2*x + 1/y(1)/cos(y(2))
];
И
function res = emdenbc(ya,yb)
res = [ ya(1) - 10
ya(3) - 1
yb(3) - 1.5
];

Решение расходится. Может другие алгоритмы нужны? Помогите, если есть похожие программы в MatLab,Mathematic,Maple,C++.