% Queremos resolver c*F''(x) + F'(x) =1 en x=(0,1)
m=input('Ingrese el numero de nodos que desea [mayor a 1 por favor]')
% Establecemos el tamaño de la grilla y el paso
h=1/m;
% Alocamos espacio para el vector de soluciones en cada punto de la
% particion, y tambien para la matriz, asi como el vector de constantes
%F=zeros(n,1); 
n=m+1;
A=zeros(n,n);
b=ones(n,1);
% Defino las constantes del problema discreto
c=-10^(-4);
c1=(c/(h^2) - 1/(2*h));
c2=-(2*c)/(h^2);
c3=(c/(h^2) + 1/(2*h));
% Completo la matriz correspondiente que es tridiagonal
for i=1:n
A(i,i)=c2;
end
for i=2:n
A(i,i-1)=c1;
end
for i=1:(n-1)
A(i,i+1)=c3;
end
A;
% Queremos resolver entonces AF=b, entonces hacemos LU
[L,U,P]=lu(A);
% Tenemos LUF=Pb
b=P*b;
% Resolvemos ahora el nuevo sistema LY=b
Y=zeros(n,1);
Y(1)=b(1);
for i=2:n
    Y(i)=b(i)-L(i,i-1)*Y(i-1);
end
%Y(n)=0;
Y;
% Ahora resolvemos el sistemita UF=Y
F=zeros(n,1);
F(n)=Y(n)/U(n,n);
for i=(n-1):-1:1
   F(i)=( Y(i) - U(i,i+1)*F(i+1) ) / U(i,i);
end
F(n)=0;
F(1)=0;
X=zeros(n,1);
X(1)=0;
for i=2:n
    X(i)=i*h;
end
X=transpose(X);
F
close
plot (X,F,'r')