code: def expand_grid_coords(dataset, dim):
      
      arrays = [dataset.interface.coords(dataset, d.name, True)
                for d in dataset.kdims]
      idx = dataset.get_dimension_index(dim)
      return cartesian_product(arrays, flat=False)[idx].T
docs: Expand the coordinates along a dimension of the gridded
      dataset into an ND-array matching the dimensionality of
      the dataset.