Raise when assigning attrs to virtual variables (default coordinate arrays)

Discussed in https://github.com/pydata/xarray/discussions/8567

^{Originally posted by matthew-brett December 21, 2023}
Sorry for the introductory question, but we (@ivanov and I) ran into this behavior while experimenting:

import numpy as np
data = np.zeros((3, 4, 5))
ds = xr.DataArray(data, dims=('i', 'j', 'k'))
print(ds['k'].attrs)

This shows {} as we might reasonably expect. But then:

ds['k'].attrs['foo'] = 'bar'
print(ds['k'].attrs)

This also gives {}, which we found surprising. We worked out why that was, after a little experimentation (the default coordinate arrays seems to get created on the fly and garbage collected immediately). But it took us a little while. Is that as intended? Is there a way of making this less confusing?

Thanks for any help.

The issue names no source file or test. Start by running the Python/NumPy reproducer with a default coordinate array, then trace how the virtual variable and its attrs are created and discarded. Done should be an explicit error when assigning attrs to such a variable, with a regression test covering the reported example.