diff --git a/docs/source/code.rst b/docs/source/code.rst
index a81bf36c..dde71fc5 100644
--- a/docs/source/code.rst
+++ b/docs/source/code.rst
@@ -15,4 +15,5 @@ Code Documentation
    vtkhandler
    unvhandler
    igeshandler
+   utilities
 
diff --git a/docs/source/utilities.rst b/docs/source/utilities.rst
new file mode 100644
index 00000000..09c21cf2
--- /dev/null
+++ b/docs/source/utilities.rst
@@ -0,0 +1,5 @@
+Utilities
+=====================
+
+.. automodule:: pygem.utilities
+	:members:
diff --git a/pygem/__init__.py b/pygem/__init__.py
index 5d940317..cbb0264a 100644
--- a/pygem/__init__.py
+++ b/pygem/__init__.py
@@ -1,5 +1,5 @@
 
-__all__ = ['affine', 'filehandler', 'freeform', 'openfhandler', 'params', 'stlhandler', 'unvhandler', 'vtkhandler', 'igeshandler']
+__all__ = ['affine', 'filehandler', 'freeform', 'openfhandler', 'params', 'stlhandler', 'unvhandler', 'vtkhandler', 'igeshandler', 'utilities']
 
 from . import affine
 from . import freeform
@@ -10,3 +10,4 @@
 from . import unvhandler
 from . import vtkhandler
 from . import igeshandler
+from . import utilities
diff --git a/pygem/utilities.py b/pygem/utilities.py
new file mode 100644
index 00000000..5b1afbbc
--- /dev/null
+++ b/pygem/utilities.py
@@ -0,0 +1,86 @@
+"""
+Auxiliary utilities for PyGeM.
+"""
+
+import vtk
+import vtk.util.numpy_support as ns
+import numpy as np
+
+# TODO: add the connectivity to the ffd control points to visualize the lattice.
+
+
+
+def write_bounding_box(parameters, outfile, write_deformed=True):
+	"""
+	Method that writes a vtk file containing the FFD lattice. This method
+	allows to visualize where the FFD control points are located before the geometrical morphing.
+	If the flag is set to original (default) the method writes out the undeformed lattice, 
+	if it is set to modified it writes out the deformed lattice.
+	
+	:param FFDParameters parameters: parameters of the Free Form Deformation.
+	:param string outfile: name of the output file.
+	:param bool write_deformed: flag to write the original or modified FFD control lattice. 
+		The default is set to deformed.
+		
+	:Example:
+	
+	>>> import pygem.utilities as util
+	>>> import pygem.params as pars
+	>>> import numpy as np
+	
+	>>> params = pars.FFDParameters()
+	>>> params.read_parameters(filename='tests/test_datasets/parameters_test_ffd_sphere.prm')
+	>>> util.write_bounding_box(params, 'tests/test_datasets/box_test_sphere.vtk')
+	"""
+
+	aux_x = np.linspace(0, parameters.lenght_box_x, parameters.n_control_points[0])
+	aux_y = np.linspace(0, parameters.lenght_box_y, parameters.n_control_points[1])
+	aux_z = np.linspace(0, parameters.lenght_box_z, parameters.n_control_points[2])
+	lattice_y_coords, lattice_x_coords, lattice_z_coords = np.meshgrid(aux_y, aux_x, aux_z)
+
+	if write_deformed == False:
+		box_points = np.array([lattice_x_coords.ravel(), lattice_y_coords.ravel(), lattice_z_coords.ravel()])
+	if write_deformed == True:
+		box_points = np.array([lattice_x_coords.ravel() + parameters.array_mu_x.ravel()*parameters.lenght_box_x, \
+		lattice_y_coords.ravel()  + parameters.array_mu_y.ravel()*parameters.lenght_box_y, \
+		lattice_z_coords.ravel()  + parameters.array_mu_z.ravel()*parameters.lenght_box_z])
+		
+	n_rows = box_points.shape[1]
+
+	box_points = np.dot(parameters.rotation_matrix,box_points) + np.transpose(np.tile(parameters.origin_box, (n_rows,1)))
+	
+	_write_vtk_box(box_points, outfile)
+	
+	
+def _write_vtk_box(box_points, filename):
+	"""
+	Method that writes a vtk file containing FFD control points.
+	
+	:param numpy.ndarray box_points: coordinates of the FFD control points.
+	:param string filename: name of the output file.
+	"""
+	# setup points and vertices
+	points = vtk.vtkPoints()
+	vertices = vtk.vtkCellArray()
+	
+	for index in range(0, box_points.shape[1]):
+		id = points.InsertNextPoint(box_points[0, index], box_points[1, index], box_points[2, index])
+		vertices.InsertNextCell(1)
+		vertices.InsertCellPoint(id)
+		
+	polydata = vtk.vtkPolyData()
+	polydata.SetPoints(points)
+	polydata.SetVerts(vertices)
+	
+	polydata.Modified()
+	writer = vtk.vtkDataSetWriter()
+	writer.SetFileName(filename)
+	
+	if vtk.VTK_MAJOR_VERSION <= 5:
+		polydata.Update()
+		writer.SetInput(polydata)
+	else:
+		writer.SetInputData(polydata)
+		
+	writer.Write()
+
diff --git a/tests/test_datasets/box_modified_test_sphere_true_version5.vtk b/tests/test_datasets/box_modified_test_sphere_true_version5.vtk
new file mode 100644
index 00000000..9b0cb110
--- /dev/null
+++ b/tests/test_datasets/box_modified_test_sphere_true_version5.vtk
@@ -0,0 +1,24 @@
+# vtk DataFile Version 3.0
+vtk output
+ASCII
+DATASET POLYDATA
+POINTS 12 float
+-20 -55 -45 -4.37166 -55 43.6327 -42.9398 31.9333 -40.9551 
+-27.3115 31.9333 47.6776 1.40315 -49.1766 -48.774 17.0315 -49.1766 39.8587 
+-21.5367 37.7568 -44.729 -5.90834 37.7568 43.9037 52.8929 -77.5023 -149.241 
+38.4346 -43.3531 36.0848 -0.133524 43.5802 -48.503 15.4948 43.5802 40.1297 
+
+VERTICES 12 24
+1 0 
+1 1 
+1 2 
+1 3 
+1 4 
+1 5 
+1 6 
+1 7 
+1 8 
+1 9 
+1 10 
+1 11 
+
diff --git a/tests/test_datasets/box_modified_test_sphere_true_version6.vtk b/tests/test_datasets/box_modified_test_sphere_true_version6.vtk
new file mode 100644
index 00000000..c547733b
--- /dev/null
+++ b/tests/test_datasets/box_modified_test_sphere_true_version6.vtk
@@ -0,0 +1,24 @@
+# vtk DataFile Version 4.0
+vtk output
+ASCII
+DATASET POLYDATA
+POINTS 12 float
+-20 -55 -45 -4.37166 -55 43.6327 -42.9398 31.9333 -40.9551 
+-27.3115 31.9333 47.6776 1.40315 -49.1766 -48.774 17.0315 -49.1766 39.8587 
+-21.5367 37.7568 -44.729 -5.90834 37.7568 43.9037 52.8929 -77.5023 -149.241 
+38.4346 -43.3531 36.0848 -0.133524 43.5802 -48.503 15.4948 43.5802 40.1297 
+
+VERTICES 12 24
+1 0 
+1 1 
+1 2 
+1 3 
+1 4 
+1 5 
+1 6 
+1 7 
+1 8 
+1 9 
+1 10 
+1 11 
+
diff --git a/tests/test_datasets/box_test_sphere_true_version5.vtk b/tests/test_datasets/box_test_sphere_true_version5.vtk
new file mode 100644
index 00000000..f0d478ff
--- /dev/null
+++ b/tests/test_datasets/box_test_sphere_true_version5.vtk
@@ -0,0 +1,24 @@
+# vtk DataFile Version 3.0
+vtk output
+ASCII
+DATASET POLYDATA
+POINTS 12 float
+-20 -55 -45 -4.37166 -55 43.6327 -42.9398 31.9333 -40.9551 
+-27.3115 31.9333 47.6776 1.40315 -49.1766 -48.774 17.0315 -49.1766 39.8587 
+-21.5367 37.7568 -44.729 -5.90834 37.7568 43.9037 22.8063 -43.3531 -52.5479 
+38.4346 -43.3531 36.0848 -0.133524 43.5802 -48.503 15.4948 43.5802 40.1297 
+
+VERTICES 12 24
+1 0 
+1 1 
+1 2 
+1 3 
+1 4 
+1 5 
+1 6 
+1 7 
+1 8 
+1 9 
+1 10 
+1 11 
+
diff --git a/tests/test_datasets/box_test_sphere_true_version6.vtk b/tests/test_datasets/box_test_sphere_true_version6.vtk
new file mode 100644
index 00000000..cc40df95
--- /dev/null
+++ b/tests/test_datasets/box_test_sphere_true_version6.vtk
@@ -0,0 +1,24 @@
+# vtk DataFile Version 4.0
+vtk output
+ASCII
+DATASET POLYDATA
+POINTS 12 float
+-20 -55 -45 -4.37166 -55 43.6327 -42.9398 31.9333 -40.9551 
+-27.3115 31.9333 47.6776 1.40315 -49.1766 -48.774 17.0315 -49.1766 39.8587 
+-21.5367 37.7568 -44.729 -5.90834 37.7568 43.9037 22.8063 -43.3531 -52.5479 
+38.4346 -43.3531 36.0848 -0.133524 43.5802 -48.503 15.4948 43.5802 40.1297 
+
+VERTICES 12 24
+1 0 
+1 1 
+1 2 
+1 3 
+1 4 
+1 5 
+1 6 
+1 7 
+1 8 
+1 9 
+1 10 
+1 11 
+
diff --git a/tests/test_utilities.py b/tests/test_utilities.py
new file mode 100644
index 00000000..8d278284
--- /dev/null
+++ b/tests/test_utilities.py
@@ -0,0 +1,68 @@
+
+from unittest import TestCase
+import unittest
+import pygem.utilities as util
+import pygem.params as pars
+import numpy as np
+import filecmp
+import os
+
+
+class TestVtkHandler(TestCase):
+
+	
+	def test_utilities_write_original_box(self):
+		params = pars.FFDParameters()
+		params.read_parameters(filename='tests/test_datasets/parameters_test_ffd_sphere.prm')
+		
+		outfilename = 'tests/test_datasets/box_test_sphere.vtk'
+		
+		util.write_bounding_box(params, outfilename, False)
+		os.remove(outfilename)
+		
+	
+	def test_utilities_write_modified_box(self):
+		params = pars.FFDParameters()
+		params.read_parameters(filename='tests/test_datasets/parameters_test_ffd_sphere.prm')
+		
+		outfilename = 'tests/test_datasets/box_test_sphere.vtk'
+		
+		util.write_bounding_box(params, outfilename)
+		os.remove(outfilename)
+
+
+	def test_utilities_check_vtk_original_box(self):
+		import vtk
+		
+		params = pars.FFDParameters()
+		params.read_parameters(filename='tests/test_datasets/parameters_test_ffd_sphere.prm')
+		
+		outfilename = 'tests/test_datasets/box_test_sphere.vtk'
+		if vtk.VTK_MAJOR_VERSION <= 5:
+			outfilename_expected = 'tests/test_datasets/box_test_sphere_true_version5.vtk'
+		else:
+			outfilename_expected = 'tests/test_datasets/box_test_sphere_true_version6.vtk'
+		
+		util.write_bounding_box(params, outfilename, False)
+		
+		self.assertTrue(filecmp.cmp(outfilename, outfilename_expected))
+		os.remove(outfilename)
+		
+		
+	def test_utilities_check_vtk_modified_box(self):
+		import vtk
+		
+		params = pars.FFDParameters()
+		params.read_parameters(filename='tests/test_datasets/parameters_test_ffd_sphere.prm')
+		
+		outfilename = 'tests/test_datasets/box_test_sphere.vtk'
+		if vtk.VTK_MAJOR_VERSION <= 5:
+			outfilename_expected = 'tests/test_datasets/box_modified_test_sphere_true_version5.vtk'
+		else:
+			outfilename_expected = 'tests/test_datasets/box_modified_test_sphere_true_version6.vtk'
+		
+		util.write_bounding_box(params, outfilename)
+		
+		self.assertTrue(filecmp.cmp(outfilename, outfilename_expected))
+		os.remove(outfilename)
+