Repository: pybamm-team/PyBaMM (python)

Issue:
plot pybamm arrays
It would be nice to automate the plotting of `pybamm.Array` objects for doing things like the example below. I think just 1D is fine for now. More complicated will usually be of solution entries, so handled via `Quickplot`, but it would be nice to easily plot parameters that are functions. I guess 2D is useful since this depends on T and sto. Thoughts?

Edit: I know calling `.entries` isn't a huge deal, but think this would be easy to implement and avoid needing to import matplotlib 

```
import pybamm
import numpy as np
import matplotlib.pyplot as plt

param = pybamm.ParameterValues(chemistry=pybamm.parameter_sets.Ecker2015)
D_fun = param["Negative electrode diffusivity [m2.s-1]"]
x = pybamm.Array(np.linspace(0,1,100))
D = D_fun(x, 300)
plt.plot(x.entries,D.entries) 
```

Pull request:
Issue 992 plot arrays
# Description
Added basic plots for 1D and 2D arrays (`pybamm.plot` and `pybamm.plot2D`). I've left these as separate methods, but I guess they could be combined and the appropriate plot called depending on what arrays you pass. On the other hand, maybe it is better to be explicit. 

Also added `pybamm.linspace` and `pybamm.meshgrid` to quickly create some common arrays.

Fixes #992 

## Type of change

Please add a line in the relevant section of [CHANGELOG.md](https://github.com/pybamm-team/PyBaMM/blob/master/CHANGELOG.md) to document the change (include PR #) - note reverse order of PR #s. If necessary, also add to the list of breaking changes.

- [ ] New feature (non-breaking change which adds functionality)
- [ ] Optimization (back-end change that speeds up the code)
- [ ] Bug fix (non-breaking change which fixes an issue)


# Key checklist:

- [ ] No style issues: `$ flake8`
- [ ] All tests pass: `$ python run-tests.py --unit`
- [ ] The documentation builds: `$ cd docs` and then `$ make clean; make html`

You can run all three at once, using `$ python run-tests.py --quick`.

## Further checks:

- [ ] Code is commented, particularly in hard-to-understand areas
- [ ] Tests added that prove fix is effective or that feature works


Patch:
diff --git a/.flake8 b/.flake8
index bad0d2ef9..a69ea9bf0 100644
--- a/.flake8
+++ b/.flake8
@@ -12,7 +12,7 @@ exclude=
     share,
     pyvenv.cfg,
     third-party,
-    sundials-5.0.0,
+    KLU_module_deps,
 ignore=
     # False positive for white space before ':' on list slice
     # black should format these correctly
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 58ce82fd8..be8432402 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -3,6 +3,7 @@
 ## Features
 
 -   Added `BackwardIndefiniteIntegral` symbol ([#1014](https://github.com/pybamm-team/PyBaMM/pull/1014))
+-   Added `plot` and `plot2D` to enable easy plotting of `pybamm.Array` objects ([#1008](https://github.com/pybamm-team/PyBaMM/pull/1008))
 -   Added SEI film resistance as an option ([#994](https://github.com/pybamm-team/PyBaMM/pull/994))
 -   Added tab, edge, and surface cooling ([#965](https://github.com/pybamm-team/PyBaMM/pull/965))
 -   Added functionality to solver to automatically discretise a 0D model ([#947](https://github.com/pybamm-team/PyBaMM/pull/947))
diff --git a/docs/index.rst b/docs/index.rst
index 0ec015f49..daa168594 100644
--- a/docs/index.rst
+++ b/docs/index.rst
@@ -66,7 +66,7 @@ API documentation
    source/solvers/index
    source/experiments/index
    source/simulation
-   source/quick_plot
+   source/plotting/index
    source/util
    source/citations
    source/parameters_cli
diff --git a/docs/source/expression_tree/array.rst b/docs/source/expression_tree/array.rst
new file mode 100644
index 000000000..a84a856a1
--- /dev/null
+++ b/docs/source/expression_tree/array.rst
@@ -0,0 +1,9 @@
+Array
+=====
+
+.. autoclass:: pybamm.Array
+  :members:
+
+.. autofunction:: pybamm.linspace
+
+.. autofunction:: pybamm.meshgrid 
diff --git a/docs/source/expression_tree/index.rst b/docs/source/expression_tree/index.rst
index 3efd97fb3..0a6f3d757 100644
--- a/docs/source/expression_tree/index.rst
+++ b/docs/source/expression_tree/index.rst
@@ -8,6 +8,7 @@ Expression Tree
   variable
   independent_variable
   scalar
+  array
   matrix
   vector
   state_vector
diff --git a/docs/source/plotting/dynamic_plot.rst b/docs/source/plotting/dynamic_plot.rst
new file mode 100644
index 000000000..b363ab222
--- /dev/null
+++ b/docs/source/plotting/dynamic_plot.rst
@@ -0,0 +1,4 @@
+Dynamic Plot
+============
+
+.. autofunction:: pybamm.dynamic_plot
diff --git a/docs/source/plotting/index.rst b/docs/source/plotting/index.rst
new file mode 100644
index 000000000..705c332e8
--- /dev/null
+++ b/docs/source/plotting/index.rst
@@ -0,0 +1,9 @@
+Plotting
+========
+
+.. toctree::
+
+  quick_plot
+  dynamic_plot
+  plot
+  plot_2D
diff --git a/docs/source/plotting/plot.rst b/docs/source/plotting/plot.rst
new file mode 100644
index 000000000..011bd34fb
--- /dev/null
+++ b/docs/source/plotting/plot.rst
@@ -0,0 +1,4 @@
+Plot
+====
+
+.. autofunction:: pybamm.plot
diff --git a/docs/source/plotting/plot_2D.rst b/docs/source/plotting/plot_2D.rst
new file mode 100644
index 000000000..cd139ad2c
--- /dev/null
+++ b/docs/source/plotting/plot_2D.rst
@@ -0,0 +1,4 @@
+Plot 2D
+=======
+
+.. autofunction:: pybamm.plot2D
diff --git a/docs/source/quick_plot.rst b/docs/source/plotting/quick_plot.rst
similarity index 100%
rename from docs/source/quick_plot.rst
rename to docs/source/plotting/quick_plot.rst
diff --git a/pybamm/__init__.py b/pybamm/__init__.py
index e3c673e2c..6b7e7f8e5 100644
--- a/pybamm/__init__.py
+++ b/pybamm/__init__.py
@@ -68,7 +68,7 @@ from .citations import Citations, citations, print_citations
 from .expression_tree.symbol import *
 from .expression_tree.binary_operators import *
 from .expression_tree.concatenations import *
-from .expression_tree.array import Array
+from .expression_tree.array import Array, linspace, meshgrid
 from .expression_tree.matrix import Matrix
 from .expression_tree.unary_operators import *
 from .expression_tree.functions import *
@@ -221,10 +221,16 @@ from .experiments.experiment import Experiment
 from . import experiments
 
 #
-# other
+# Plotting
 #
-from .quick_plot import QuickPlot, dynamic_plot, ax_min, ax_max
+from .plotting.quick_plot import QuickPlot
+from .plotting.plot import plot
+from .plotting.plot2D import plot2D
+from .plotting.dynamic_plot import dynamic_plot
 
+#
+# Simulation
+#
 from .simulation import Simulation, load_sim, is_notebook
 
 #
diff --git a/pybamm/expression_tree/array.py b/pybamm/expression_tree/array.py
index 6acdd0795..e2d9869b9 100644
--- a/pybamm/expression_tree/array.py
+++ b/pybamm/expression_tree/array.py
@@ -101,3 +101,24 @@ class Array(pybamm.Symbol):
     def _base_evaluate(self, t=None, y=None, y_dot=None, inputs=None):
         """ See :meth:`pybamm.Symbol._base_evaluate()`. """
         return self._entries
+
+
+def linspace(start, stop, num=50, **kwargs):
+    """
+    Creates a linearly spaced array by calling `numpy.linspace` with keyword
+    arguments 'kwargs'. For a list of 'kwargs' see the
+    `numpy linspace documentation <https://tinyurl.com/yc4ne47x>`_
+    """
+    return pybamm.Array(np.linspace(start, stop, num, **kwargs))
+
+
+def meshgrid(x, y, **kwargs):
+    """
+    Return coordinate matrices as from coordinate vectors by calling
+    `numpy.meshgrid` with keyword arguments 'kwargs'. For a list of 'kwargs'
+    see the `numpy meshgrid documentation <https://tinyurl.com/y8azewrj>`_
+    """
+    [X, Y] = np.meshgrid(x.entries, y.entries)
+    X = pybamm.Array(X)
+    Y = pybamm.Array(Y)
+    return X, Y
diff --git a/pybamm/plotting/__init__.py b/pybamm/plotting/__init__.py
new file mode 100644
index 000000000..e69de29bb
diff --git a/pybamm/plotting/dynamic_plot.py b/pybamm/plotting/dynamic_plot.py
new file mode 100644
index 000000000..eb4d468be
--- /dev/null
+++ b/pybamm/plotting/dynamic_plot.py
@@ -0,0 +1,22 @@
+#
+# Method for creating a dynamic plot
+#
+import pybamm
+
+
+def dynamic_plot(*args, **kwargs):
+    """
+    Creates a :class:`pybamm.QuickPlot` object (with arguments 'args' and keyword
+    arguments 'kwargs') and then calls :meth:`pybamm.QuickPlot.dynamic_plot`.
+    The key-word argument 'testing' is passed to the 'dynamic_plot' method, not the
+    `QuickPlot` class.
+
+    Returns
+    -------
+    plot : :class:`pybamm.QuickPlot`
+        The 'QuickPlot' object that was created
+    """
+    kwargs_for_class = {k: v for k, v in kwargs.items() if k != "testing"}
+    plot = pybamm.QuickPlot(*args, **kwargs_for_class)
+    plot.dynamic_plot(kwargs.get("testing", False))
+    return plot
diff --git a/pybamm/plotting/plot.py b/pybamm/plotting/plot.py
new file mode 100644
index 000000000..d4176d838
--- /dev/null
+++ b/pybamm/plotting/plot.py
@@ -0,0 +1,44 @@
+#
+# Method for creating a 1D plot of pybamm arrays
+#
+import pybamm
+from .quick_plot import ax_min, ax_max
+
+
+def plot(x, y, xlabel=None, ylabel=None, title=None, testing=False, **kwargs):
+    """
+    Generate a simple 1D plot. Calls `matplotlib.pyplot.plot` with keyword
+    arguments 'kwargs'. For a list of 'kwargs' see the
+    `matplotlib plot documentation <https://tinyurl.com/ycblw9bx>`_
+
+    Parameters
+    ----------
+    x : :class:`pybamm.Array`
+        The array to plot on the x axis
+    y : :class:`pybamm.Array`
+        The array to plot on the y axis
+    xlabel : str, optional
+        The label for the x axis
+    ylabel : str, optional
+        The label for the y axis
+    testing : bool, optional
+        Whether to actually make the plot (turned off for unit tests)
+
+    """
+    import matplotlib.pyplot as plt
+
+    if not isinstance(x, pybamm.Array):
+        raise TypeError("x must be 'pybamm.Array'")
+    if not isinstance(y, pybamm.Array):
+        raise TypeError("y must be 'pybamm.Array'")
+
+    plt.plot(x.entries, y.entries, **kwargs)
+    plt.ylim([ax_min(y.entries), ax_max(y.entries)])
+    plt.xlabel(xlabel)
+    plt.ylabel(ylabel)
+    plt.title(title)
+
+    if not testing:  # pragma: no cover
+        plt.show()
+
+    return
diff --git a/pybamm/plotting/plot2D.py b/pybamm/plotting/plot2D.py
new file mode 100644
index 000000000..45f0d2b10
--- /dev/null
+++ b/pybamm/plotting/plot2D.py
@@ -0,0 +1,66 @@
+#
+# Method for creating a filled contour plot of pybamm arrays
+#
+import pybamm
+from .quick_plot import ax_min, ax_max
+
+
+def plot2D(x, y, z, xlabel=None, ylabel=None, title=None, testing=False, **kwargs):
+    """
+    Generate a simple 2D plot. Calls `matplotlib.pyplot.contourf` with keyword
+    arguments 'kwargs'.  For a list of 'kwargs' see the
+    `matplotlib contourf documentation <https://tinyurl.com/y8mnadtn>`_
+
+    Parameters
+    ----------
+    x : :class:`pybamm.Array`
+        The array to plot on the x axis. Can be of shape (M, N) or (N, 1)
+    y : :class:`pybamm.Array`
+        The array to plot on the y axis. Can be of shape (M, N)  or (M, 1)
+    z : :class:`pybamm.Array`
+        The array to plot on the z axis. Is of shape (M, N)
+    xlabel : str, optional
+        The label for the x axis
+    ylabel : str, optional
+        The label for the y axis
+    title : str, optional
+        The title for the plot
+    testing : bool, optional
+        Whether to actually make the plot (turned off for unit tests)
+
+    """
+    import matplotlib.pyplot as plt
+
+    if not isinstance(x, pybamm.Array):
+        raise TypeError("x must be 'pybamm.Array'")
+    if not isinstance(y, pybamm.Array):
+        raise TypeError("y must be 'pybamm.Array'")
+    if not isinstance(z, pybamm.Array):
+        raise TypeError("z must be 'pybamm.Array'")
+
+    # Get correct entries of x and y depending on shape
+    if x.shape == y.shape == z.shape:
+        x_entries = x.entries
+        y_entries = y.entries
+    else:
+        x_entries = x.entries[:, 0]
+        y_entries = y.entries[:, 0]
+
+    plt.contourf(
+        x_entries,
+        y_entries,
+        z.entries,
+        vmin=ax_min(z.entries),
+        vmax=ax_max(z.entries),
+        cmap="coolwarm",
+        **kwargs
+    )
+    plt.xlabel(xlabel)
+    plt.ylabel(ylabel)
+    plt.title(title)
+    plt.colorbar()
+
+    if not testing:  # pragma: no cover
+        plt.show()
+
+    return
diff --git a/pybamm/quick_plot.py b/pybamm/plotting/quick_plot.py
similarity index 98%
rename from pybamm/quick_plot.py
rename to pybamm/plotting/quick_plot.py
index 7f569840d..238cbfa77 100644
--- a/pybamm/quick_plot.py
+++ b/pybamm/plotting/quick_plot.py
@@ -44,24 +44,6 @@ def split_long_string(title, max_words=4):
         return first_line + "\n" + second_line
 
 
-def dynamic_plot(*args, **kwargs):
-    """
-    Creates a :class:`pybamm.QuickPlot` object (with arguments 'args' and keyword
-    arguments 'kwargs') and then calls :meth:`pybamm.QuickPlot.dynamic_plot`.
-    The key-word argument 'testing' is passed to the 'dynamic_plot' method, not the
-    `QuickPlot' class.
-
-    Returns
-    -------
-    plot : :class:`pybamm.QuickPlot`
-        The 'QuickPlot' object that was created
-    """
-    kwargs_for_class = {k: v for k, v in kwargs.items() if k != "testing"}
-    plot = pybamm.QuickPlot(*args, **kwargs_for_class)
-    plot.dynamic_plot(kwargs.get("testing", False))
-    return plot
-
-
 class QuickPlot(object):
     """
     Generates a quick plot of a subset of key outputs of the model so that the model


Test patch:
diff --git a/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_compare_outputs.py b/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_compare_outputs.py
index 440a8e97f..30eecc79a 100644
--- a/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_compare_outputs.py
+++ b/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_compare_outputs.py
@@ -99,7 +99,7 @@ class TestCompareOutputs(unittest.TestCase):
             solutions = []
             t_eval = np.linspace(0, 3600, 100)
             for model in models:
-                solution = pybamm.CasadiSolver(dt_max=0.01).solve(model, t_eval)
+                solution = pybamm.CasadiSolver().solve(model, t_eval)
                 solutions.append(solution)
 
             # compare outputs
diff --git a/tests/unit/test_expression_tree/test_array.py b/tests/unit/test_expression_tree/test_array.py
new file mode 100644
index 000000000..550abd84c
--- /dev/null
+++ b/tests/unit/test_expression_tree/test_array.py
@@ -0,0 +1,38 @@
+#
+# Tests for the Array class
+#
+import pybamm
+import numpy as np
+
+import unittest
+
+
+class TestArray(unittest.TestCase):
+    def test_name(self):
+        arr = pybamm.Array(np.array([1, 2, 3]))
+        self.assertEqual(arr.name, "Array of shape (3, 1)")
+
+    def test_linspace(self):
+        x = np.linspace(0, 1, 100)[:, np.newaxis]
+        y = pybamm.linspace(0, 1, 100)
+        np.testing.assert_array_equal(x, y.entries)
+
+    def test_meshgrid(self):
+        a = np.linspace(0, 5)
+        b = np.linspace(0, 3)
+        A, B = np.meshgrid(a, b)
+        c = pybamm.linspace(0, 5)
+        d = pybamm.linspace(0, 3)
+        C, D = pybamm.meshgrid(c, d)
+        np.testing.assert_array_equal(A, C.entries)
+        np.testing.assert_array_equal(B, D.entries)
+
+
+if __name__ == "__main__":
+    print("Add -v for more debug output")
+    import sys
+
+    if "-v" in sys.argv:
+        debug = True
+    pybamm.settings.debug_mode = True
+    unittest.main()
diff --git a/tests/unit/test_expression_tree/test_matrix.py b/tests/unit/test_expression_tree/test_matrix.py
index fa3aa26bb..0d2300ed5 100644
--- a/tests/unit/test_expression_tree/test_matrix.py
+++ b/tests/unit/test_expression_tree/test_matrix.py
@@ -32,12 +32,6 @@ class TestMatrix(unittest.TestCase):
         )
 
 
-class TestArray(unittest.TestCase):
-    def test_name(self):
-        arr = pybamm.Array(np.array([1, 2, 3]))
-        self.assertEqual(arr.name, "Array of shape (3, 1)")
-
-
 if __name__ == "__main__":
     print("Add -v for more debug output")
     import sys
diff --git a/tests/unit/test_plot.py b/tests/unit/test_plot.py
new file mode 100644
index 000000000..445095dc0
--- /dev/null
+++ b/tests/unit/test_plot.py
@@ -0,0 +1,47 @@
+import pybamm
+import unittest
+import numpy as np
+
+
+class TestPlot(unittest.TestCase):
+    def test_plot(self):
+        x = pybamm.Array(np.array([0, 3, 10]))
+        y = pybamm.Array(np.array([6, 16, 78]))
+        pybamm.plot(x, y, xlabel="x", ylabel="y", title="title", testing=True)
+
+    def test_plot_fail(self):
+        x = pybamm.Array(np.array([0]))
+        with self.assertRaisesRegex(TypeError, "x must be 'pybamm.Array'"):
+            pybamm.plot("bad", x)
+        with self.assertRaisesRegex(TypeError, "y must be 'pybamm.Array'"):
+            pybamm.plot(x, "bad")
+
+    def test_plot2D(self):
+        x = pybamm.Array(np.array([0, 3, 10]))
+        y = pybamm.Array(np.array([6, 16, 78]))
+        X, Y = pybamm.meshgrid(x, y)
+
+        # plot with array directly
+        pybamm.plot2D(x, y, Y, xlabel="x", ylabel="y", title="title", testing=True)
+
+        # plot with meshgrid
+        pybamm.plot2D(X, Y, Y, xlabel="x", ylabel="y", title="title", testing=True)
+
+    def test_plot2D_fail(self):
+        x = pybamm.Array(np.array([0]))
+        with self.assertRaisesRegex(TypeError, "x must be 'pybamm.Array'"):
+            pybamm.plot2D("bad", x, x)
+        with self.assertRaisesRegex(TypeError, "y must be 'pybamm.Array'"):
+            pybamm.plot2D(x, "bad", x)
+        with self.assertRaisesRegex(TypeError, "z must be 'pybamm.Array'"):
+            pybamm.plot2D(x, x, "bad")
+
+
+if __name__ == "__main__":
+    print("Add -v for more debug output")
+    import sys
+
+    if "-v" in sys.argv:
+        debug = True
+    pybamm.settings.debug_mode = True
+    unittest.main()