diff --git a/report/main.tex b/report/main.tex
index 4761044..9568085 100644
--- a/report/main.tex
+++ b/report/main.tex
@@ -26,6 +26,15 @@ Scopo del progetto è verificare numericamente i risultati
 
 Nell'analisi consideriamo i grafi di Erdo''s-Reiny (Figura)
 
+Consideriamo un grafo non diretto e pesato $ G = (V, E, W)$.
+
+Studiamo i grafi di Erdos-Reiny e di tipo Sensors. (Metti figure)
+
+\section{Esperimento 1}
+
+Di seguito consideriamo 
+
+
 \printbibliography
 
 \end{document}
diff --git a/report/references.bib b/report/references.bib
index 735d5f5..e69de29 100644
--- a/report/references.bib
+++ b/report/references.bib
@@ -1,7 +0,0 @@
-@book{bringhurst2004,
-  author    = {Robert Bringhurst},
-  title     = {The Elements of Typographic Style},
-  year      = {2004},
-  publisher = {Hartley \& Marks},
-  address   = {Vancouver}
-}
diff --git a/src/afgl/main.py b/src/afgl/main.py
index a54cbc6..578a77c 100644
--- a/src/afgl/main.py
+++ b/src/afgl/main.py
@@ -3,12 +3,22 @@ import sys
 
 import matplotlib.pyplot as plt
 import numpy as np
+import numpy.linalg as LA
+
+# Requires latex installed
 from pygsp import graphs, plotting
 
+from afgl.util.build_T_matrix import build_T_matrix
+from afgl.util.lanczos import lanczos
 
-def test_plot():
+
+def plot_setup():
     plotting.BACKEND = "matplotlib"
     plt.style.use(["science"])
+    # TODO match font with document
+
+
+def test_plot():
     fig, ax = plt.subplots(figsize=(3.3, 2.5))
     # G = graphs.ErdosRenyi()
     G = graphs.Sensor()
@@ -23,8 +33,50 @@ def test_plot():
     plt.savefig("./out/test.pdf", bbox_inches="tight")
 
 
+def g_extended(t):
+    return np.sin(0.5 * np.pi * np.cos(np.pi * t) ** 2)
+
+
+"""
+Evaluates the function sin(0.5*pi*cos(pi*t)^2)chi_[-1/2,1/2] where chi_I is the
+characteristic function of I, as defined in example 1 of the paper.
+"""
+
+
+def g(T):
+    Ind = ((T >= -1 / 2) & (T <= 1 / 2)).astype(int)
+    Val = g_extended(T)
+    # Perform element wise multiplication to resemble applying indicating
+    # function to all matrix.
+    return np.multiply(Val, Ind)
+
+
+def example_1():
+    N = 500
+    M = 300
+    p = 0.04
+
+    G = graphs.ErdosRenyi(N, p)
+
+    G.compute_laplacian("combinatorial")
+    L = G.L
+
+    s = np.random.randint(1, 10000, N)
+
+    [V, alp, beta] = lanczos(L, s, M)
+
+    T = build_T_matrix(alp, beta)
+
+    e_1 = np.zeros(M)
+    e_1[0] = 1
+    y = LA.norm(s) * (g(T) * e_1)
+    G_L_s = V @ y
+    return G_L_s
+
+
 def run():
-    test_plot()
+    plot_setup()
+    example_1()
 
 
 if __name__ == "__main__":
diff --git a/src/afgl/arnoldi.py b/src/afgl/util/arnoldi.py
similarity index 100%
rename from src/afgl/arnoldi.py
rename to src/afgl/util/arnoldi.py
diff --git a/src/afgl/util/build_T_matrix.py b/src/afgl/util/build_T_matrix.py
new file mode 100644
index 0000000..d274384
--- /dev/null
+++ b/src/afgl/util/build_T_matrix.py
@@ -0,0 +1,5 @@
+import numpy as np
+
+
+def build_T_matrix(alp, beta):
+    return np.diag(alp) + np.diag(beta, -1) + np.diag(beta, 1)
diff --git a/src/afgl/lanczos.py b/src/afgl/util/lanczos.py
similarity index 83%
rename from src/afgl/lanczos.py
rename to src/afgl/util/lanczos.py
index 3797a16..82d5cea 100644
--- a/src/afgl/lanczos.py
+++ b/src/afgl/util/lanczos.py
@@ -2,10 +2,12 @@ import numpy as np
 import numpy.linalg as LA
 
 """
+Classic Lanczos method (without re-orthogonalization)
+
 Arguments
-L Real valued NxN symmetric matrix
-s vector of size N
-M natural number indicating basis size
+L : Real valued NxN symmetric matrix
+s : vector of size N
+M : natural number indicating basis size
 
 Returns
 -------
diff --git a/tests/test_main.py b/tests/test_main.py
index d7270dd..3c69326 100644
--- a/tests/test_main.py
+++ b/tests/test_main.py
@@ -1,6 +1,7 @@
 import numpy as np
 import numpy.linalg as LA
-from afgl.lanczos import lanczos
+from afgl.util.build_T_matrix import build_T_matrix
+from afgl.util.lanczos import lanczos
 
 """
 Todo: better test case
@@ -18,7 +19,7 @@ def test_lanczos_return_correct_solution():
     s = np.random.randint(1, 10, N)
     [V, alp, beta] = lanczos(L, s, M)
 
-    T = np.diag(alp) + np.diag(beta, -1) + np.diag(beta, 1)
+    T = build_T_matrix(alp, beta)
 
     x = LA.solve(L, s)
     e_1 = np.zeros(M)