Repository: dask/dask

Issue:
Task dependencies do not match chunking when using `from_array`
If I create an array chunked along the first axis like so:
```python
shape = (2,2)
chunks = (1,-1)
x = da.random.random(size=shape, chunks=chunks)
x.visualize()
```
![image](https://user-images.githubusercontent.com/3805136/64438270-d3956a80-d095-11e9-83cd-c700a42ad965.png)
The task graph reflects that each chunk is independent. So far so good. 

But if I use the same data and chunking to construct a dask array using `from_array`, the topology of the resulting task graph is different:
``` python
x_ = x.compute()
y = da.from_array(x_, chunks=chunks, name='test')
y.visualize()
```
![image](https://user-images.githubusercontent.com/3805136/64438546-4dc5ef00-d096-11e9-9db2-6b5660419d29.png)

Now both leaves of the task graph share a dependency on the original array _which is independent of the chunking scheme_. 

I feel like the task graph created by `da.random.random` is "right" and the task graph created by `from_array` is "wrong" with respect to my chunking scheme. However, `from_array` seems natural for users, and this can lead to dead workers when using the distributed scheduler, as I describe in another issue (See https://github.com/dask/distributed/issues/3032).

How to fix this? Maybe just make `from_array` call `stack`, `concatenate`, etc under the hood? This guarantees that the outputs of these array creation routines are more uniform. I'd be willing to look into a PR along these lines, but I want to be sure I'm on the right track.

Discussion:
mrocklin: In principle I agree that supporting a special case here makes sense.  

I wonder if we might run into trouble with memmapped arrays though?  Are those a different type and so won't be caught due to the strict type check, or do we need to be more careful here?
dcherian: > Are those a different type and so won't be caught due to the strict type check, or do we need to be more careful here?

Yes. I just checked. Memmaps have type `numpy.memmap` so they fail the strict type check. The behaviour for memmaps is unchanged.
TomAugspurger: This seems worth trying out, and I think basing it off array.chunk-size makes sense.

Could you add a short section to array-creation.rst explaining this behavior?
mrocklin: > Could you add a short section to array-creation.rst explaining this behavior?

This might also be overly nuanced for the main docs.  I don't expect that most users will want to go this deeply.  Maybe an entry in the docstring instead?



Patch:
diff --git a/dask/array/core.py b/dask/array/core.py
index 404e5cecf..767785073 100644
--- a/dask/array/core.py
+++ b/dask/array/core.py
@@ -2831,6 +2831,14 @@ def from_array(
     >>> import dask.base
     >>> token = dask.base.tokenize(x)  # doctest: +SKIP
     >>> a = da.from_array('myarray-' + token)  # doctest: +SKIP
+
+    Numpy ndarrays are eagerly sliced and then embedded in the graph.
+
+    >>> import dask.array
+    >>> a = dask.array.from_array(np.array([[1, 2], [3, 4]]), chunks=(1,1))
+    >>> a.dask[a.name, 0, 0][0]
+    array([1])
+
     """
     if isinstance(x, Array):
         raise ValueError(
@@ -2866,17 +2874,24 @@ def from_array(
     if lock is True:
         lock = SerializableLock()
 
+    is_ndarray = type(x) is np.ndarray
+    is_single_block = all(len(c) == 1 for c in chunks)
     # Always use the getter for h5py etc. Not using isinstance(x, np.ndarray)
     # because np.matrix is a subclass of np.ndarray.
-    if type(x) is np.ndarray and all(len(c) == 1 for c in chunks):
+    if is_ndarray and not is_single_block and not lock:
+        # eagerly slice numpy arrays to prevent memory blowup
+        # GH5367, GH5601
+        slices = slices_from_chunks(chunks)
+        keys = product([name], *(range(len(bds)) for bds in chunks))
+        values = [x[slc] for slc in slices]
+        dsk = dict(zip(keys, values))
+
+    elif is_ndarray and is_single_block:
         # No slicing needed
         dsk = {(name,) + (0,) * x.ndim: x}
     else:
         if getitem is None:
-            if type(x) is np.ndarray and not lock:
-                # simpler and cleaner, but missing all the nuances of getter
-                getitem = operator.getitem
-            elif fancy:
+            if fancy:
                 getitem = getter
             else:
                 getitem = getter_nofancy
diff --git a/dask/array/creation.py b/dask/array/creation.py
index 8d4e1f448..4e890122f 100644
--- a/dask/array/creation.py
+++ b/dask/array/creation.py
@@ -691,6 +691,7 @@ def triu(m, k=0):
     token = tokenize(m, k)
     name = "triu-" + token
 
+    triu_is_empty = True
     dsk = {}
     for i in range(rdim):
         for j in range(hdim):
@@ -701,9 +702,13 @@ def triu(m, k=0):
                 )
             elif chunk * (j - i - 1) < k <= chunk * (j - i + 1):
                 dsk[(name, i, j)] = (np.triu, (m.name, i, j), k - (chunk * (j - i)))
+                triu_is_empty = False
             else:
                 dsk[(name, i, j)] = (m.name, i, j)
-    graph = HighLevelGraph.from_collections(name, dsk, dependencies=[m])
+                triu_is_empty = False
+    graph = HighLevelGraph.from_collections(
+        name, dsk, dependencies=[] if triu_is_empty else [m]
+    )
     return Array(graph, name, shape=m.shape, chunks=m.chunks, meta=m)
 
 
@@ -744,19 +749,24 @@ def tril(m, k=0):
     token = tokenize(m, k)
     name = "tril-" + token
 
+    tril_is_empty = True
     dsk = {}
     for i in range(rdim):
         for j in range(hdim):
             if chunk * (j - i + 1) < k:
                 dsk[(name, i, j)] = (m.name, i, j)
+                tril_is_empty = False
             elif chunk * (j - i - 1) < k <= chunk * (j - i + 1):
                 dsk[(name, i, j)] = (np.tril, (m.name, i, j), k - (chunk * (j - i)))
+                tril_is_empty = False
             else:
                 dsk[(name, i, j)] = (
                     partial(zeros_like_safe, shape=(m.chunks[0][i], m.chunks[1][j])),
                     m._meta,
                 )
-    graph = HighLevelGraph.from_collections(name, dsk, dependencies=[m])
+    graph = HighLevelGraph.from_collections(
+        name, dsk, dependencies=[] if tril_is_empty else [m]
+    )
     return Array(graph, name, shape=m.shape, chunks=m.chunks, meta=m)
 
 
diff --git a/dask/array/linalg.py b/dask/array/linalg.py
index 5e69999f0..06ff71453 100644
--- a/dask/array/linalg.py
+++ b/dask/array/linalg.py
@@ -139,7 +139,7 @@ def tsqr(data, compute_svd=False, _max_vchunk_size=None):
         numblocks={data.name: numblocks},
     )
     layers[name_qr_st1] = dsk_qr_st1
-    dependencies[name_qr_st1] = data.__dask_layers__()
+    dependencies[name_qr_st1] = set(data.__dask_layers__())
 
     # Block qr[0]
     name_q_st1 = "getitem" + token + "-q1"
@@ -248,13 +248,13 @@ def tsqr(data, compute_svd=False, _max_vchunk_size=None):
             )
         )
         layers[name_q_st2] = dsk_q_st2
-        dependencies[name_q_st2] = q_inner.__dask_layers__()
+        dependencies[name_q_st2] = set(q_inner.__dask_layers__())
 
         # R: R_inner
         name_r_st2 = "r-inner" + token
         dsk_r_st2 = {(name_r_st2, 0, 0): (r_inner.name, 0, 0)}
         layers[name_r_st2] = dsk_r_st2
-        dependencies[name_r_st2] = r_inner.__dask_layers__()
+        dependencies[name_r_st2] = set(r_inner.__dask_layers__())
 
         # Q: Block qr[0] (*) Q_inner
         name_q_st3 = "dot" + token + "-q3"
@@ -301,7 +301,8 @@ def tsqr(data, compute_svd=False, _max_vchunk_size=None):
         layers[name_q_st2_aux] = dsk_q_st2_aux
         dependencies[name_q_st2_aux] = {name_qr_st2}
 
-        if not any(np.isnan(c) for cs in data.chunks for c in cs):
+        chucks_are_all_known = not any(np.isnan(c) for cs in data.chunks for c in cs)
+        if chucks_are_all_known:
             # when chunks are all known...
             # obtain slices on q from in-core compute (e.g.: (slice(10, 20), slice(0, 5)))
             q2_block_sizes = [min(e, n) for e in data.chunks[0]]
@@ -349,7 +350,7 @@ def tsqr(data, compute_svd=False, _max_vchunk_size=None):
                 dsk_n, dsk_q2_shapes, dsk_q2_cumsum, dsk_block_slices
             )
 
-            deps = {data.name, name_q2bs, name_q2cs}
+            deps = {data.name}
             block_slices = [(name_blockslice, i) for i in range(numblocks[0])]
 
         layers["q-blocksizes" + token] = dsk_q_blockslices
@@ -362,7 +363,10 @@ def tsqr(data, compute_svd=False, _max_vchunk_size=None):
             for i, b in enumerate(block_slices)
         )
         layers[name_q_st2] = dsk_q_st2
-        dependencies[name_q_st2] = {name_q_st2_aux, "q-blocksizes" + token}
+        if chucks_are_all_known:
+            dependencies[name_q_st2] = {name_q_st2_aux}
+        else:
+            dependencies[name_q_st2] = {name_q_st2_aux, "q-blocksizes" + token}
 
         # Q: Block qr[0] (*) In-core qr[0]
         name_q_st3 = "dot" + token + "-q3"
@@ -565,11 +569,14 @@ def sfqr(data, name=None):
     name_A_1 = prefix + "A_1"
     name_A_rest = prefix + "A_rest"
     layers[name_A_1] = {(name_A_1, 0, 0): (data.name, 0, 0)}
-    dependencies[name_A_1] = data.__dask_layers__()
+    dependencies[name_A_1] = set(data.__dask_layers__())
     layers[name_A_rest] = {
         (name_A_rest, 0, idx): (data.name, 0, 1 + idx) for idx in range(nc - 1)
     }
-    dependencies[name_A_rest] = data.__dask_layers__()
+    if len(layers[name_A_rest]) > 0:
+        dependencies[name_A_rest] = set(data.__dask_layers__())
+    else:
+        dependencies[name_A_rest] = set()
 
     # Q R_1 = A_1
     name_Q_R1 = prefix + "Q_R_1"
@@ -1300,7 +1307,7 @@ def lstsq(a, b):
             (np.sqrt, (np.linalg.eigvals, (np.dot, (rt.name, 0, 0), (r.name, 0, 0)))),
         )
     }
-    graph = HighLevelGraph.from_collections(sname, sdsk, dependencies=[rt])
+    graph = HighLevelGraph.from_collections(sname, sdsk, dependencies=[rt, r])
     _, _, _, ss = np.linalg.lstsq(
         np.array([[1, 0], [1, 2]], dtype=a.dtype),
         np.array([0, 1], dtype=b.dtype),
diff --git a/dask/array/utils.py b/dask/array/utils.py
index 9b25ab634..7237fb46a 100644
--- a/dask/array/utils.py
+++ b/dask/array/utils.py
@@ -178,6 +178,7 @@ def _check_dsk(dsk):
     if not isinstance(dsk, HighLevelGraph):
         return
 
+    dsk.validate()
     assert all(isinstance(k, (tuple, str)) for k in dsk.layers)
     freqs = frequencies(concat(dsk.dicts.values()))
     non_one = {k: v for k, v in freqs.items() if v != 1}
diff --git a/dask/blockwise.py b/dask/blockwise.py
index 872b0d379..b634563d0 100644
--- a/dask/blockwise.py
+++ b/dask/blockwise.py
@@ -6,7 +6,7 @@ import numpy as np
 
 import tlz as toolz
 
-from .core import reverse_dict
+from .core import reverse_dict, keys_in_tasks
 from .delayed import unpack_collections
 from .highlevelgraph import HighLevelGraph
 from .optimization import SubgraphCallable, fuse
@@ -636,7 +636,14 @@ def _optimize_blockwise(full_graph, keys=()):
             # Merge these Blockwise layers into one
             new_layer = rewrite_blockwise([layers[l] for l in blockwise_layers])
             out[layer] = new_layer
-            dependencies[layer] = {k for k, v in new_layer.indices if v is not None}
+
+            new_deps = set()
+            for k, v in new_layer.indices:
+                if v is None:
+                    new_deps |= keys_in_tasks(full_graph.dependencies, [k])
+                else:
+                    new_deps.add(k)
+            dependencies[layer] = new_deps
         else:
             out[layer] = layers[layer]
             dependencies[layer] = full_graph.dependencies.get(layer, set())
diff --git a/dask/core.py b/dask/core.py
index 2d4f0f988..533bba86d 100644
--- a/dask/core.py
+++ b/dask/core.py
@@ -156,6 +156,41 @@ def get(dsk, out, cache=None):
     return result
 
 
+def keys_in_tasks(keys, tasks, as_list=False):
+    """Returns the keys in `keys` that are also in `tasks`
+
+    Examples
+    --------
+    >>> dsk = {'x': 1,
+    ...        'y': (inc, 'x'),
+    ...        'z': (add, 'x', 'y'),
+    ...        'w': (inc, 'z'),
+    ...        'a': (add, (inc, 'x'), 1)}
+
+    >>> keys_in_tasks(dsk, ['x', 'y', 'j'])  # doctest: +SKIP
+    {'x', 'y'}
+    """
+    ret = []
+    while tasks:
+        work = []
+        for w in tasks:
+            typ = type(w)
+            if typ is tuple and w and callable(w[0]):  # istask(w)
+                work.extend(w[1:])
+            elif typ is list:
+                work.extend(w)
+            elif typ is dict:
+                work.extend(w.values())
+            else:
+                try:
+                    if w in keys:
+                        ret.append(w)
+                except TypeError:  # not hashable
+                    pass
+        tasks = work
+    return ret if as_list else set(ret)
+
+
 def get_dependencies(dsk, key=None, task=no_default, as_list=False):
     """Get the immediate tasks on which this task depends
 
@@ -192,28 +227,7 @@ def get_dependencies(dsk, key=None, task=no_default, as_list=False):
     else:
         raise ValueError("Provide either key or task")
 
-    result = []
-    work = [arg]
-
-    while work:
-        new_work = []
-        for w in work:
-            typ = type(w)
-            if typ is tuple and w and callable(w[0]):  # istask(w)
-                new_work.extend(w[1:])
-            elif typ is list:
-                new_work.extend(w)
-            elif typ is dict:
-                new_work.extend(w.values())
-            else:
-                try:
-                    if w in dsk:
-                        result.append(w)
-                except TypeError:  # not hashable
-                    pass
-        work = new_work
-
-    return result if as_list else set(result)
+    return keys_in_tasks(dsk, [arg], as_list=as_list)
 
 
 def get_deps(dsk):
diff --git a/dask/dataframe/core.py b/dask/dataframe/core.py
index 60c3bab6b..b8a134c5f 100644
--- a/dask/dataframe/core.py
+++ b/dask/dataframe/core.py
@@ -4743,6 +4743,7 @@ def elemwise(op, *args, **kwargs):
     dfs = [df for df in dasks if isinstance(df, _Frame)]
 
     # Clean up dask arrays if present
+    deps = dasks.copy()
     for i, a in enumerate(dasks):
         if not isinstance(a, Array):
             continue
@@ -4785,7 +4786,7 @@ def elemwise(op, *args, **kwargs):
     # adjust the key length of Scalar
     dsk = partitionwise_graph(op, _name, *args, **kwargs)
 
-    graph = HighLevelGraph.from_collections(_name, dsk, dependencies=dasks)
+    graph = HighLevelGraph.from_collections(_name, dsk, dependencies=deps)
 
     if meta is no_default:
         if len(dfs) >= 2 and not all(hasattr(d, "npartitions") for d in dasks):
diff --git a/dask/dataframe/indexing.py b/dask/dataframe/indexing.py
index 9fd095d2a..0380af426 100644
--- a/dask/dataframe/indexing.py
+++ b/dask/dataframe/indexing.py
@@ -164,10 +164,11 @@ class _LocIndexer(_IndexerBase):
                 divisions.append(sorted(indexer)[0])
             # append maximum value of the last division
             divisions.append(sorted(items[-1][1])[-1])
+            graph = HighLevelGraph.from_collections(name, dsk, dependencies=[self.obj])
         else:
             divisions = [None, None]
             dsk = {(name, 0): meta.head(0)}
-        graph = HighLevelGraph.from_collections(name, dsk, dependencies=[self.obj])
+            graph = HighLevelGraph.from_collections(name, dsk)
         return new_dd_object(graph, name, meta=meta, divisions=divisions)
 
     def _loc_element(self, iindexer, cindexer):
diff --git a/dask/dataframe/utils.py b/dask/dataframe/utils.py
index 8285ad62a..eca5b2517 100644
--- a/dask/dataframe/utils.py
+++ b/dask/dataframe/utils.py
@@ -715,7 +715,10 @@ def index_summary(idx, name=None):
 def _check_dask(dsk, check_names=True, check_dtypes=True, result=None):
     import dask.dataframe as dd
 
-    if hasattr(dsk, "dask"):
+    if hasattr(dsk, "__dask_graph__"):
+        graph = dsk.__dask_graph__()
+        if hasattr(graph, "validate"):
+            graph.validate()
         if result is None:
             result = dsk.compute(scheduler="sync")
         if isinstance(dsk, dd.Index):
diff --git a/dask/delayed.py b/dask/delayed.py
index 5c9714961..a8bfbec15 100644
--- a/dask/delayed.py
+++ b/dask/delayed.py
@@ -470,13 +470,20 @@ class Delayed(DaskMethodsMixin, OperatorMethodMixin):
     Equivalent to the output from a single key in a dask graph.
     """
 
-    __slots__ = ("_key", "dask", "_length")
+    __slots__ = ("_key", "dask", "_length", "_dask_layers")
 
     def __init__(self, key, dsk, length=None):
         self._key = key
         self.dask = dsk
         self._length = length
 
+        # `__dask_layers__` of an one-layered HLG is the layer name
+        # TODO: what is the name of a multi-layered HLG?
+        if isinstance(dsk, HighLevelGraph) and len(dsk.layers) == 1:
+            self._dask_layers = set(iter(dsk.layers))
+        else:
+            self._dask_layers = (self.key,)
+
     def __dask_graph__(self):
         return self.dask
 
@@ -484,7 +491,10 @@ class Delayed(DaskMethodsMixin, OperatorMethodMixin):
         return [self.key]
 
     def __dask_layers__(self):
-        return (self.key,)
+        if hasattr(self, "_dask_layers"):
+            return self._dask_layers
+        else:
+            return (self.key,)
 
     def __dask_tokenize__(self):
         return self.key
diff --git a/dask/highlevelgraph.py b/dask/highlevelgraph.py
index 9fd70f765..a0f772699 100644
--- a/dask/highlevelgraph.py
+++ b/dask/highlevelgraph.py
@@ -4,7 +4,24 @@ import tlz as toolz
 
 from .utils import ignoring
 from .base import is_dask_collection
-from .core import reverse_dict
+from .core import reverse_dict, keys_in_tasks
+
+
+def compute_layer_dependencies(layers):
+    """Returns the dependencies between layers"""
+
+    def _find_layer_containing_key(key):
+        for k, v in layers.items():
+            if key in v:
+                return k
+        raise RuntimeError(f"{repr(key)} not found")
+
+    all_keys = set(key for layer in layers.values() for key in layer)
+    ret = {k: set() for k in layers.keys()}
+    for k, v in layers.items():
+        for key in keys_in_tasks(all_keys.difference(v.keys()), v.values()):
+            ret[k].add(_find_layer_containing_key(key))
+    return ret
 
 
 class HighLevelGraph(Mapping):
@@ -228,6 +245,26 @@ class HighLevelGraph(Mapping):
                 if dep not in self.dependencies:
                     raise ValueError(f"{repr(dep)} not found in dependencies")
 
+        # Re-calculate all layer dependencies
+        dependencies = compute_layer_dependencies(self.layers)
+
+        # Check keys
+        dep_key1 = set(self.dependencies.keys())
+        dep_key2 = set(dependencies.keys())
+        if dep_key1 != dep_key2:
+            raise ValueError(
+                f"incorrect dependencies keys {repr(dep_key1)} "
+                f"expected {repr(dep_key2)}"
+            )
+
+        # Check values
+        for k in dep_key1:
+            if self.dependencies[k] != dependencies[k]:
+                raise ValueError(
+                    f"incorrect dependencies[{repr(k)}]: {repr(self.dependencies[k])} "
+                    f"expected {repr(dependencies[k])}"
+                )
+
 
 def to_graphviz(
     hg,
diff --git a/docs/source/develop.rst b/docs/source/develop.rst
index 0af2963c2..4bee27e9a 100644
--- a/docs/source/develop.rst
+++ b/docs/source/develop.rst
@@ -124,7 +124,7 @@ language support, testing, documentation, and style.
 Python Versions
 ~~~~~~~~~~~~~~~
 
-Dask supports Python versions 3.5, 3.6, and 3.7.
+Dask supports Python versions 3.6, 3.7, and 3.8.
 Name changes are handled by the :file:`dask/compatibility.py` file.
 
 Test
@@ -288,7 +288,7 @@ Documentation is maintained in the RestructuredText markup language (``.rst``
 files) in ``dask/docs/source``.  The documentation consists both of prose
 and API documentation.
 
-To build the documentation locally, clone this repository and install 
+To build the documentation locally, clone this repository and install
 the necessary requirements using ``pip`` or ``conda``::
 
   git clone https://github.com/dask/dask.git
@@ -312,4 +312,4 @@ The resulting HTML files end up in the ``build/html`` directory.
 You can now make edits to rst files and run ``make html`` again to update
 the affected pages.
 
-.. _Sphinx: https://www.sphinx-doc.org/
\ No newline at end of file
+.. _Sphinx: https://www.sphinx-doc.org/


Test patch:
diff --git a/dask/array/tests/test_array_core.py b/dask/array/tests/test_array_core.py
index 9d7700e7b..54d496aa1 100644
--- a/dask/array/tests/test_array_core.py
+++ b/dask/array/tests/test_array_core.py
@@ -2387,7 +2387,7 @@ def test_from_array_ndarray_getitem():
     x = np.array([[1, 2], [3, 4]])
     dx = da.from_array(x, chunks=(1, 2))
     assert_eq(x, dx)
-    assert dx.dask[dx.name, 0, 0][0] == operator.getitem
+    assert (dx.dask[dx.name, 0, 0] == np.array([[1, 2]])).all()
 
 
 @pytest.mark.parametrize("x", [[1, 2], (1, 2), memoryview(b"abc")])
@@ -2399,8 +2399,7 @@ def test_from_array_list(x):
 
     dx = da.from_array(x, chunks=1)
     assert_eq(np.array(x), dx)
-    assert dx.dask[dx.name, 0][0] == operator.getitem
-    assert isinstance(dx.dask[dx.name.replace("array", "array-original")], np.ndarray)
+    assert dx.dask[dx.name, 0][0] == x[0]
 
 
 @pytest.mark.parametrize("type_", [t for t in np.ScalarType if t is not memoryview])
@@ -3295,7 +3294,7 @@ def test_from_array_names():
     d = da.from_array(x, chunks=2)
 
     names = countby(key_split, d.dask)
-    assert set(names.values()) == set([1, 5])
+    assert set(names.values()) == set([5])
 
 
 @pytest.mark.parametrize(
diff --git a/dask/array/tests/test_optimization.py b/dask/array/tests/test_optimization.py
index bef4722e4..aa632604d 100644
--- a/dask/array/tests/test_optimization.py
+++ b/dask/array/tests/test_optimization.py
@@ -264,7 +264,8 @@ def test_dont_fuse_numpy_arrays():
 
 
 def test_minimize_data_transfer():
-    x = np.ones(100)
+    zarr = pytest.importorskip("zarr")
+    x = zarr.ones((100,))
     y = da.from_array(x, chunks=25)
     z = y + 1
     dsk = z.__dask_optimize__(z.dask, z.__dask_keys__())