% Create tau objects.
[x, y] = tau('LegendreP', [0 60], 10);

% Specify problem, conditions and exact solution.
ode = {'diff(y1) -  10*y2  +  10*y1               =  0'; ...
       'diff(y2) -  28*y1  + yo3*y1 + yo1*y3 + y2 =  yo1*yo3'; ...
       'diff(y3) + 8/3*y3  - yo2*y1 - yo1*y2      = -yo1*yo2'};

conditions = {'y1(0)=-8';'y2(0)=-7';'y3(0)=29'};

% Solve the problem.
a = tausysnewton(x, y, .........................................................% Tau variables.
                 ode, ....................% Linearized system of rdinary differential equations.
                 conditions, ......................................................% Conditions.
                 'pieces', 200, .............................................% Number of pieces.
                 'iter', 1, ...........................................% Show newton iterations.
                 'phase', 1, .............................................% Plot ode phase plan.
                 'step', 0.001); ....................................% Step for results ploting.