diff --git a/bladex/blade.py b/bladex/blade.py
index e94cb23..488f144 100644
--- a/bladex/blade.py
+++ b/bladex/blade.py
@@ -322,7 +322,7 @@ def rotate(self, deg_angle=None, rad_angle=None, axis='x'):
         in the 3D Euclidean space about the specified axis, which is
         -- by default -- the x axis.
 
-        :math: when the axis of rotation is the x-axis `R(\\theta)` is defined
+        when the axis of rotation is the x-axis :math: `R(\\theta)` is defined
         by:
 
         .. math::
@@ -401,12 +401,12 @@ def rotate(self, deg_angle=None, rad_angle=None, axis='x'):
             self.blade_coordinates_down[i][2] = new_coord_matrix_down[2]
 
     def scale(self, factor):
-        '''
+        """
         Scale the blade coordinates by a specified factor.
 
         :param float factor: scaling factor
+        """
 
-        '''
         scaling_matrix = np.array([factor, 0, 0, 0, factor,
             0, 0, 0, factor]).reshape((3, 3))
 
diff --git a/bladex/customprofile.py b/bladex/customprofile.py
index 7720ef7..23c0e97 100644
--- a/bladex/customprofile.py
+++ b/bladex/customprofile.py
@@ -9,6 +9,7 @@
 
 import numpy as np
 from .profilebase import ProfileBase
+from scipy.optimize import newton
 
 
 class CustomProfile(ProfileBase):
@@ -89,7 +90,7 @@ def _check_parameters(self):
             raise ValueError(
                 'object "thickness_perc" refers to an empty array.')
         if self.chord_len is None:
-            raise ValueError('object "chorf_len" refers to an empty array.')
+            raise ValueError('object "chord_len" refers to an empty array.')
         if self.camber_max is None:
             raise ValueError('object "camber_max" refers to an empty array.')
         if self.thickness_max is None:
@@ -127,13 +128,13 @@ def _check_parameters(self):
                 'thickness_perc and chord_perc must have same shape.')
 
     def _compute_orth_camber_coordinates(self):
-        '''
+        """
         Compute the coordinates of points on upper and lower profile on the
         line orthogonal to the camber line.
 
         :return: x and y coordinates of section points on line orthogonal to
         camber line
-        '''
+        """
         # Compute the angular coefficient of the camber line
         n_pos = self.chord_percentage.shape[0]
         m = np.zeros(n_pos)
@@ -157,7 +158,7 @@ def _compute_orth_camber_coordinates(self):
             xup_tmp[1], xdown_tmp[1] = xup_tmp[2]-1e-16, xdown_tmp[2]-1e-16
             yup_tmp[1], ydown_tmp[1] = yup_tmp[2]-1e-16, ydown_tmp[2]-1e-16
 
-        return xup_tmp, xdown_tmp, yup_tmp, ydown_tmp
+        return [xup_tmp, xdown_tmp, yup_tmp, ydown_tmp]
 
 
 
@@ -169,7 +170,8 @@ def generate_coordinates(self):
         and camber.
         """
 
-        self.xup_coordinates, self.xdown_coordinates, self.yup_coordinates, self.ydown_coordinates = self._compute_orth_camber_coordinates()
+        [self.xup_coordinates, self.xdown_coordinates, self.yup_coordinates,
+            self.ydown_coordinates] = self._compute_orth_camber_coordinates()
 
         self.ydown_coordinates = self.ydown_curve(
                 self.chord_len*(self.chord_percentage).reshape(-1,1)).reshape(
@@ -186,10 +188,10 @@ def generate_coordinates(self):
 
 
     def adimensionalize(self):
-        '''
+        """
         Rescale coordinates of upper and lower profiles of section such that
         coordinates on x axis are between 0 and 1.
-        '''
+        """
         factor = abs(self.xup_coordinates[-1]-self.xup_coordinates[0])
         self.yup_coordinates *= 1/factor
         self.xdown_coordinates *= 1/factor
@@ -197,12 +199,12 @@ def adimensionalize(self):
         self.xup_coordinates *= 1/factor
 
     def generate_parameters(self):
-        '''
+        """
         Method that generates the parameters of a general airfoil profile
         (chord length, chord percentages, camber max, thickness max, camber and
         thickness percentages), starting from the upper and lower
         coordinates of the section profile.
-        '''
+        """
         n_pos = self.xup_coordinates.shape[0]
 
         self.chord_len = abs(np.max(self.xup_coordinates)-
@@ -220,11 +222,12 @@ def generate_parameters(self):
                     self.chord_percentage[i-1])*self.camber_max/self.chord_len
         m_angle = np.arctan(m)
 
-        from scipy.optimize import newton
-
+        # generating temporary profile coordinates orthogonal to the camber
+        # line
         ind_horizontal_camber = (np.sin(m_angle)==0)
         def eq_to_solve(x):
-            spline_curve = self.ydown_curve(x.reshape(-1,1)).reshape(x.shape[0],)
+            spline_curve = self.ydown_curve(x.reshape(-1,1)).reshape(
+                    x.shape[0],)
             line_orth_camber = (camber[~ind_horizontal_camber] +
                     np.cos(m_angle[~ind_horizontal_camber])/
                     np.sin(m_angle[~ind_horizontal_camber])*(self.chord_len*
@@ -235,7 +238,8 @@ def eq_to_solve(x):
         xdown_tmp[~ind_horizontal_camber] = newton(eq_to_solve,
                 xdown_tmp[~ind_horizontal_camber])
         xup_tmp = 2*self.chord_len*self.chord_percentage - xdown_tmp
-        ydown_tmp = self.ydown_curve(xdown_tmp.reshape(-1,1)).reshape(xdown_tmp.shape[0],)
+        ydown_tmp = self.ydown_curve(xdown_tmp.reshape(-1,1)).reshape(
+                xdown_tmp.shape[0],)
         yup_tmp = 2*self.camber_max*self.camber_percentage - ydown_tmp
         if xup_tmp[1]<self.xup_coordinates[0]:
             xup_tmp[1], xdown_tmp[1] = xup_tmp[2], xdown_tmp[2]
diff --git a/bladex/profilebase.py b/bladex/profilebase.py
index 2ffa797..a595a47 100644
--- a/bladex/profilebase.py
+++ b/bladex/profilebase.py
@@ -5,6 +5,7 @@
 import numpy as np
 import matplotlib.pyplot as plt
 from .ndinterpolator import reconstruct_f
+from scipy.interpolate import RBFInterpolator
 
 
 class ProfileBase(object):
@@ -276,7 +277,7 @@ def deform_camber_line(self, percent_change, n_interpolated_points=None):
 
     @property
     def yup_curve(self):
-        '''
+        """
         Return the spline function corresponding to the upper profile
         of the airfoil
 
@@ -285,15 +286,14 @@ def yup_curve(self):
 
         .. todo::
             generalize the interpolation function
-        '''
-        from scipy.interpolate import RBFInterpolator
+        """
         spline = RBFInterpolator(self.xup_coordinates.reshape(-1,1),
                self.yup_coordinates.reshape(-1,1))
         return spline
 
     @property
     def ydown_curve(self):
-        '''
+        """
         Return the spline function corresponding to the lower profile
         of the airfoil
 
@@ -302,8 +302,7 @@ def ydown_curve(self):
 
         .. todo::
             generalize the interpolation function
-        '''
-        from scipy.interpolate import RBFInterpolator
+        """
         spline = RBFInterpolator(self.xdown_coordinates.reshape(-1,1),
                self.ydown_coordinates.reshape(-1,1))
         return spline
@@ -537,7 +536,7 @@ def scale(self, factor, translate=True):
 
         :param float factor: the scaling factor
         """
-        if translate==True:
+        if translate:
             ref_point = self.reference_point
             self.translate(-ref_point)
             self.xup_coordinates *= factor
diff --git a/bladex/reversepropeller.py b/bladex/reversepropeller.py
index 05cce5f..f735130 100644
--- a/bladex/reversepropeller.py
+++ b/bladex/reversepropeller.py
@@ -84,7 +84,7 @@ def optimized_path(coords, start=None):
 
 
 def point_inside_polygon(x, y, poly, include_edges=True):
-    '''
+    """
     Test if point (x,y) is inside polygon poly.
 
     poly is N-vertices polygon defined as
@@ -94,7 +94,7 @@ def point_inside_polygon(x, y, poly, include_edges=True):
     Geometrical idea: point is inside polygon if horisontal beam
     to the right from point crosses polygon even number of times.
     Works fine for non-convex polygons.
-    '''
+    """
     n = len(poly)
     inside = False
 
diff --git a/tests/test_ndinterpolator.py b/tests/test_ndinterpolator.py
index d8f6e94..1fd2f11 100644
--- a/tests/test_ndinterpolator.py
+++ b/tests/test_ndinterpolator.py
@@ -44,39 +44,39 @@ def test_wrong_basis(self):
     def test_gaussian_evaluation(self):
         rbf = nd.RBF(basis='gaussian_spline', radius=1.)
         result = rbf.basis(X=1., r=1.)
-        assert result == 0.36787944117144233
+        np.testing.assert_almost_equal(result, 0.36787944117144233)
 
     def test_biharmonic_evaluation(self):
         rbf = nd.RBF(basis='multi_quadratic_biharmonic_spline', radius=1.)
         result = rbf.basis(X=1., r=1.)
-        assert result == 1.4142135623730951
+        np.testing.assert_almost_equal(result, 1.4142135623730951)
 
     def test_inv_biharmonic_evaluation(self):
         rbf = nd.RBF(basis='inv_multi_quadratic_biharmonic_spline', radius=1.)
         result = rbf.basis(X=1., r=1.)
-        assert result == 0.7071067811865475
+        np.testing.assert_almost_equal(result, 0.7071067811865475)
 
     def test_thin_plate_evaluation(self):
         rbf = nd.RBF(basis='thin_plate_spline', radius=1.)
         result = rbf.basis(X=1., r=0.5)
-        assert result == 2.772588722239781
+        np.testing.assert_almost_equal(result, 2.772588722239781)
 
     def test_wendland_evaluation(self):
         rbf = nd.RBF(basis='beckert_wendland_c2_basis', radius=1.)
         result = rbf.basis(X=1., r=2.)
-        assert result == 0.1875
+        np.testing.assert_almost_equal(result, 0.1875)
 
     def test_wendland_outside_cutoff(self):
         rbf = nd.RBF(basis='beckert_wendland_c2_basis', radius=1.)
         result = rbf.basis(X=2., r=1.)
-        assert result == 0.0
+        np.testing.assert_almost_equal(result, 0.0)
 
     def test_weight_matrix(self):
         x = np.arange(10)
         rbf = nd.RBF(basis='beckert_wendland_c2_basis', radius=1.)
         weights_matrix = rbf.weights_matrix(X1=x, X2=x)
         expected = np.diag(np.ones(10))
-        np.testing.assert_array_equal(weights_matrix, expected)
+        np.testing.assert_array_almost_equal(weights_matrix, expected)
 
     def test_reconstruct_f_gaussian(self):
         x, y, xx, yy = sample_data()
@@ -91,7 +91,7 @@ def test_reconstruct_f_gaussian(self):
         # for argmin(yy) nearest to point y=20.25, where y=x*x
         idx = (np.abs(xx - 4.5)).argmin()
         idx2 = (np.abs(yy - 20.25)).argmin()
-        assert idx == idx2
+        np.testing.assert_array_almost_equal(idx, idx2)
 
     def test_reconstruct_f_biharmonic(self):
         x, y, xx, yy = sample_data()
@@ -106,7 +106,7 @@ def test_reconstruct_f_biharmonic(self):
         # for argmin(yy) nearest to point y=20.25, where y=x*x
         idx = (np.abs(xx - 4.5)).argmin()
         idx2 = (np.abs(yy - 20.25)).argmin()
-        assert idx == idx2
+        np.testing.assert_array_almost_equal(idx, idx2)
 
     def test_reconstruct_f_inv_biharmonic(self):
         x, y, xx, yy = sample_data()
@@ -121,7 +121,7 @@ def test_reconstruct_f_inv_biharmonic(self):
         # for argmin(yy) nearest to point y=20.25, where y=x*x
         idx = (np.abs(xx - 4.5)).argmin()
         idx2 = (np.abs(yy - 20.25)).argmin()
-        assert idx == idx2
+        np.testing.assert_array_almost_equal(idx, idx2)
 
     def test_reconstruct_f_plate(self):
         x, y, xx, yy = sample_data()
@@ -136,7 +136,7 @@ def test_reconstruct_f_plate(self):
         # for argmin(yy) nearest to point y=20.25, where y=x*x
         idx = (np.abs(xx - 4.5)).argmin()
         idx2 = (np.abs(yy - 20.25)).argmin()
-        assert idx == idx2
+        np.testing.assert_array_almost_equal(idx, idx2)
 
     def test_reconstruct_f_wendland(self):
         x, y, xx, yy = sample_data()
@@ -151,7 +151,7 @@ def test_reconstruct_f_wendland(self):
         # for argmin(yy) nearest to point y=20.25, where y=x*x
         idx = (np.abs(xx - 4.5)).argmin()
         idx2 = (np.abs(yy - 20.25)).argmin()
-        assert idx == idx2
+        np.testing.assert_array_almost_equal(idx, idx2)
 
     def test_reconstruct_f_scalar(self):
         x = np.arange(10)