diff --git a/.DS_Store b/.DS_Store
index 054e7c7..669318a 100644
Binary files a/.DS_Store and b/.DS_Store differ
diff --git a/Exercise_03/calore_CN.m b/Exercise_03/calore_CN.m
new file mode 100644
index 0000000..2cb2331
--- /dev/null
+++ b/Exercise_03/calore_CN.m
@@ -0,0 +1,60 @@
+function [x,t,u] = calore_CN(L,N,T,K,c1,c2,f,u0)
+% ----  Risoluzione dell'equazione del calore ----
+%   u_t - u_xx = f nell'intervallo [-L,L] 
+%   con condizioni al bordo di Dirichlet
+%   e condizioni iniziali.
+% -----------------------------------------------
+% Sintassi:
+%   [x,t,u]=calore_template(L,N,T,K,c1,c2,fun,u0)
+%
+% Input:
+%   L semiampiezza intervallo spaziale (-L,L) 
+%	N numero di sottointervalli in (-L,L)
+%   T estremo finale intervallo temporale (0,T)
+%	K numero di sottointervalli in (0,T)
+%   c1 funzione che descrive la condizione di Dirichlet in x=-L
+%   c2 funzione che descrive la condizione di Dirichlet in x=L
+%   f funzione che descrive il termine noto dell'equazione
+%   u0 funzione che descrive la condizione iniziale in t=0
+%
+% Output:
+%   x vettore dei nodi spaziali
+%   t vettore dei nodi temporali
+%   u soluzione numerica 
+%     del problema
+
+% Calcolo passo di discretizzazione in spazio e tempo
+h=2*L/N;
+tau=T/K;
+% Inizializzazione del vettore t
+t=linspace(0,T,K+1)';
+% Inizializzazione del vettore x
+x=linspace(-L,L,N+1);
+
+% Inizializzazione della matrice soluzione u
+u=zeros(N+1,K+1);
+% Condizione iniziale
+u(:,1)=u0(x);
+
+% Condizioni al bordo
+u(1,:)=c1(t);
+u(end,:)=c2(t);
+
+% Costruzione della matrice A
+e=ones(N-1,1);
+A=spdiags([-e,2*e,-e],[-1,0,1],N-1,N-1)/(h^2);
+I=speye(N-1,N-1);
+
+% Ciclo iterativo
+for k=1:K
+    % Assemblaggio termine noto
+    F1=f(x(2:end-1),t(k));
+    F2=f(x(2:end-1),t(k+1))
+    % Correzione del termine noto con le condizioni al bordo
+    F1(1)=F1(1) + c1(t(k))/(h^2);
+    F1(end)=F1(end) + c2(t(k))/(h^2);
+    F2(1)=F2(1) + c1(t(k+1))/(h^2);
+    F2(end)=F2(end) + c2(t(k+1))/(h^2);
+    % Risoluzione del problema
+    u(2:end-1,k+1) = ((I + 0.5*tau*A)\((I - 0.5*tau*A)*u(2:end-1,k)))' + 0.5*tau*F1 + 0.5*tau*F2;
+end
diff --git a/Exercise_03/image.png b/Exercise_03/image.png
new file mode 100644
index 0000000..a820c4f
--- /dev/null
+++ b/Exercise_03/image.png
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:33d48866166fd1a5e429d5b52f83a55a1ce3f166656cb0a55754f10d7d38026f
+size 25398