Repository: ross-rotordynamics/ross
Pull request #307 (2019-09-25 13:33:00)

Add BallBearing and RollerBearing Element
This PR adds two new classes:
- BallBearingElement
- RollerBearingElement

These classes create a bearing element for ball and roller bearings models. They are instantiated based os some geometric and constructive parameters like: number of rotating elements, size (diameter for ballbearing; length for rollerbearing) of the rotating elements, contact angle, and static loading force. 

The direct stiffness coefficientes are calculated using these parameters.

I also left an opened option to instantiate the direct damping coefficients. But if it's set as None, they are calculated based on the stiffness coefficient as the literature suggests.
However the cross-coupling coefficients are set to zero, and the coefficients are not dependent on speed.

```
    def __init__(
        self,
        n,
        n_balls,
        d_balls,
        fs,
        alpha,
        cxx=None,
        cyy=None,
        tag=None,
    ):

        Kb = 13.0e6
        kyy = (
            Kb * n_balls ** (2./3) * d_balls ** (1./3) *
            fs ** (1./3) * (np.cos(alpha)) ** (5./3)
        )

        nb = [8, 12, 16]
        ratio = [0.46, 0.64, 0.73]
        dict_ratio = dict(zip(nb, ratio))

        if n_balls in dict_ratio.keys():
            kxx = dict_ratio[n_balls] * kyy
        else:
            f = interpolate.interp1d(nb, ratio, kind="quadratic")
            kxx = f(n_balls)

        if cxx is None:
            cxx = 1.25e-5 * kxx
        if cyy is None:
            cyy = 1.25e-5 * kyy

        super().__init__(
            n=n,
            w=None,
            kxx=kxx,
            kxy=0.0,
            kyx=0.0,
            kyy=kyy,
            cxx=cxx,
            cxy=0.0,
            cyx=0.0,
            cyy=cyy,
            tag=tag,
        )
```


Patch:
diff --git a/.travis.yml b/.travis.yml
index 98f2b2ff..91699779 100644
--- a/.travis.yml
+++ b/.travis.yml
@@ -1,21 +1,17 @@
 sudo: required
 dist: xenial
 language: python
+
 matrix:
   include:
-  - python: 3.6
-  - python: 3.7
+    - python: 3.6
+    - python: 3.7
+
 install:
-- pip install -e .[dev]
+  - pip install -e .[dev]
+
 script:
-- pytest ross --cov=ross
+  - pytest ross --cov=ross
+
 after_success:
-- codecov
-deploy:
-  provider: pypi
-  user: raphaelts
-  password:
-    secure: GgnV7C40a+CXKa+y35qSb9SMcceVM42LZrEfqh983m8fEpP+B6LyM5bWlvAsP2rLIkjiZl0zU5ZzgzAmV8buSUqwOmVnagjeaaicrOpWCeuc+1A8Uz0amoII8cZQJCvyYJScxiYhxap8eOqwXFZ2SKnVPzxaToF1T5NaJxoESIriuXD6mW0UOUaXaEX8qKB37iFVIyFy/YYdmdBG6atRnL2VEBzC/Gi4GTFim1tSRkEE9MX6TFzF7v9UqXcrdFcCdQt5KtKWcFTMlV8QTPiT+ZZqPXbjliwF78z/rOpDe1Cv0quMbGLEBvZjxUHJFNfyA4+s6hqGgC+3VXKNZ0f0nIa6mjt1rvd/UJM3105iaGG51xjCHyw+CpEPASjB3i15wMwevKV2di/ecGA0hZSXVdycENymf9TArjEEU/bXJnByQY/ejFJScueBmYnE3U8kn87eT2Dtz3GbipHiHqrhNLyCbe+PdpsZH8y7cF9BWeUTxRQlp8+Y/67jM2bCeGtuxiLksucipkssg5a6ek59SPRGb50hDjQEwk/EdDlON/4Exu6ABBVNBskjXiSX3HzIcHa+jz3Xr84SblPbFzn4ymd706/kWhk73wiRaYGkPgV6l7fdckCO23K0bAiCZWSojf2WgbA1rKbs7/yzK3HJbEu+TUVq9QBwj4RLnFw2P7M=
-  on:
-    tags: true
-    python: 3.7
+  - codecov
diff --git a/CONTRIBUTING.rst b/CONTRIBUTING.rst
index dcd60e04..2c439eff 100644
--- a/CONTRIBUTING.rst
+++ b/CONTRIBUTING.rst
@@ -132,15 +132,18 @@ The following blog posts have some good information on how to write commit messa
 
 Step 5: Push changes to the main repo
 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+For contributors
+++++++++++++++++
 To create a Pull Request (PR), refer to `the github PR guide <https://help.github.com/articles/about-pull-requests/>`_.
 
-Making new releases
--------------------
-To make a new release we need only to create a tag using git and push to GitHub:
+For core developers
++++++++++++++++++++
+If there are only a few, unrelated commits:
 
-    $ git tag <version number>
+::
 
-    $ git push upstream --tags
+    git fetch upstream
+    git rebase upstream/master
+    git push upstream my-feature-branch:master
 
-Pushing the new tag to the GitHub repository will start a new build on Travis CI. If all the tests succeed, Travis will
-upload the new package to PyPI (see the deploy command on .travis.yml).
diff --git a/ross/__init__.py b/ross/__init__.py
index e29a25ea..b9c92489 100644
--- a/ross/__init__.py
+++ b/ross/__init__.py
@@ -1,4 +1,4 @@
-__version__ = "0.1.6"
+__version__ = "0.1.2"
 
 from .bearing_seal_element import *
 from .disk_element import *
diff --git a/ross/bearing_seal_element.py b/ross/bearing_seal_element.py
index 180a35e1..1ebde0a1 100644
--- a/ross/bearing_seal_element.py
+++ b/ross/bearing_seal_element.py
@@ -1105,6 +1105,635 @@ class SealElement(BearingElement):
         )
 
 
+class BallBearingElement(BearingElement):
+    """A bearing element for ball bearings.
+    This class will create a bearing element based on some geometric and
+    constructive parameters of ball bearings. The main difference is that
+    cross-coupling stiffness and damping are not modeled in this case.    
+
+    Parameters
+    ----------
+    n: int
+        Node which the bearing will be located in
+    n_balls: float
+        Number of steel spheres in the bearing
+    d_balls: float
+        Diameter of the steel sphere
+    fs: float,optional
+        Static bearing loading force
+    alpha: float, optional
+        Contact angle between the steel sphere and the inner / outer raceway
+        (defaults to cxx)
+    cxx: float, optional
+        Direct stiffness in the x direction.
+        Default is None)
+    cyy: float, optional
+        Direct damping in the y direction.
+        Defaults is None
+    tag : str, optional
+        A tag to name the element
+        Default is None
+    Examples
+    --------
+    >>> n = 0
+    >>> n_balls= 8
+    >>> d_balls = 0.03
+    >>> fs = 500.0
+    >>> alpha = np.pi / 6
+    >>> tag = "ballbearing"
+    >>> bearing = BallBearingElement(n=n, n_balls=n_balls, d_balls=d_balls,
+    ...                              fs=fs, alpha=alpha, tag=tag)
+    >>> bearing.K(0)
+    array([[4.64168838e+07, 0.00000000e+00],
+           [0.00000000e+00, 1.00906269e+08]])
+    """
+    def __init__(
+        self,
+        n,
+        n_balls,
+        d_balls,
+        fs,
+        alpha,
+        cxx=None,
+        cyy=None,
+        tag=None,
+    ):
+
+        Kb = 13.0e6
+        kyy = (
+            Kb * n_balls ** (2./3) * d_balls ** (1./3) *
+            fs ** (1./3) * (np.cos(alpha)) ** (5./3)
+        )
+
+        nb = [8, 12, 16]
+        ratio = [0.46, 0.64, 0.73]
+        dict_ratio = dict(zip(nb, ratio))
+
+        if n_balls in dict_ratio.keys():
+            kxx = dict_ratio[n_balls] * kyy
+        else:
+            f = interpolate.interp1d(nb, ratio, kind="quadratic")
+            kxx = f(n_balls)
+
+        if cxx is None:
+            cxx = 1.25e-5 * kxx
+        if cyy is None:
+            cyy = 1.25e-5 * kyy
+
+        super().__init__(
+            n=n,
+            w=None,
+            kxx=kxx,
+            kxy=0.0,
+            kyx=0.0,
+            kyy=kyy,
+            cxx=cxx,
+            cxy=0.0,
+            cyx=0.0,
+            cyy=cyy,
+            tag=tag,
+        )
+
+        self.color = "#77ACA2"
+
+    def patch(self, position, mean, ax, **kwargs):
+        """Bearing element patch.
+        Patch that will be used to draw the bearing element.
+
+        Parameters
+        ----------
+        position : tuple
+            Position (z, y) in which the patch will be drawn.
+        mean : float
+            Mean value of shaft element outer diameters
+        ax : matplotlib axes, optional
+            Axes in which the plot will be drawn.
+
+        Returns
+        -------
+        """
+        default_values = dict(lw=1.0, alpha=1.0, c="k")
+        for k, v in default_values.items():
+            kwargs.setdefault(k, v)
+
+        # geometric factors
+        zpos, ypos = position
+        coils = 6  # number of points to generate spring
+        step = mean / 5  # spring step
+        a = 2.2  # range(1.8 to 2.8)
+        b = 1.4  # range(1.2 to 1.5)
+        c = mean / 2.9  # defines width
+        d = (coils / a) * step  # modifies damper width
+        n = 5  # number of ground lines
+
+        zs0 = zpos - c
+        zs1 = zpos + c
+        ys0 = ypos + mean / 2
+
+        # plot bottom base
+        x_bot = [zpos, zpos, zs0, zs1]
+        yl_bot = [ypos, ys0, ys0, ys0]
+        yu_bot = [-y for y in yl_bot]
+        ax.add_line(mlines.Line2D(x_bot, yl_bot, **kwargs))
+        ax.add_line(mlines.Line2D(x_bot, yu_bot, **kwargs))
+
+        # plot top base
+        x_top = [zpos, zpos, zs0, zs1]
+        yl_top = [
+            b * ys0 + (coils + 1) * step,
+            ys0 + (coils + 1) * step,
+            ys0 + (coils + 1) * step,
+            ys0 + (coils + 1) * step,
+        ]
+        yu_top = [-y for y in yl_top]
+
+        ax.add_line(mlines.Line2D(x_top, yl_top, **kwargs))
+        ax.add_line(mlines.Line2D(x_top, yu_top, **kwargs))
+
+        # plot ground
+        zl_g = [zs0 - step, zs1 + step]
+        yl_g = [yl_top[0], yl_top[0]]
+        yu_g = [-y for y in yl_g]
+        ax.add_line(mlines.Line2D(zl_g, yl_g, **kwargs))
+        ax.add_line(mlines.Line2D(zl_g, yu_g, **kwargs))
+
+        step2 = (zl_g[1] - zl_g[0]) / n
+        for i in range(n + 1):
+            zl_g2 = [(zs0 - step) + step2 * (i), (zs0 - step) + step2 * (i + 1)]
+            yl_g2 = [yl_g[0], 1.1 * yl_g[0]]
+            yu_g2 = [-y for y in yl_g2]
+            ax.add_line(mlines.Line2D(zl_g2, yl_g2, **kwargs))
+            ax.add_line(mlines.Line2D(zl_g2, yu_g2, **kwargs))
+
+        # plot spring
+        z_spring = np.array([zs0, zs0, zs0, zs0])
+        yl_spring = np.array(
+            [ys0, ys0 + step, ys0 + coils * step, ys0 + (coils + 1) * step]
+        )
+
+        for i in range(coils):
+            z_spring = np.insert(z_spring, i + 2, zs0 - (-1) ** i * step)
+            yl_spring = np.insert(yl_spring, i + 2, ys0 + (i + 1) * step)
+        yu_spring = [-y for y in yl_spring]
+
+        ax.add_line(mlines.Line2D(z_spring, yl_spring, **kwargs))
+        ax.add_line(mlines.Line2D(z_spring, yu_spring, **kwargs))
+
+        # plot damper
+        z_damper1 = [zs1, zs1]
+        yl_damper1 = [ys0, ys0 + d]
+        yu_damper1 = [-y for y in yl_damper1]
+
+        ax.add_line(mlines.Line2D(z_damper1, yl_damper1, **kwargs))
+        ax.add_line(mlines.Line2D(z_damper1, yu_damper1, **kwargs))
+
+        z_damper2 = [zs1 - step, zs1 - step, zs1 + step, zs1 + step]
+        yl_damper2 = [ys0 + 2 * d, ys0 + d, ys0 + d, ys0 + 2 * d]
+        yu_damper2 = [-y for y in yl_damper2]
+
+        ax.add_line(mlines.Line2D(z_damper2, yl_damper2, **kwargs))
+        ax.add_line(mlines.Line2D(z_damper2, yu_damper2, **kwargs))
+
+        z_damper3 = [zs1 - step, zs1 + step, zs1, zs1]
+        yl_damper3 = [
+            ys0 + 1.5 * d,
+            ys0 + 1.5 * d,
+            ys0 + 1.5 * d,
+            ys0 + (coils + 1) * step,
+        ]
+        yu_damper3 = [-y for y in yl_damper3]
+
+        ax.add_line(mlines.Line2D(z_damper3, yl_damper3, **kwargs))
+        ax.add_line(mlines.Line2D(z_damper3, yu_damper3, **kwargs))
+
+    def bokeh_patch(self, position, mean, bk_ax, **kwargs):
+        """Bearing element patch.
+        Patch that will be used to draw the bearing element.
+
+        Parameters
+        ----------
+        position : tuple
+            Position (z, y) in which the patch will be drawn.
+        mean : float
+            Mean value of shaft element outer diameters
+        bk_ax : bokeh plotting axes, optional
+            Axes in which the plot will be drawn.
+
+        Returns
+        -------
+        """
+        default_values = dict(line_width=3, line_alpha=1, color=bokeh_colors[1])
+        for k, v in default_values.items():
+            kwargs.setdefault(k, v)
+
+        # geometric factors
+        zpos, ypos = position
+        coils = 6  # number of points to generate spring
+        step = mean / 5  # spring step
+        a = 2.2  # range(1.8 to 2.8)
+        b = 1.4  # range(1.2 to 1.5)
+        c = mean / 2.9  # defines width
+        d = (coils / a) * step  # modifies damper width
+        n = 5  # number of ground lines
+
+        zs0 = zpos - c
+        zs1 = zpos + c
+        ys0 = ypos + mean / 2
+
+        # plot bottom base
+        x_bot = [zpos, zpos, zs0, zs1]
+        yl_bot = [ypos, ys0, ys0, ys0]
+        yu_bot = [-y for y in yl_bot]
+        bk_ax.line(x=x_bot, y=yl_bot, **kwargs)
+        bk_ax.line(x=x_bot, y=yu_bot, **kwargs)
+
+        # plot top base
+        x_top = [zpos, zpos, zs0, zs1]
+        yl_top = [
+            b * ys0 + (coils + 1) * step,
+            ys0 + (coils + 1) * step,
+            ys0 + (coils + 1) * step,
+            ys0 + (coils + 1) * step,
+        ]
+        yu_top = [-y for y in yl_top]
+
+        bk_ax.line(x=x_top, y=yl_top, legend="Bearing", **kwargs)
+        bk_ax.line(x=x_top, y=yu_top, legend="Bearing", **kwargs)
+
+        # plot ground
+        if self.n_link is None:
+            zl_g = [zs0 - step, zs1 + step]
+            yl_g = [yl_top[0], yl_top[0]]
+            yu_g = [-y for y in yl_g]
+            bk_ax.line(x=zl_g, y=yl_g, **kwargs)
+            bk_ax.line(x=zl_g, y=yu_g, **kwargs)
+
+            step2 = (zl_g[1] - zl_g[0]) / n
+            for i in range(n + 1):
+                zl_g2 = [(zs0 - step) + step2 * (i), (zs0 - step) + step2 * (i + 1)]
+                yl_g2 = [yl_g[0], 1.1 * yl_g[0]]
+                yu_g2 = [-y for y in yl_g2]
+                bk_ax.line(x=zl_g2, y=yl_g2, **kwargs)
+                bk_ax.line(x=zl_g2, y=yu_g2, **kwargs)
+
+        # plot spring
+        z_spring = np.array([zs0, zs0, zs0, zs0])
+        yl_spring = np.array(
+            [ys0, ys0 + step, ys0 + coils * step, ys0 + (coils + 1) * step]
+        )
+
+        for i in range(coils):
+            z_spring = np.insert(z_spring, i + 2, zs0 - (-1) ** i * step)
+            yl_spring = np.insert(yl_spring, i + 2, ys0 + (i + 1) * step)
+        yu_spring = [-y for y in yl_spring]
+
+        bk_ax.line(x=z_spring, y=yl_spring, **kwargs)
+        bk_ax.line(x=z_spring, y=yu_spring, **kwargs)
+
+        # plot damper
+        z_damper1 = [zs1, zs1]
+        yl_damper1 = [ys0, ys0 + d]
+        yu_damper1 = [-y for y in yl_damper1]
+
+        bk_ax.line(x=z_damper1, y=yl_damper1, **kwargs)
+        bk_ax.line(x=z_damper1, y=yu_damper1, **kwargs)
+
+        z_damper2 = [zs1 - step, zs1 - step, zs1 + step, zs1 + step]
+        yl_damper2 = [ys0 + 2 * d, ys0 + d, ys0 + d, ys0 + 2 * d]
+        yu_damper2 = [-y for y in yl_damper2]
+
+        bk_ax.line(x=z_damper2, y=yl_damper2, **kwargs)
+        bk_ax.line(x=z_damper2, y=yu_damper2, **kwargs)
+
+        z_damper3 = [zs1 - step, zs1 + step, zs1, zs1]
+        yl_damper3 = [
+            ys0 + 1.5 * d,
+            ys0 + 1.5 * d,
+            ys0 + 1.5 * d,
+            ys0 + (coils + 1) * step,
+        ]
+        yu_damper3 = [-y for y in yl_damper3]
+
+        bk_ax.line(x=z_damper3, y=yl_damper3, **kwargs)
+        bk_ax.line(x=z_damper3, y=yu_damper3, **kwargs)
+
+
+class RollerBearingElement(BearingElement):
+    """A bearing element for roller bearings.
+    This class will create a bearing element based on some geometric and
+    constructive parameters of roller bearings. The main difference is that
+    cross-coupling stiffness and damping are not modeled in this case.    
+
+    Parameters
+    ----------
+    n: int
+        Node which the bearing will be located in
+    n_rollers: float
+        Number of steel spheres in the bearing
+    l_rollers: float
+        Diameter of the steel sphere
+    fs: float,optional
+        Static bearing loading force
+    alpha: float, optional
+        Contact angle between the steel sphere and the inner / outer raceway
+        (defaults to cxx)
+    cxx: float, optional
+        Direct stiffness in the x direction.
+        Default is None)
+    cyy: float, optional
+        Direct damping in the y direction.
+        Defaults is None
+    tag : str, optional
+        A tag to name the element
+        Default is None
+    Examples
+    --------
+    >>> n = 0
+    >>> n_rollers= 8
+    >>> l_rollers = 0.03
+    >>> fs = 500.0
+    >>> alpha = np.pi / 6
+    >>> tag = "rollerbearing"
+    >>> bearing = RollerBearingElement(n=n, n_rollers=n_rollers, l_rollers=l_rollers,
+    ...                            fs=fs, alpha=alpha, tag=tag)
+    >>> bearing.K(0)
+    array([[2.72821927e+08, 0.00000000e+00],
+           [0.00000000e+00, 5.56779444e+08]])
+    """
+    def __init__(
+        self,
+        n,
+        n_rollers,
+        l_rollers,
+        fs,
+        alpha,
+        cxx=None,
+        cyy=None,
+        tag=None,
+    ):
+
+        Kb = 1.0e9
+        kyy = (
+            Kb * n_rollers ** 0.9 * l_rollers ** 0.8 *
+            fs ** 0.1 * (np.cos(alpha)) ** 1.9
+        )
+
+        nr = [8, 12, 16]
+        ratio = [0.49, 0.66, 0.74]
+        dict_ratio = dict(zip(nr, ratio))
+
+        if n_rollers in dict_ratio.keys():
+            kxx = dict_ratio[n_rollers] * kyy
+        else:
+            f = interpolate.interp1d(nr, ratio, kind="quadratic")
+            kxx = f(n_rollers)
+
+        if cxx is None:
+            cxx = 1.25e-5 * kxx
+        if cyy is None:
+            cyy = 1.25e-5 * kyy
+
+        super().__init__(
+            n=n,
+            w=None,
+            kxx=kxx,
+            kxy=0.0,
+            kyx=0.0,
+            kyy=kyy,
+            cxx=cxx,
+            cxy=0.0,
+            cyx=0.0,
+            cyy=cyy,
+            tag=tag,
+        )
+
+        self.color = "#77ACA2"
+
+    def patch(self, position, mean, ax, **kwargs):
+        """Bearing element patch.
+        Patch that will be used to draw the bearing element.
+
+        Parameters
+        ----------
+        position : tuple
+            Position (z, y) in which the patch will be drawn.
+        mean : float
+            Mean value of shaft element outer diameters
+        ax : matplotlib axes, optional
+            Axes in which the plot will be drawn.
+
+        Returns
+        -------
+
+        Examples
+        --------
+        """
+        default_values = dict(lw=1.0, alpha=1.0, c="k")
+        for k, v in default_values.items():
+            kwargs.setdefault(k, v)
+
+        # geometric factors
+        zpos, ypos = position
+        coils = 6  # number of points to generate spring
+        step = mean / 5  # spring step
+        a = 2.2  # range(1.8 to 2.8)
+        b = 1.4  # range(1.2 to 1.5)
+        c = mean / 2.9  # defines width
+        d = (coils / a) * step  # modifies damper width
+        n = 5  # number of ground lines
+
+        zs0 = zpos - c
+        zs1 = zpos + c
+        ys0 = ypos + mean / 2
+
+        # plot bottom base
+        x_bot = [zpos, zpos, zs0, zs1]
+        yl_bot = [ypos, ys0, ys0, ys0]
+        yu_bot = [-y for y in yl_bot]
+        ax.add_line(mlines.Line2D(x_bot, yl_bot, **kwargs))
+        ax.add_line(mlines.Line2D(x_bot, yu_bot, **kwargs))
+
+        # plot top base
+        x_top = [zpos, zpos, zs0, zs1]
+        yl_top = [
+            b * ys0 + (coils + 1) * step,
+            ys0 + (coils + 1) * step,
+            ys0 + (coils + 1) * step,
+            ys0 + (coils + 1) * step,
+        ]
+        yu_top = [-y for y in yl_top]
+
+        ax.add_line(mlines.Line2D(x_top, yl_top, **kwargs))
+        ax.add_line(mlines.Line2D(x_top, yu_top, **kwargs))
+
+        # plot ground
+        zl_g = [zs0 - step, zs1 + step]
+        yl_g = [yl_top[0], yl_top[0]]
+        yu_g = [-y for y in yl_g]
+        ax.add_line(mlines.Line2D(zl_g, yl_g, **kwargs))
+        ax.add_line(mlines.Line2D(zl_g, yu_g, **kwargs))
+
+        step2 = (zl_g[1] - zl_g[0]) / n
+        for i in range(n + 1):
+            zl_g2 = [(zs0 - step) + step2 * (i), (zs0 - step) + step2 * (i + 1)]
+            yl_g2 = [yl_g[0], 1.1 * yl_g[0]]
+            yu_g2 = [-y for y in yl_g2]
+            ax.add_line(mlines.Line2D(zl_g2, yl_g2, **kwargs))
+            ax.add_line(mlines.Line2D(zl_g2, yu_g2, **kwargs))
+
+        # plot spring
+        z_spring = np.array([zs0, zs0, zs0, zs0])
+        yl_spring = np.array(
+            [ys0, ys0 + step, ys0 + coils * step, ys0 + (coils + 1) * step]
+        )
+
+        for i in range(coils):
+            z_spring = np.insert(z_spring, i + 2, zs0 - (-1) ** i * step)
+            yl_spring = np.insert(yl_spring, i + 2, ys0 + (i + 1) * step)
+        yu_spring = [-y for y in yl_spring]
+
+        ax.add_line(mlines.Line2D(z_spring, yl_spring, **kwargs))
+        ax.add_line(mlines.Line2D(z_spring, yu_spring, **kwargs))
+
+        # plot damper
+        z_damper1 = [zs1, zs1]
+        yl_damper1 = [ys0, ys0 + d]
+        yu_damper1 = [-y for y in yl_damper1]
+
+        ax.add_line(mlines.Line2D(z_damper1, yl_damper1, **kwargs))
+        ax.add_line(mlines.Line2D(z_damper1, yu_damper1, **kwargs))
+
+        z_damper2 = [zs1 - step, zs1 - step, zs1 + step, zs1 + step]
+        yl_damper2 = [ys0 + 2 * d, ys0 + d, ys0 + d, ys0 + 2 * d]
+        yu_damper2 = [-y for y in yl_damper2]
+
+        ax.add_line(mlines.Line2D(z_damper2, yl_damper2, **kwargs))
+        ax.add_line(mlines.Line2D(z_damper2, yu_damper2, **kwargs))
+
+        z_damper3 = [zs1 - step, zs1 + step, zs1, zs1]
+        yl_damper3 = [
+            ys0 + 1.5 * d,
+            ys0 + 1.5 * d,
+            ys0 + 1.5 * d,
+            ys0 + (coils + 1) * step,
+        ]
+        yu_damper3 = [-y for y in yl_damper3]
+
+        ax.add_line(mlines.Line2D(z_damper3, yl_damper3, **kwargs))
+        ax.add_line(mlines.Line2D(z_damper3, yu_damper3, **kwargs))
+
+    def bokeh_patch(self, position, mean, bk_ax, **kwargs):
+        """Bearing element patch.
+        Patch that will be used to draw the bearing element.
+
+        Parameters
+        ----------
+        position : tuple
+            Position (z, y) in which the patch will be drawn.
+        mean : float
+            Mean value of shaft element outer diameters
+        bk_ax : bokeh plotting axes, optional
+            Axes in which the plot will be drawn.
+
+        Returns
+        -------
+        """
+        default_values = dict(line_width=3, line_alpha=1, color=bokeh_colors[1])
+        for k, v in default_values.items():
+            kwargs.setdefault(k, v)
+
+        # geometric factors
+        zpos, ypos = position
+        coils = 6  # number of points to generate spring
+        step = mean / 5  # spring step
+        a = 2.2  # range(1.8 to 2.8)
+        b = 1.4  # range(1.2 to 1.5)
+        c = mean / 2.9  # defines width
+        d = (coils / a) * step  # modifies damper width
+        n = 5  # number of ground lines
+
+        zs0 = zpos - c
+        zs1 = zpos + c
+        ys0 = ypos + mean / 2
+
+        # plot bottom base
+        x_bot = [zpos, zpos, zs0, zs1]
+        yl_bot = [ypos, ys0, ys0, ys0]
+        yu_bot = [-y for y in yl_bot]
+        bk_ax.line(x=x_bot, y=yl_bot, **kwargs)
+        bk_ax.line(x=x_bot, y=yu_bot, **kwargs)
+
+        # plot top base
+        x_top = [zpos, zpos, zs0, zs1]
+        yl_top = [
+            b * ys0 + (coils + 1) * step,
+            ys0 + (coils + 1) * step,
+            ys0 + (coils + 1) * step,
+            ys0 + (coils + 1) * step,
+        ]
+        yu_top = [-y for y in yl_top]
+
+        bk_ax.line(x=x_top, y=yl_top, legend="Bearing", **kwargs)
+        bk_ax.line(x=x_top, y=yu_top, legend="Bearing", **kwargs)
+
+        # plot ground
+        if self.n_link is None:
+            zl_g = [zs0 - step, zs1 + step]
+            yl_g = [yl_top[0], yl_top[0]]
+            yu_g = [-y for y in yl_g]
+            bk_ax.line(x=zl_g, y=yl_g, **kwargs)
+            bk_ax.line(x=zl_g, y=yu_g, **kwargs)
+
+            step2 = (zl_g[1] - zl_g[0]) / n
+            for i in range(n + 1):
+                zl_g2 = [(zs0 - step) + step2 * (i), (zs0 - step) + step2 * (i + 1)]
+                yl_g2 = [yl_g[0], 1.1 * yl_g[0]]
+                yu_g2 = [-y for y in yl_g2]
+                bk_ax.line(x=zl_g2, y=yl_g2, **kwargs)
+                bk_ax.line(x=zl_g2, y=yu_g2, **kwargs)
+
+        # plot spring
+        z_spring = np.array([zs0, zs0, zs0, zs0])
+        yl_spring = np.array(
+            [ys0, ys0 + step, ys0 + coils * step, ys0 + (coils + 1) * step]
+        )
+
+        for i in range(coils):
+            z_spring = np.insert(z_spring, i + 2, zs0 - (-1) ** i * step)
+            yl_spring = np.insert(yl_spring, i + 2, ys0 + (i + 1) * step)
+        yu_spring = [-y for y in yl_spring]
+
+        bk_ax.line(x=z_spring, y=yl_spring, **kwargs)
+        bk_ax.line(x=z_spring, y=yu_spring, **kwargs)
+
+        # plot damper
+        z_damper1 = [zs1, zs1]
+        yl_damper1 = [ys0, ys0 + d]
+        yu_damper1 = [-y for y in yl_damper1]
+
+        bk_ax.line(x=z_damper1, y=yl_damper1, **kwargs)
+        bk_ax.line(x=z_damper1, y=yu_damper1, **kwargs)
+
+        z_damper2 = [zs1 - step, zs1 - step, zs1 + step, zs1 + step]
+        yl_damper2 = [ys0 + 2 * d, ys0 + d, ys0 + d, ys0 + 2 * d]
+        yu_damper2 = [-y for y in yl_damper2]
+
+        bk_ax.line(x=z_damper2, y=yl_damper2, **kwargs)
+        bk_ax.line(x=z_damper2, y=yu_damper2, **kwargs)
+
+        z_damper3 = [zs1 - step, zs1 + step, zs1, zs1]
+        yl_damper3 = [
+            ys0 + 1.5 * d,
+            ys0 + 1.5 * d,
+            ys0 + 1.5 * d,
+            ys0 + (coils + 1) * step,
+        ]
+        yu_damper3 = [-y for y in yl_damper3]
+
+        bk_ax.line(x=z_damper3, y=yl_damper3, **kwargs)
+        bk_ax.line(x=z_damper3, y=yu_damper3, **kwargs)
+
+
 def bearing_example():
     """This function returns an instance of a simple bearing.
     The purpose is to make available a simple model
diff --git a/ross/element.py b/ross/element.py
index 40d34c7a..4a4a82c6 100644
--- a/ross/element.py
+++ b/ross/element.py
@@ -7,8 +7,8 @@ import toml
 
 class Element(ABC):
     """Element class.
-    This class is a general abstract class to be implemented in other files, in order to
-    create specific elements for the user.
+    This class is a general class to be called for other files which
+    create specific elements for the user
     """
 
     def __init__(self, n):
@@ -17,161 +17,36 @@ class Element(ABC):
         pass
 
     def save(self, file_name):
-        """Saves the element in a file.
-        Parameters
-        ----------
-        file_name: string
-            The name of the file the element will be saved in.
-        Returns
-        -------
-        None
-        Examples
-        --------
-        >>> # Example using DiskElement
-        >>> from ross.disk_element import disk_example
-        >>> disk = disk_example()
-        >>> disk.save('DiskElement.toml')
-        """
         pass
 
     @staticmethod
-    def load(file_name):
-        """Loads elements saved in a file.
-        Parameters
-        ----------
-        file_name: string
-            The name of the file to be loaded.
-        Returns
-        -------
-        The element object.
-        Examples
-        --------
-        >>> # Example using BearingElement
-        >>> from ross.bearing_seal_element import bearing_example
-        >>> from ross.bearing_seal_element import BearingElement
-        >>> bearing1 = bearing_example()
-        >>> bearing1.save('BearingElement.toml')
-        >>> list_of_bearings = BearingElement.load('BearingElement.toml')
-        >>> bearing1 == list_of_bearings[0]
-        True
-        """
+    def load(self, file_name):
         pass
 
     @abstractmethod
     def M(self):
-        """Mass matrix.
-        Parameters
-        ----------
-
-        Returns
-        -------
-        A matrix of floats.
-
-        Examples
-        --------
-        >>> # Example using BearingElement
-        >>> from ross.bearing_seal_element import bearing_example
-        >>> bearing = bearing_example()
-        >>> bearing.M()
-        array([[0., 0.],
-               [0., 0.]])
-        """
+        """Mass matrix."""
         pass
 
     @abstractmethod
     def C(self, frequency):
-        """Frequency dependent damping coefficients matrix.
-        Parameters
-        ----------
-        frequency: float
-            The frequency in which the coefficients depend on.
-
-        Returns
-        -------
-        A matrix of floats.
-
-        Examples
-        --------
-        >>> # Example using BearingElement
-        >>> from ross.bearing_seal_element import bearing_example
-        >>> bearing = bearing_example()
-        >>> bearing.C(0)
-        array([[200.,   0.],
-               [  0., 150.]])
-        """
+        """Frequency dependent damping coefficients matrix."""
         pass
 
     @abstractmethod
     def K(self, frequency):
-        """Frequency dependent stiffness coefficients matrix.
-        Parameters
-        ----------
-        frequency: float
-            The frequency in which the coefficients depend on.
-
-        Returns
-        -------
-        A matrix of floats.
-
-        Examples
-        --------
-        >>> # Example using BearingElement
-        >>> from ross.bearing_seal_element import bearing_example
-        >>> bearing = bearing_example()
-        >>> bearing.K(0)
-        array([[1000000.,       0.],
-               [      0.,  800000.]])
-        """
+        """Frequency dependent stiffness coefficients matrix."""
         pass
 
     @abstractmethod
     def G(self):
-        """Gyroscopic matrix.
-        Parameters
-        ----------
-
-        Returns
-        -------
-        A matrix of floats.
-
-        Examples
-        --------
-        >>> # Example using BearingElement
-        >>> from ross.bearing_seal_element import bearing_example
-        >>> bearing = bearing_example()
-        >>> bearing.G()
-        array([[0., 0.],
-               [0., 0.]])
-        """
+        """Gyroscopic matrix."""
         pass
 
     def summary(self):
         """A summary for the element.
 
         A pandas series with the element properties as variables.
-        Parameters
-        ----------
-
-        Returns
-        -------
-        A pandas series.
-
-        Examples
-        --------
-        >>> # Example using DiskElement
-        >>> from ross.disk_element import disk_example
-        >>> disk = disk_example()
-        >>> disk.summary()
-        n                  0
-        n_l                0
-        n_r                0
-        m            32.5897
-        Id          0.178089
-        Ip          0.329564
-        tag             None
-        color        #bc625b
-        type     DiskElement
-        dtype: object
         """
         attributes = self.__dict__
         attributes["type"] = self.__class__.__name__
@@ -179,26 +54,6 @@ class Element(ABC):
 
     @staticmethod
     def load_data(file_name):
-        """Loads elements data saved in a toml format.
-        Parameters
-        ----------
-        file_name: string
-            The name of the file to be loaded.
-
-        Returns
-        -------
-        dict
-            The element parameters presented as a dictionary.
-        Examples
-        --------
-        >>> # Example using BearingElement
-        >>> from ross.bearing_seal_element import bearing_example
-        >>> from ross.bearing_seal_element import BearingElement
-        >>> bearing = bearing_example()
-        >>> bearing.save('BearingElement.toml')
-        >>> BearingElement.load_data('BearingElement.toml') # doctest: +ELLIPSIS
-        {'BearingElement': {'0': {'n': 0, 'kxx': [1000000.0, 1000000.0,...
-        """
         try:
             with open(file_name, "r") as f:
                 data = toml.load(f)
@@ -212,26 +67,6 @@ class Element(ABC):
 
     @staticmethod
     def dump_data(data, file_name):
-        """Dumps element data in a toml file.
-        Parameters
-        ----------
-        data: dict
-            The data that should be dumped.
-        file_name: string
-            The name of the file the data will be dumped in.
-        Returns
-        -------
-
-        Examples
-        --------
-        >>> # Example using BearingElement
-        >>> from ross.bearing_seal_element import bearing_example
-        >>> from ross.bearing_seal_element import BearingElement
-        >>> bearing = bearing_example()
-        >>> bearing.save('BearingElement.toml')
-        >>> data = BearingElement.load_data('BearingElement.toml')
-        >>> BearingElement.dump_data(data, 'bearing_data.toml')
-        """
         with open(file_name, "w") as f:
             toml.dump(data, f)
 
@@ -239,40 +74,24 @@ class Element(ABC):
     def dof_mapping(self):
         """Should return a dictionary with a mapping between degree of freedom
         and its index.
-        Parameters
-        ----------
-
-        Returns
-        -------
-        A dictionary of degrees of freedom.
 
-        Examples
-        --------
-        >>> # Example using BearingElement
-        >>> from ross.bearing_seal_element import bearing_example
-        >>> bearing = bearing_example()
-        >>> bearing.dof_mapping()
-        {'x_0': 0, 'y_0': 1}
+        Example considering a shaft element:
+        def dof_mapping(self):
+            return dict(
+                x_0=0,
+                y_0=1,
+                alpha_0=2,
+                beta_0=3,
+                x_1=4,
+                y_1=5,
+                alpha_1=6,
+                beta_1=7,
+            )
         """
         pass
 
     def dof_local_index(self):
-        """Get the local index for a element specific degree of freedom.
-        Parameters
-        ----------
-
-        Returns
-        -------
-        A named tuple containing the local index.
-
-        Examples
-        --------
-        >>> # Example using BearingElement
-        >>> from ross.bearing_seal_element import bearing_example
-        >>> bearing = bearing_example()
-        >>> bearing.dof_local_index()
-        LocalIndex(x_0=0, y_0=1)
-        """
+        """Get the local index for a element specific degree of freedom."""
         dof_mapping = self.dof_mapping()
         dof_tuple = namedtuple("LocalIndex", dof_mapping)
         local_index = dof_tuple(**dof_mapping)
@@ -280,22 +99,7 @@ class Element(ABC):
         return local_index
 
     def dof_global_index(self):
-        """Get the global index for a element specific degree of freedom.
-        Parameters
-        ----------
-
-        Returns
-        -------
-        A named tuple containing the global index.
-
-        Examples
-        --------
-        >>> # Example using BearingElement
-        >>> from ross.bearing_seal_element import bearing_example
-        >>> bearing = bearing_example()
-        >>> bearing.dof_global_index()
-        GlobalIndex(x_0=0, y_0=1)
-        """
+        """Get the global index for a element specific degree of freedom."""
         dof_mapping = self.dof_mapping()
         global_dof_mapping = {}
         for k, v in dof_mapping.items():
diff --git a/ross/results.py b/ross/results.py
index 4907efe7..02ae3e64 100644
--- a/ross/results.py
+++ b/ross/results.py
@@ -68,34 +68,6 @@ class Results(np.ndarray):
 
 
 class CampbellResults:
-    """Class used to store results and provide plots for Campbell Diagram.
-
-    It's possible to visualize multiples harmonics in a single plot to check
-    other speeds which also excite a specific natural frequency.
-    Two options for plooting are available: Matplotlib and Bokeh. The user
-    chooses between them using the attribute plot_type. The default is bokeh
-
-    Parameters
-    ----------
-    speed_range : array
-        Array with the speed range in rad/s.
-    wd : array
-        Array with the damped natural frequencies
-    log_dec : array
-        Array with the Logarithmic decrement
-    whirl_values : array
-        Array with the whirl values (0, 0.5 or 1)
-
-    Returns
-    -------
-    ax : matplotlib axes
-        Returns the matplotlib axes object with the plot
-        if plot_type == "matplotlib"
-    bk_ax : bokeh axes
-        Returns the bokeh axes object with the plot
-        if plot_type == "bokeh"
-    """
-
     def __init__(self, speed_range, wd, log_dec, whirl_values):
         self.speed_range = speed_range
         self.wd = wd
@@ -103,31 +75,6 @@ class CampbellResults:
         self.whirl_values = whirl_values
 
     def _plot_matplotlib(self, harmonics=[1], fig=None, ax=None, **kwargs):
-        """
-        Method to create Campbell Diagram figure using Matplotlib library.
-
-        Parameters
-        ----------
-        harmonics: list, optional
-            List withe the harmonics to be plotted.
-            The default is to plot 1x.
-        fig : matplotlib figure, optional
-            Figure in which the plot will be drawn
-            Default is None
-        ax : matplotlib plotting axes, optional
-            Axes which the plot will take to draw.
-            Default is None
-        kwargs : optional
-            Additional key word arguments can be passed to change
-            the plot (e.g. linestyle='--')
-
-        Returns
-        -------
-        fig : matplotlib figure
-            A figure with the Campbell Diagram plot
-        ax : matplotlib plotting axes
-            The axes from Campbell Diagram plot
-        """
         if fig is None and ax is None:
             fig, ax = plt.subplots()
 
@@ -226,23 +173,6 @@ class CampbellResults:
         return fig, ax
 
     def _plot_bokeh(self, harmonics=[1], **kwargs):
-        """
-        Method to create Campbell Diagram figure using Bokeh library.
-
-        Parameters
-        ----------
-        harmonics: list, optional
-            List withe the harmonics to be plotted.
-            The default is to plot 1x.
-        kwargs : optional
-            Additional key word arguments can be passed to change
-            the plot (e.g. linestyle='--')
-
-        Returns
-        -------
-        camp : Bokeh axes
-            The bokeh axes object with the plot
-        """
         wd = self.wd
         num_frequencies = wd.shape[1]
         log_dec = self.log_dec
@@ -386,28 +316,22 @@ class CampbellResults:
 
     def plot(self, *args, plot_type="bokeh", **kwargs):
         """Plot campbell results.
-
         Parameters
         ----------
-        args: optional
-            harmonics : list, optional
-                List with the harmonics to be plotted.
-                The default is to plot 1x.
+        harmonics: list, optional
+            List withe the harmonics to be plotted.
+            The default is to plot 1x.
         plot_type: str
             Matplotlib or bokeh.
-            The default is bokeh
-        kwargs : optional
-            Additional key word arguments can be passed to change
-            the plot (e.g. linestyle='--')
-
+            The default is matplotlib
+        fig : matplotlib figure, optional
+            Figure to insert axes with log_dec colorbar.
+        ax : matplotlib axes, optional
+            Axes in which the plot will be drawn.
         Returns
         -------
         ax : matplotlib axes
-            Returns the matplotlib axes object with the plot
-            if plot_type == "matplotlib"
-        bk_ax : bokeh axes
-            Returns the bokeh axes object with the plot
-            if plot_type == "bokeh"
+            Returns the axes object with the plot.
         """
         if plot_type == "matplotlib":
             return self._plot_matplotlib(*args, **kwargs)
@@ -418,34 +342,6 @@ class CampbellResults:
 
 
 class FrequencyResponseResults:
-    """Class used to store results and provide plots for Frequency Response.
-
-    Two options for plooting are available: Matplotlib and Bokeh. The user
-    chooses between them using the attribute plot_type. The default is bokeh
-
-    Parameters
-    ----------
-    freq_resp : array
-        Array with the transfer matrix
-    speed_range : array
-        Array with the speed range in rad/s.
-    magnitude : array
-        Array with the frequencies, magnitude (dB) of the frequency
-        response for each pair input/output
-    phase : array
-        Array with the frequencies, phase of the frequency
-        response for each pair input/output
-
-    Returns
-    -------
-    ax : matplotlib axes
-        Returns the matplotlib axes object with the plot
-        if plot_type == "matplotlib"
-    bk_ax : bokeh axes
-        Returns the bokeh axes object with the plot
-        if plot_type == "bokeh"
-    """
-
     def __init__(self, freq_resp, speed_range, magnitude, phase):
         self.freq_resp = freq_resp
         self.speed_range = speed_range
@@ -455,8 +351,7 @@ class FrequencyResponseResults:
     def plot_magnitude_matplotlib(self, inp, out, ax=None, units="mic-pk-pk", **kwargs):
         """Plot frequency response.
         This method plots the frequency response magnitude given an output and
-        an input using Matplotlib.
-
+        an input.
         Parameters
         ----------
         inp : int
@@ -464,19 +359,19 @@ class FrequencyResponseResults:
         out : int
             Output.
         ax : matplotlib.axes, optional
-            Matplotlib axes to plot the magnitude.
+            Matplotlib axes where the phase will be plotted.
             If None creates a new.
-        units : str
-            Unit system
-            Default is "mic-pk-pk"
         kwargs : optional
             Additional key word arguments can be passed to change
             the plot (e.g. linestyle='--')
-
         Returns
         -------
         ax : matplotlib.axes
-            Matplotlib axes with magnitude plot.
+            Matplotlib axes with amplitude plot.
+        mag_plot : bokeh plot axes
+            Bokeh plot axes with amplitude plot.
+        Examples
+        --------
         """
         if ax is None:
             ax = plt.gca()
@@ -496,28 +391,30 @@ class FrequencyResponseResults:
 
         return ax
 
-    def plot_magnitude_bokeh(self, inp, out, units="mic-pk-pk", **kwargs):
+    def plot_magnitude_bokeh(self, inp, out, units="m", **kwargs):
         """Plot frequency response.
         This method plots the frequency response magnitude given an output and
-        an input using Bokeh.
-
+        an input.
         Parameters
         ----------
         inp : int
             Input.
         out : int
             Output.
-        units : str
-            Unit system
-            Default is "mic-pk-pk"
+        ax : matplotlib.axes, optional
+            Matplotlib axes where the phase will be plotted.
+            If None creates a new.
         kwargs : optional
             Additional key word arguments can be passed to change
             the plot (e.g. linestyle='--')
-
         Returns
         -------
+        ax : matplotlib.axes
+            Matplotlib axes with amplitude plot.
         mag_plot : bokeh plot axes
-            Bokeh plot axes with magnitude plot.
+            Bokeh plot axes with amplitude plot.
+        Examples
+        --------
         """
         frequency_range = self.speed_range
         mag = self.magnitude
@@ -525,7 +422,7 @@ class FrequencyResponseResults:
         if units == "m":
             y_axis_label = "Amplitude (m)"
         elif units == "mic-pk-pk":
-            y_axis_label = "Amplitude ($\mu$ pk-pk)"
+            y_axis_label = "Amplitude (\mu pk-pk)"
         else:
             y_axis_label = "Amplitude (dB)"
 
@@ -556,8 +453,7 @@ class FrequencyResponseResults:
     def plot_phase_matplotlib(self, inp, out, ax=None, **kwargs):
         """Plot frequency response.
         This method plots the frequency response phase given an output and
-        an input using Matplotlib.
-
+        an input.
         Parameters
         ----------
         inp : int
@@ -570,11 +466,14 @@ class FrequencyResponseResults:
         kwargs : optional
             Additional key word arguments can be passed to change
             the plot (e.g. linestyle='--')
-
         Returns
         -------
         ax : matplotlib.axes
             Matplotlib axes with phase plot.
+        phase_plot : bokeh plot axes
+            Bokeh plot axes with phase plot.
+        Examples
+        --------
         """
         if ax is None:
             ax = plt.gca()
@@ -596,22 +495,27 @@ class FrequencyResponseResults:
     def plot_phase_bokeh(self, inp, out, **kwargs):
         """Plot frequency response.
         This method plots the frequency response phase given an output and
-        an input using bokeh.
-
+        an input.
         Parameters
         ----------
         inp : int
             Input.
         out : int
             Output.
+        ax : matplotlib.axes, optional
+            Matplotlib axes where the phase will be plotted.
+            If None creates a new.
         kwargs : optional
             Additional key word arguments can be passed to change
             the plot (e.g. linestyle='--')
-
         Returns
         -------
+        ax : matplotlib.axes
+            Matplotlib axes with phase plot.
         phase_plot : bokeh plot axes
             Bokeh plot axes with phase plot.
+        Examples
+        --------
         """
         frequency_range = self.speed_range
         phase = self.phase
@@ -643,30 +547,34 @@ class FrequencyResponseResults:
     def _plot_matplotlib(self, inp, out, ax0=None, ax1=None, **kwargs):
         """Plot frequency response.
         This method plots the frequency response given
-        an output and an input using Matplotib.
-
+        an output and an input.
         Parameters
         ----------
         inp : int
             Input.
         out : int
             Output.
-        ax0 : matplotlib.axes, optional
-            Matplotlib axes where the magnitude will be plotted.
+        ax0 : matplotlib.axes, bokeh plot axes optional
+            Matplotlib and bokeh plot axes where the amplitude will be plotted.
             If None creates a new.
-        ax1 : matplotlib.axes, optional
-            Matplotlib axes where the phase will be plotted.
+        ax1 : matplotlib.axes, bokeh plot axes optional
+            Matplotlib and bokeh plot axes where the phase will be plotted.
             If None creates a new.
         kwargs : optional
             Additional key word arguments can be passed to change
             the plot (e.g. linestyle='--')
-
         Returns
         -------
         ax0 : matplotlib.axes
             Matplotlib axes with amplitude plot.
         ax1 : matplotlib.axes
             Matplotlib axes with phase plot.
+        bk_ax0 : bokeh plot axes
+            Bokeh plot axes with amplitude plot
+        bk_ax1 : bokeh plot axes
+            Bokeh plot axes with phase plot
+        Examples
+        --------
         """
         if ax0 is None and ax1 is None:
             fig, (ax0, ax1) = plt.subplots(2)
@@ -682,28 +590,34 @@ class FrequencyResponseResults:
     def _plot_bokeh(self, inp, out, ax0=None, ax1=None, **kwargs):
         """Plot frequency response.
         This method plots the frequency response given
-        an output and an input using Bokeh.
-
+        an output and an input.
         Parameters
         ----------
         inp : int
             Input.
         out : int
             Output.
-        ax0 : bokeh axes, optional
-            Bokeh plot axes where the magnitude will be plotted.
+        ax0 : matplotlib.axes, bokeh plot axes optional
+            Matplotlib and bokeh plot axes where the amplitude will be plotted.
             If None creates a new.
-        ax1 : bokeh axes, optional
-            Bokeh plot axes where the phase will be plotted.
+        ax1 : matplotlib.axes, bokeh plot axes optional
+            Matplotlib and bokeh plot axes where the phase will be plotted.
             If None creates a new.
         kwargs : optional
             Additional key word arguments can be passed to change
             the plot (e.g. linestyle='--')
-
         Returns
         -------
-        grid_plots : bokeh column
-            Bokeh column with magnitude and phase plots.
+        ax0 : matplotlib.axes
+            Matplotlib axes with amplitude plot.
+        ax1 : matplotlib.axes
+            Matplotlib axes with phase plot.
+        bk_ax0 : bokeh plot axes
+            Bokeh plot axes with amplitude plot
+        bk_ax1 : bokeh plot axes
+            Bokeh plot axes with phase plot
+        Examples
+        --------
         """
         # bokeh plot axes
         bk_ax0 = self.plot_magnitude_bokeh(inp, out, ax=ax0)
@@ -719,33 +633,36 @@ class FrequencyResponseResults:
         """Plot frequency response.
         This method plots the frequency response given
         an output and an input.
-
         Parameters
         ----------
         inp : int
             Input.
         out : int
             Output.
-        args : optional
-            Additional bokeh plot axes or matplolib.axes
         plot_type: str
             Matplotlib or bokeh.
-            The default is bokeh
+            The default is matplotlib
+        ax0 : matplotlib.axes, bokeh plot axes optional
+            Matplotlib and bokeh plot axes where the amplitude will be plotted.
+            If None creates a new.
+        ax1 : matplotlib.axes, bokeh plot axes optional
+            Matplotlib and bokeh plot axes where the phase will be plotted.
+            If None creates a new.
         kwargs : optional
             Additional key word arguments can be passed to change
             the plot (e.g. linestyle='--')
-
         Returns
         -------
         ax0 : matplotlib.axes
             Matplotlib axes with amplitude plot.
-            if plot_type == "matplotlib"
         ax1 : matplotlib.axes
             Matplotlib axes with phase plot.
-            if plot_type == "matplotlib"
-        grid_plots : bokeh column
-            Bokeh column with amplitude and phase plot
-            if plot_type == "bokeh"
+        bk_ax0 : bokeh plot axes
+            Bokeh plot axes with amplitude plot
+        bk_ax1 : bokeh plot axes
+            Bokeh plot axes with phase plot
+        Examples
+        --------
         """
         if plot_type == "matplotlib":
             return self._plot_matplotlib(inp, out, *args, **kwargs)
@@ -758,7 +675,6 @@ class FrequencyResponseResults:
         """Plot frequency response.
         This method plots the frequency response given
         an output and an input.
-
         Parameters
         ----------
         outs : list
@@ -768,7 +684,6 @@ class FrequencyResponseResults:
         ax : array with matplotlib.axes, optional
             Matplotlib axes array created with plt.subplots.
             It needs to have a shape of (2*inputs, outputs).
-
         Returns
         -------
         ax : array with matplotlib.axes, optional
@@ -797,36 +712,6 @@ class FrequencyResponseResults:
 
 
 class ForcedResponseResults:
-    """Class used to store results and provide plots for Unbalance and Forced
-    Response analysis.
-
-    Two options for plooting are available: Matplotlib and Bokeh. The user
-    chooses between them using the attribute plot_type. The default is bokeh
-
-    Parameters
-    ----------
-    force_resp : array
-        Array with the force response for each node for each frequency
-    speed_range : array
-        Array with the frequencies
-    magnitude : array
-        Magnitude (dB) of the frequency response for node for each frequency
-    phase : array
-        Phase of the frequency response for node for each frequency
-
-    Returns
-    -------
-    ax0 : matplotlib.axes
-        Matplotlib axes with magnitude plot.
-        if plot_type == "matplotlib"
-    ax1 : matplotlib.axes
-        Matplotlib axes with phase plot.
-        if plot_type == "matplotlib"
-    grid_plots : bokeh column
-        Bokeh colum with magnitude and phase plot
-        if plot_type == "bokeh"
-    """
-
     def __init__(self, forced_resp, speed_range, magnitude, phase):
         self.forced_resp = forced_resp
         self.speed_range = speed_range
@@ -836,18 +721,17 @@ class ForcedResponseResults:
     def plot_magnitude_matplotlib(self, dof, ax=None, units="m", **kwargs):
         """Plot frequency response.
         This method plots the frequency response magnitude given an output and
-        an input using Matplotlib.
+        an input.
 
         Parameters
         ----------
         dof : int
             Degree of freedom.
         ax : matplotlib.axes, optional
-            Matplotlib axes where the magnitude will be plotted.
+            Matplotlib axes where the phase will be plotted.
             If None creates a new.
         units : str
             Units to plot the magnitude ('m' or 'mic-pk-pk')
-            Default is 'm'
         kwargs : optional
             Additional key word arguments can be passed to change
             the plot (e.g. linestyle='--')
@@ -855,7 +739,10 @@ class ForcedResponseResults:
         Returns
         -------
         ax : matplotlib.axes
-            Matplotlib axes with magnitude plot.
+            Matplotlib axes with phase plot.
+
+        Examples
+        --------
         """
         if ax is None:
             ax = plt.gca()
@@ -883,7 +770,7 @@ class ForcedResponseResults:
     def plot_magnitude_bokeh(self, dof, units="m", **kwargs):
         """Plot frequency response.
         This method plots the frequency response magnitude given an output and
-        an input using Bokeh.
+        an input.
 
         Parameters
         ----------
@@ -891,7 +778,6 @@ class ForcedResponseResults:
             Degree of freedom.
         units : str
             Units to plot the magnitude ('m' or 'mic-pk-pk')
-            Default is 'm'
         kwargs : optional
             Additional key word arguments can be passed to change
             the plot (e.g. linestyle='--')
@@ -900,6 +786,9 @@ class ForcedResponseResults:
         -------
         mag_plot : bokeh axes
             bokeh axes with magnitude plot
+
+        Examples
+        --------
         """
         frequency_range = self.speed_range
         mag = self.magnitude
@@ -938,7 +827,7 @@ class ForcedResponseResults:
     def plot_phase_matplotlib(self, dof, ax=None, **kwargs):
         """Plot frequency response.
         This method plots the frequency response phase given an output and
-        an input using Matplotlib.
+        an input.
 
         Parameters
         ----------
@@ -955,6 +844,9 @@ class ForcedResponseResults:
         -------
         ax : matplotlib.axes
             Matplotlib axes with phase plot.
+
+        Examples
+        --------
         """
         if ax is None:
             ax = plt.gca()
@@ -977,7 +869,7 @@ class ForcedResponseResults:
     def plot_phase_bokeh(self, dof, **kwargs):
         """Plot frequency response.
         This method plots the frequency response phase given an output and
-        an input using Bokeh.
+        an input.
 
         Parameters
         ----------
@@ -991,6 +883,9 @@ class ForcedResponseResults:
         -------
         phase_plot : bokeh axes
             Bokeh axes with phase plot
+
+        Examples
+        --------
         """
         frequency_range = self.speed_range
         phase = self.phase
@@ -1019,31 +914,7 @@ class ForcedResponseResults:
         return phase_plot
 
     def _plot_matplotlib(self, dof, ax0=None, ax1=None, **kwargs):
-        """Plot frequency response.
-        This method plots the frequency response magnitude and phase given
-        an output and an input using Matplotlib.
-
-        Parameters
-        ----------
-        dof : int
-            Degree of freedom.
-        ax0 : matplotlib.axes, optional
-            Matplotlib axes where the magnitude will be plotted.
-            If None creates a new.
-        ax1 : matplotlib.axes, optional
-            Matplotlib axes where the phase will be plotted.
-            If None creates a new.            
-        kwargs : optional
-            Additional key word arguments can be passed to change
-            the plot (e.g. linestyle='--')
 
-        Returns
-        -------
-        ax0 : matplotlib.axes
-            Matplotlib axes with magnitude plot.
-        ax1 : matplotlib.axes
-            Matplotlib axes with phase plot.
-        """
         if ax0 is None and ax1 is None:
             fig, (ax0, ax1) = plt.subplots(2)
 
@@ -1059,23 +930,6 @@ class ForcedResponseResults:
         return ax0, ax1
 
     def _plot_bokeh(self, dof, **kwargs):
-        """Plot frequency response.
-        This method plots the frequency response magnitude and phase given
-        an output and an input using Bokeh.
-
-        Parameters
-        ----------
-        dof : int
-            Degree of freedom.
-        kwargs : optional
-            Additional key word arguments can be passed to change
-            the plot (e.g. linestyle='--')
-
-        Returns
-        -------
-        grid_plots : bokeh column
-            Bokeh colum with magnitude and phase plot
-        """
         # bokeh plot axes
         bk_ax0 = self.plot_magnitude_bokeh(dof, **kwargs)
         bk_ax1 = self.plot_phase_bokeh(dof, **kwargs)
@@ -1088,7 +942,8 @@ class ForcedResponseResults:
 
     def plot(self, dof, plot_type="bokeh", **kwargs):
         """Plot frequency response.
-        This method plots the frequency response given an output and an input.
+        This method plots the frequency response given
+        an output and an input.
 
         Parameters
         ----------
@@ -1096,22 +951,16 @@ class ForcedResponseResults:
             Degree of freedom.
         plot_type: str
             Matplotlib or bokeh.
-            The default is bokeh
+            The default is matplotlib
         kwargs : optional
             Additional key word arguments can be passed to change
             the plot (e.g. linestyle='--')
 
         Returns
         -------
-        ax0 : matplotlib.axes
-            Matplotlib axes with magnitude plot.
-            if plot_type == "matplotlib"
-        ax1 : matplotlib.axes
-            Matplotlib axes with phase plot.
-            if plot_type == "matplotlib"
-        grid_plots : bokeh column
-            Bokeh colum with magnitude and phase plot
-            if plot_type == "bokeh"
+
+        Examples
+        --------
         """
         if plot_type == "matplotlib":
             return self._plot_matplotlib(dof, **kwargs)
@@ -1128,8 +977,8 @@ class ModeShapeResults:
 
     Parameters
     ----------
-    modes : array
-        Array of eigenvectors
+    modes : list
+        list of eigenvectors
     ndof : int
         Number of degrees of freedom of the system
     nodes : list
@@ -1149,7 +998,10 @@ class ModeShapeResults:
 
     Returns
     -------
-    A graphic object
+        A graphic object
+
+    Example
+    -------
     """
 
     def __init__(
@@ -1158,7 +1010,7 @@ class ModeShapeResults:
         ndof,
         nodes,
         nodes_pos,
-        shaft_elements_length,
+        elements_length,
         w,
         wd,
         log_dec,
@@ -1168,7 +1020,7 @@ class ModeShapeResults:
         self.ndof = ndof
         self.nodes = nodes
         self.nodes_pos = nodes_pos
-        self.shaft_elements_length = shaft_elements_length
+        self.elements_length = elements_length
         self.w = w
         self.wd = wd
         self.log_dec = log_dec
@@ -1202,11 +1054,14 @@ class ModeShapeResults:
             axial location of each orbit
         nn : int
             number of points to plot lines
+
+        Example
+        -------
         """
         evec0 = self.modes[:, mode]
         nodes = self.nodes
         nodes_pos = self.nodes_pos
-        shaft_elements_length = self.shaft_elements_length
+        elements_length = self.elements_length
 
         modex = evec0[0::4]
         modey = evec0[1::4]
@@ -1254,7 +1109,7 @@ class ModeShapeResults:
         N3 = 3 * zeta ** 2 - 2 * zeta ** 3
         N4 = -zeta ** 2 + zeta ** 3
 
-        for Le, n in zip(shaft_elements_length, nodes):
+        for Le, n in zip(elements_length, nodes):
             node_pos = nodes_pos[n]
             Nx = np.hstack((N1, Le * N2, N3, Le * N4))
             Ny = np.hstack((N1, -Le * N2, N3, -Le * N4))
@@ -1273,7 +1128,7 @@ class ModeShapeResults:
 
     def plot(self, mode=None, evec=None, fig=None, ax=None):
         """
-        Method that plots the mode shapes.
+        Method that calculate the arrays describing the mode shapes.
 
         Parameters
         ----------
@@ -1293,6 +1148,9 @@ class ModeShapeResults:
             Returns the figure object with the plot.
         ax : matplotlib axes
             Returns the axes object with the plot.
+
+        Example
+        -------
         """
         if ax is None:
             fig = plt.figure()
@@ -1342,46 +1200,6 @@ class ModeShapeResults:
 
 
 class StaticResults:
-    """Class used to store results and provide plots for Static Analysis.
-
-    This class plots free-body diagram, deformed shaft, shearing
-    force diagram and bending moment diagram.
-
-    Parameters
-    ----------
-    disp_y : array
-        shaft displacement in y direction
-    Vx : array
-        shearing force array
-    Bm : array
-        bending moment array
-    df_shaft : dataframe
-        shaft dataframe
-    df_disks : dataframe
-        disks dataframe
-    df_bearings : dataframe
-        bearing dataframe
-    nodes : list
-        list of nodes numbers
-    nodes_pos : list
-        list of nodes positions
-    Vx_axis : array
-        X axis for displaying shearing force
-    force_data : dict
-        A dictionary containing the information about:
-        Static displacement vector,
-        Shearing force vector,
-        Bending moment vector,
-        Shaft total weight,
-        Disks forces,
-        Bearings reaction forces
-
-    Returns
-    -------
-    fig : bokeh figures
-        Bokeh figure with Static Analysis plots depending on which method
-        is called.
-    """
     def __init__(
         self,
         disp_y,
@@ -1407,85 +1225,87 @@ class StaticResults:
         self.Vx_axis = Vx_axis
         self.force_data = force_data
 
-    def plot_deformation(self):
-        """Plot the shaft static deformation.
-
+    def plot(self):
+        """Plot static analysis graphs.
         This method plots:
-            deformed shaft
+            free-body diagram,
+            deformed shaft,
+            shearing force diagram,
+            bending moment diagram.
 
         Parameters
         ----------
 
         Returns
         -------
-        fig : bokeh figure
-            Bokeh figure with static deformation plot
+        grid_plots : bokeh.gridplot
+        --------
         """
         source = ColumnDataSource(
             data=dict(
-                x=self.nodes_pos, y0=self.disp_y * 1000, y1=[0] * len(self.nodes_pos)
+                x0=self.nodes_pos, y0=self.disp_y * 1000, y1=[0] * len(self.nodes_pos)
             )
         )
 
+        TOOLS = "pan,wheel_zoom,box_zoom,reset,save,box_select,hover"
         TOOLTIPS = [
-            ("Shaft lenght:", "@x"),
+            ("Shaft lenght:", "@x0"),
             ("Underformed:", "@y1"),
             ("Displacement:", "@y0"),
         ]
 
         # create displacement plot
-        fig = figure(
-            tools="pan, wheel_zoom, box_zoom, reset, save, box_select, hover",
+        disp_graph = figure(
+            tools=TOOLS,
             tooltips=TOOLTIPS,
-            width=640,
-            height=480,
+            width=900,
+            height=450,
             title="Static Analysis",
-            x_axis_label="Shaft lenght",
-            y_axis_label="Lateral displacement",
+            x_axis_label="shaft lenght",
+            y_axis_label="lateral displacement",
         )
-        fig.xaxis.axis_label_text_font_size = "14pt"
-        fig.yaxis.axis_label_text_font_size = "14pt"
-        fig.title.text_font_size = "14pt"
+        disp_graph.xaxis.axis_label_text_font_size = "14pt"
+        disp_graph.yaxis.axis_label_text_font_size = "14pt"
 
         interpolated = interpolate.interp1d(
-            source.data["x"], source.data["y0"], kind="cubic"
+            source.data["x0"], source.data["y0"], kind="cubic"
         )
         xnew = np.linspace(
-            source.data["x"][0],
-            source.data["x"][-1],
+            source.data["x0"][0],
+            source.data["x0"][-1],
             num=len(self.nodes_pos) * 20,
             endpoint=True,
         )
 
         ynew = interpolated(xnew)
-        auxsource = ColumnDataSource(data=dict(x=xnew, y0=ynew, y1=[0] * len(xnew)))
+        auxsource = ColumnDataSource(data=dict(x0=xnew, y0=ynew, y1=[0] * len(xnew)))
 
-        fig.line(
-            "x",
+        disp_graph.line(
+            "x0",
             "y0",
             source=auxsource,
             legend="Deformed shaft",
             line_width=3,
             line_color=bokeh_colors[9],
         )
-        fig.circle(
-            "x",
+        disp_graph.circle(
+            "x0",
             "y0",
             source=source,
             legend="Deformed shaft",
             size=8,
             fill_color=bokeh_colors[9],
         )
-        fig.line(
-            "x",
+        disp_graph.line(
+            "x0",
             "y1",
             source=source,
             legend="underformed shaft",
             line_width=3,
             line_color=bokeh_colors[0],
         )
-        fig.circle(
-            "x",
+        disp_graph.circle(
+            "x0",
             "y1",
             source=source,
             legend="underformed shaft",
@@ -1493,50 +1313,28 @@ class StaticResults:
             fill_color=bokeh_colors[0],
         )
 
-        return fig
-
-    def plot_free_body_diagram(self):
-        """Plot the rotor free-body diagram.
-
-        This method plots:
-            free-body diagram.
-
-        Parameters
-        ----------
-
-        Returns
-        -------
-        fig : bokeh figure
-            Bokeh figure with the free-body diagram plot
-        """
-        source = ColumnDataSource(
-            data=dict(
-                x=self.nodes_pos, y0=self.disp_y * 1000, y1=[0] * len(self.nodes_pos)
-            )
-        )
-
+        # create a new plot for free body diagram (FDB)
         y_start = 5.0
         sh_weight = sum(self.df_shaft["m"].values) * 9.8065
 
         shaft_end = self.nodes_pos[-1]
-        fig = figure(
-            tools="pan, wheel_zoom, box_zoom, reset, save, box_select, hover",
-            width=640,
-            height=480,
+        FBD = figure(
+            tools=TOOLS,
+            width=900,
+            height=450,
             title="Free-Body Diagram",
             x_axis_label="shaft lenght",
             x_range=[-0.1 * shaft_end, 1.1 * shaft_end],
             y_range=[-3 * y_start, 3 * y_start],
         )
-        fig.yaxis.visible = False
-        fig.xaxis.axis_label_text_font_size = "14pt"
-        fig.title.text_font_size = "14pt"
+        FBD.yaxis.visible = False
+        FBD.xaxis.axis_label_text_font_size = "14pt"
 
-        fig.line("x", "y1", source=source, line_width=5, line_color=bokeh_colors[0])
+        FBD.line("x0", "y1", source=source, line_width=5, line_color=bokeh_colors[0])
 
-        # fig - plot arrows indicating shaft weight distribution
+        # FBD - plot arrows indicating shaft weight distribution
         text = str("%.1f" % sh_weight)
-        fig.line(
+        FBD.line(
             x=self.nodes_pos,
             y=[y_start] * len(self.nodes_pos),
             line_width=2,
@@ -1547,7 +1345,7 @@ class StaticResults:
         fin = self.nodes_pos[-1]
         arrows_list = np.arange(ini, 1.01 * fin, fin / 5.0)
         for node in arrows_list:
-            fig.add_layout(
+            FBD.add_layout(
                 Arrow(
                     end=NormalHead(
                         fill_color=bokeh_colors[2],
@@ -1563,7 +1361,7 @@ class StaticResults:
                 )
             )
 
-        fig.add_layout(
+        FBD.add_layout(
             Label(
                 x=self.nodes_pos[0],
                 y=y_start,
@@ -1576,11 +1374,11 @@ class StaticResults:
             )
         )
 
-        # fig - calculate the reaction force of bearings and plot arrows
+        # FBD - calculate the reaction force of bearings and plot arrows
         for i, node in enumerate(self.df_bearings["n"]):
             Fb = -self.disp_y[node] * self.df_bearings.loc[i, "kyy"].coefficient[0]
             text = str("%.1f" % Fb)
-            fig.add_layout(
+            FBD.add_layout(
                 Arrow(
                     end=NormalHead(
                         fill_color=bokeh_colors[6],
@@ -1595,7 +1393,7 @@ class StaticResults:
                     y_end=0,
                 )
             )
-            fig.add_layout(
+            FBD.add_layout(
                 Label(
                     x=self.nodes_pos[node],
                     y=-2 * y_start,
@@ -1609,12 +1407,12 @@ class StaticResults:
                 )
             )
 
-        # fig - plot arrows indicating disk weight
+        # FBD - plot arrows indicating disk weight
         if len(self.df_disks) != 0:
             for i, node in enumerate(self.df_disks["n"]):
                 Fd = self.df_disks.loc[i, "m"] * 9.8065
                 text = str("%.1f" % Fd)
-                fig.add_layout(
+                FBD.add_layout(
                     Arrow(
                         end=NormalHead(
                             fill_color=bokeh_colors[9],
@@ -1629,7 +1427,7 @@ class StaticResults:
                         y_end=0,
                     )
                 )
-                fig.add_layout(
+                FBD.add_layout(
                     Label(
                         x=self.nodes_pos[node],
                         y=2 * y_start,
@@ -1643,44 +1441,27 @@ class StaticResults:
                     )
                 )
 
-        return fig
-
-    def plot_shearing_force(self):
-        """Plot the rotor shearing force diagram.
-
-        This method plots:
-            shearing force diagram.
-
-        Parameters
-        ----------
-
-        Returns
-        -------
-        fig : bokeh figure
-            Bokeh figure with the shearing force diagram plot
-        """
-        shaft_end = self.nodes_pos[-1]
+        # Shearing Force Diagram plot (SF)
         source_SF = ColumnDataSource(data=dict(x=self.Vx_axis, y=self.Vx))
         TOOLTIPS_SF = [("Shearing Force:", "@y")]
-        fig = figure(
-            tools="pan, wheel_zoom, box_zoom, reset, save, box_select, hover",
+        SF = figure(
+            tools=TOOLS,
             tooltips=TOOLTIPS_SF,
-            width=640,
-            height=480,
+            width=900,
+            height=450,
             title="Shearing Force Diagram",
             x_axis_label="Shaft lenght",
             y_axis_label="Force",
             x_range=[-0.1 * shaft_end, 1.1 * shaft_end],
         )
-        fig.xaxis.axis_label_text_font_size = "14pt"
-        fig.yaxis.axis_label_text_font_size = "14pt"
-        fig.title.text_font_size = "14pt"
+        SF.xaxis.axis_label_text_font_size = "14pt"
+        SF.yaxis.axis_label_text_font_size = "14pt"
 
-        fig.line("x", "y", source=source_SF, line_width=4, line_color=bokeh_colors[0])
-        fig.circle("x", "y", source=source_SF, size=8, fill_color=bokeh_colors[0])
+        SF.line("x", "y", source=source_SF, line_width=4, line_color=bokeh_colors[0])
+        SF.circle("x", "y", source=source_SF, size=8, fill_color=bokeh_colors[0])
 
-        # fig - plot centerline
-        fig.line(
+        # SF - plot centerline
+        SF.line(
             [-0.1 * shaft_end, 1.1 * shaft_end],
             [0, 0],
             line_width=3,
@@ -1688,38 +1469,21 @@ class StaticResults:
             line_color=bokeh_colors[0],
         )
 
-        return fig
-
-    def plot_bending_moment(self):
-        """Plot the rotor bending moment diagram.
-
-        This method plots:
-            bending moment diagram.
-
-        Parameters
-        ----------
-
-        Returns
-        -------
-        fig : bokeh figure
-            Bokeh figure with the bending moment diagram plot
-        """
-        shaft_end = self.nodes_pos[-1]
+        # Bending Moment Diagram plot (BM)
         source_BM = ColumnDataSource(data=dict(x=self.nodes_pos, y=self.Bm))
         TOOLTIPS_BM = [("Bending Moment:", "@y")]
-        fig = figure(
-            tools="pan, wheel_zoom, box_zoom, reset, save, box_select, hover",
+        BM = figure(
+            tools=TOOLS,
             tooltips=TOOLTIPS_BM,
-            width=640,
-            height=480,
+            width=900,
+            height=450,
             title="Bending Moment Diagram",
             x_axis_label="Shaft lenght",
             y_axis_label="Bending Moment",
             x_range=[-0.1 * shaft_end, 1.1 * shaft_end],
         )
-        fig.xaxis.axis_label_text_font_size = "14pt"
-        fig.yaxis.axis_label_text_font_size = "14pt"
-        fig.title.text_font_size = "14pt"
+        BM.xaxis.axis_label_text_font_size = "14pt"
+        BM.yaxis.axis_label_text_font_size = "14pt"
 
         i = 0
         while True:
@@ -1735,14 +1499,14 @@ class StaticResults:
 
             ynew_BM = interpolated_BM(xnew_BM)
             auxsource_BM = ColumnDataSource(data=dict(x=xnew_BM, y=ynew_BM))
-            fig.line(
+            BM.line(
                 "x", "y", source=auxsource_BM, line_width=4, line_color=bokeh_colors[0]
             )
             i += 2
-        fig.circle("x", "y", source=source_BM, size=8, fill_color=bokeh_colors[0])
+        BM.circle("x", "y", source=source_BM, size=8, fill_color=bokeh_colors[0])
 
-        # fig - plot centerline
-        fig.line(
+        # BM - plot centerline
+        BM.line(
             [-0.1 * shaft_end, 1.1 * shaft_end],
             [0, 0],
             line_width=3,
@@ -1750,31 +1514,12 @@ class StaticResults:
             line_color=bokeh_colors[0],
         )
 
-        return fig
+        grid_plots = gridplot([[FBD, SF], [disp_graph, BM]])
 
+        return grid_plots
 
-class ConvergenceResults:
-    """Class used to store results and provide plots for Convergence Analysis.
-
-    This class plots:
-        Natural Frequency vs Number of Elements
-        Relative Error vs Number of Elements
-
-    Parameters
-    ----------
-    el_num : array
-        Array with number of elements in each iteraction
-    eigv_arr : array
-        Array with the n'th natural frequency in each iteraction
-    error_arr : array
-        Array with the relative error in each iteraction
-
-    Returns
-    -------
-    plot : bokeh.gridplot
-        Bokeh column with Convergence Analysis plots
-    """
 
+class ConvergenceResults:
     def __init__(self, el_num, eigv_arr, error_arr):
         self.el_num = el_num
         self.eigv_arr = eigv_arr
@@ -1790,8 +1535,9 @@ class ConvergenceResults:
 
         Returns
         -------
-        plot : bokeh.gridplot
-            Bokeh column with Convergence Analysis plots
+        plot : bokeh.figure
+            Bokeh plot showing the results
+        --------
         """
         source = ColumnDataSource(
             data=dict(x0=self.el_num, y0=self.eigv_arr, y1=self.error_arr)
@@ -1842,33 +1588,6 @@ class ConvergenceResults:
 
 
 class TimeResponseResults:
-    """Class used to store results and provide plots for Time Response
-    Analysis.
-
-    This class takes the results from time response analysis and creates a
-    plot given a force and a time.
-
-    Parameters
-    ----------
-    t : array
-        Time values for the output.
-    yout : array
-        System response.
-    xout : array
-        Time evolution of the state vector.
-    dof : int
-        Degree of freedom
-
-    Returns
-    -------
-    ax : matplotlib.axes
-        Matplotlib axes with time response plot.
-        if plot_type == "matplotlib"
-    bk_ax : bokeh axes
-        Bokeh axes with time response plot
-        if plot_type == "bokeh"
-    """
-
     def __init__(self, t, yout, xout, dof):
         self.t = t
         self.yout = yout
@@ -1876,22 +1595,7 @@ class TimeResponseResults:
         self.dof = dof
 
     def _plot_matplotlib(self, ax=None):
-        """Plot time response.
-
-        This function will take a rotor object and plot its time response
-        using Matplotlib
-
-        Parameters
-        ----------
-        ax : matplotlib.axes
-            Matplotlib axes where time response will be plotted.
-            if None, creates a new one
 
-        Returns
-        -------
-        ax : matplotlib.axes
-            Matplotlib axes with time response plot.
-        """
         if ax is None:
             ax = plt.gca()
 
@@ -1917,20 +1621,7 @@ class TimeResponseResults:
         )
 
     def _plot_bokeh(self):
-        """Plot time response.
 
-        This function will take a rotor object and plot its time response
-        using Bokeh
-
-        Parameters
-        ----------
-
-        Returns
-        -------
-        bk_ax : bokeh axes
-            Bokeh axes with time response plot
-            if plot_type == "bokeh"
-        """
         if self.dof % 4 == 0:
             obs_dof = "x"
             amp = "m"
@@ -1968,24 +1659,24 @@ class TimeResponseResults:
         """Plot time response.
 
         This function will take a rotor object and plot its time response
+        given a force and a time.
 
         Parameters
         ----------
+        dof : int
+            Degree of freedom.
         plot_type: str
             Matplotlib or bokeh.
-            The default is bokeh
+            The default is matplotlib
         kwargs : optional
             Additional key word arguments can be passed to change
             the plot (e.g. linestyle='--')
 
         Returns
         -------
-        ax : matplotlib.axes
-            Matplotlib axes with time response plot.
-            if plot_type == "matplotlib"
-        bk_ax : bokeh axes
-            Bokeh axes with time response plot
-            if plot_type == "bokeh"
+
+        Examples
+        --------
         """
         if plot_type == "matplotlib":
             return self._plot_matplotlib(**kwargs)
diff --git a/ross/rotor_assembly.py b/ross/rotor_assembly.py
index 3a8b1883..e4fc8e06 100644
--- a/ross/rotor_assembly.py
+++ b/ross/rotor_assembly.py
@@ -203,7 +203,7 @@ class Rotor(object):
             if brg.__class__.__name__ == "SealElement" and brg.tag is None:
                 brg.tag = "Seal " + str(i)
 
-        self.shaft_elements = sorted(shaft_elements, key=lambda el: el.n)
+        self.shaft_elements = shaft_elements
         self.bearing_seal_elements = bearing_seal_elements
         self.disk_elements = disk_elements
         self.point_mass_elements = point_mass_elements
@@ -309,10 +309,12 @@ class Rotor(object):
         nodes_o_d.append(df_shaft["o_d"].iloc[-1])
         self.nodes_o_d = nodes_o_d
 
-        shaft_elements_length = list(df_shaft.groupby("n_l")["L"].min())
-        self.shaft_elements_length = shaft_elements_length
+        nodes_le = list(df_shaft.groupby("n_l")["L"].min())
+        nodes_le.append(df_shaft["L"].iloc[-1])
+        self.nodes_le = nodes_le
 
         self.nodes = list(range(len(self.nodes_pos)))
+        self.elements_length = [sh_el.L for sh_el in self.shaft_elements]
         self.L = nodes_pos[-1]
 
         # rotor mass can also be calculated with self.M()[::4, ::4].sum()
@@ -1845,7 +1847,7 @@ class Rotor(object):
             ndof=self.ndof,
             nodes=self.nodes,
             nodes_pos=self.nodes_pos,
-            shaft_elements_length=self.shaft_elements_length,
+            elements_length=self.elements_length,
             w=self.w,
             wd=self.wd,
             log_dec=self.log_dec,
diff --git a/ross/shaft_element.py b/ross/shaft_element.py
index ae0c8374..a47a397e 100644
--- a/ross/shaft_element.py
+++ b/ross/shaft_element.py
@@ -279,7 +279,7 @@ class ShaftElement(Element):
         >>> shaft2 = ShaftElement.load("ShaftElement.toml")
         >>> shaft2
         [ShaftElement(L=0.25, i_d=0.0, o_d=0.05, material='Steel', n=None)]
-        """
+        """ 
         shaft_elements = []
         with open("ShaftElement.toml", "r") as f:
             shaft_elements_dict = toml.load(f)
@@ -671,7 +671,7 @@ class ShaftElement(Element):
             bk_color = "yellow"
             legend = "Shaft - Slenderness Ratio < 1.6"
         else:
-            bk_color = self.material.color
+            bk_color = bokeh_colors[2]
             legend = "Shaft"
 
         source_u = ColumnDataSource(
@@ -1708,7 +1708,7 @@ class ShaftTaperedElement(Element):
             bk_color = "yellow"
             legend = "Shaft - Slenderness Ratio < 1.6"
         else:
-            bk_color = self.material.color
+            bk_color = bokeh_colors[2]
             legend = "Shaft"
 
         # bokeh plot - plot the shaft


Test patch:
diff --git a/ross/tests/test_bearing_seal_element.py b/ross/tests/test_bearing_seal_element.py
index d1f9c46c..51f73dc7 100644
--- a/ross/tests/test_bearing_seal_element.py
+++ b/ross/tests/test_bearing_seal_element.py
@@ -4,7 +4,11 @@ import os
 from numpy.testing import assert_allclose
 import math
 
-from ross.bearing_seal_element import BearingElement
+from ross.bearing_seal_element import (
+        BearingElement,
+        BallBearingElement,
+        RollerBearingElement
+)
 
 
 @pytest.fixture
@@ -228,3 +232,61 @@ def test_bearing_link_matrices():
 
     assert_allclose(b0.K(0), M)
     assert_allclose(b0.C(0), M)
+
+
+def test_ball_bearing_element():
+    n = 0
+    n_balls = 8
+    d_balls = 0.03
+    fs = 500.0
+    alpha = np.pi / 6
+    tag = "ballbearing"
+    ballbearing = BallBearingElement(
+            n=n,
+            n_balls=n_balls,
+            d_balls=d_balls,
+            fs=fs,
+            alpha=alpha,
+            tag=tag
+    )
+
+    M = np.zeros((2, 2))
+    K = np.array([[4.64168838e+07, 0.00000000e+00],
+                  [0.00000000e+00, 1.00906269e+08]])
+    C = np.array([[ 580.2110481 ,    0.        ],
+                  [   0.        , 1261.32836543]])
+    G = np.zeros((2, 2))
+
+    assert_allclose(ballbearing.M(), M)
+    assert_allclose(ballbearing.K(0), K)
+    assert_allclose(ballbearing.C(0), C)
+    assert_allclose(ballbearing.G(), G)
+
+
+def test_roller_bearing_element():
+    n = 0
+    n_rollers = 8
+    l_rollers = 0.03
+    fs = 500.0
+    alpha = np.pi / 6
+    tag = "rollerbearing"
+    rollerbearing = RollerBearingElement(
+            n=n,
+            n_rollers=n_rollers,
+            l_rollers=l_rollers,
+            fs=fs,
+            alpha=alpha,
+            tag=tag
+    )
+
+    M = np.zeros((2, 2))
+    K = np.array([[2.72821927e+08, 0.00000000e+00],
+                  [0.00000000e+00, 5.56779444e+08]])
+    C = np.array([[3410.27409251,    0.        ],
+                  [   0.        , 6959.74304593]])
+    G = np.zeros((2, 2))
+
+    assert_allclose(rollerbearing.M(), M)
+    assert_allclose(rollerbearing.K(0), K)
+    assert_allclose(rollerbearing.C(0), C)
+    assert_allclose(rollerbearing.G(), G)
diff --git a/ross/tests/test_rotor_assembly.py b/ross/tests/test_rotor_assembly.py
index b3200959..99cd3856 100644
--- a/ross/tests/test_rotor_assembly.py
+++ b/ross/tests/test_rotor_assembly.py
@@ -316,8 +316,8 @@ def test_rotor_attributes(rotor1, rotor3, rotor3_odd):
     assert rotor1.nodes_pos[0] == 0
     assert rotor1.nodes_pos[1] == 0.25
     assert rotor1.nodes_pos[-1] == 0.5
-    assert len(rotor3.shaft_elements_length) == 6
-    assert len(rotor3_odd.shaft_elements_length) == 7
+    assert len(rotor3.elements_length) == 6
+    assert len(rotor3_odd.elements_length) == 7
 
 
 def test_evects_sorted_rotor3(rotor3):