function values=expl(Smax,M,T,N,r,sig,X)
dS=Smax/M;
dt=T/N;
matf=zeros(M+1,N+1);
vetS=linspace(0,Smax,M+1)';
matf(:,N+1)=max(X-vetS,0);
matf(1,:)=X;
matf(M+1,:)=0;
veti=0:N;
vetj=0:M;
a=0.5*(-r*dt*vetj+sig^2*dt*(vetj.^2))/(1+r*dt);
b=(1-sig^2*dt*(vetj.^2))/(1+r*dt);
c=0.5*(r*dt*vetj+sig^2*dt*(vetj.^2))/(1+r*dt);
for i=N:-1:1
   for j=2:M
      matf(j,i) = a(j)*matf(j-1,i+1) + b(j)*matf(j,i+1)+ ...
         c(j)*matf(j+1,i+1);
      matf(j,i) = max(matf(j,i),X-(j-1)*dS);
   end
end
values = flipud([vetS, matf(:,1)]);