function [y1,T,y0] = prob8(t,y,jaco)
%
% [alpha,T,y0] = prob8
%
% or
%            f = prob8(t,y)
%
% or
%           f' = prob8(t,y,1)
%
% or
%            y = prob8(t)
%
% Problem 8  in:
%
%  L.Brugnano, G.Gurioli, F.Iavernaro, M.Vikerpuur. 
%  FDE-testset: comparing Matlab© codes for solving fractional
%  differential equations of Caputo type
%  Fractal Fract. 2025, 9(5), 312; https://doi.org/10.3390/fractalfract9050312
%
% 
alfa     = 1.25;
if nargin==0    % problem data
    y1   = alfa;
    T    = 1;
    y0   = [0 0;0 0];
elseif nargin==1 % solution
    y1   = [t(:).^(3+alfa) t(:).^(4+alfa)];
elseif nargin==2 % vector field
%    y1   = [gamma(4+alfa)/6*t.^3-t.^(8+2*alfa)+y(:,2).^2   ...
%            gamma(5+alfa)/4*t.^4+t.^(3+alfa)-y(:,1)];
   if length(t)>1    

      y1 = [gamma(4+alfa)/6*t.^3-t.^(8+2*alfa)+y(:,2).^2   ...
            gamma(5+alfa)/24*t.^4+t.^(3+alfa)-y(:,1)];
    
   else
   
      y1 = [gamma(4+alfa)/6*t^3-t^(8+2*alfa)+y(2)^2;   
            gamma(5+alfa)/24*t^4+t^(3+alfa)-y(1)];    
       
   end    
else             % Jacobian
    y1   = [0 2*y(2); -1 0]; 
end
return

           