Repository: SciTools/iris (language: python)

Issue:
Investigate segfaults with NetCDF4 `>1.6.0`
## 📰 Custom Issue
<!-- Provide a clear description of what the issue is, and we'll try our best to help 😀 -->

A bad interaction between Iris, Dask and NetCDF4 is causing segfaults when loading. This is apparently fixed if Dask is set to single-threaded mode - see https://github.com/Unidata/netcdf4-python/issues/1192#issuecomment-1267390240

This workaround removes a large part of the benefit of using Dask in Iris, so we need to find out if we can instead fix via changes to the Iris loading code.

Before the next release ([`v3.4`](https://github.com/orgs/SciTools/projects/3)), 1 of the below needs to happen:

* Permanent fix - change Iris' loading code to avoid the segfaults
* Temporary workaround - pin NetCDF4 `<=1.6.0`
* Temporary workaround, with big performance regressions - set Dask to single-threaded within Iris
Netcdf thread safety
Ultimately, to address #5016 

For now, **trial code for adding thread-safety and file-locking**.
Both in saver and loader.

This trial code is rather rough, and not very thought-through.

**At present**, all of the tests are locking up + not completing.
See notes below..

### learning so far

#### netcdf library is not considered thread-safe
Actually, it strictly never was, but it got worse with 1.6.1 which now releases the GIL during file operations
according to at least one [authorative comment elsewhere](https://github.com/Unidata/netcdf4-python/issues/1192#issuecomment-1255347911)

So in addition to any per-file locking that might be considered desirable (e.g. to control access to a file when doing lazy arithmetic between different netcdf variables), we **_do_** need a global lock to control all our access to the netcdf library, if any of it can be multi-threaded.  
... and with Dask, presumably, access **_can_** always be multi-threaded, since data chunks are processed by separate tasks --> different worker threads.

#### process-based scheduling
[like this](https://github.com/SciTools/iris/blob/4c430bb732f9224963e23d40cd9dc0b693175ca1/lib/iris/tests/integration/test_netcdf.py#L45)
Rather to my surprise, lots of our existing code breaks if you simply switch to process-based Dask workers (as are used by default for Bag computations).
I had rather thought that this would "just work"  :tm:  But it's now clear that a simple "da.store" operation is not serialisable, unless you give it a serialisable lock.  Just as we already saw in the [lazy saving investigation](https://github.com/SciTools/iris/blob/ed654db19a32d596dd5996cb7ca7e7e9ba51cf5b/lib/iris/fileformats/netcdf/saver.py#L2594-L2603).

#### single-threaded
[Like this](https://github.com/SciTools/iris/blob/4c430bb732f9224963e23d40cd9dc0b693175ca1/lib/iris/tests/integration/test_netcdf.py#L44-L46)
Solves the above problems, and probably all others, but it's not a practical solution.
Has already proved useful in tracking down where code uses locks wrong ,I.E. failing to release before re-acquiring --> simply hangs when single-threaded.

#### open/closing of netcdf Dataset objects
This is somewhat simpler than I had feared.
For the reading process, a netCDF4.Dataset is opened by creating a CFReader, from which cubes are extracted, and is **_held_** open :  in the new code, that is [_along with_ the netcdf library global access lock](https://github.com/SciTools/iris/blob/4c430bb732f9224963e23d40cd9dc0b693175ca1/lib/iris/fileformats/cf.py#L1075-L1092) 
During the process, the load code and 'actions' refer to to the dataset and its Variables, but once all the cubes have been produced, the CFReader is discarded : Then the dataset is closed, and the only remaining references are in the [NetcdfDataProxy](https://github.com/SciTools/iris/blob/5c9d9b30015003dc76f65c90db8e780c80a33062/lib/iris/fileformats/netcdf/loader.py#L58-L82) objects used to build lazy arrays -- which contain a filepath + variable-name, and thus don't refer to the 'original' Dataset object.

The only tricky part was that the existing code is relying on [`CFReader.__del__` to close the dataset](https://github.com/SciTools/iris/blob/5c9d9b30015003dc76f65c90db8e780c80a33062/lib/iris/fileformats/cf.py#L1297-L1300).
Which is in fact unreliable, and may not always happen when you expected/wanted it.  Effectively it depends on how the garbage collector works, which is not a constant.  I also found that this is not a technical detail, but really happens.
<details>
<summary> really simple example.. click to expand </summary>

```
    >>> class MyC:
    ...   def __init__(self, x):
    ...     self.x = x
    ...   def __str__(self):
    ...     return f"<MyC id={id(self)}(x={self.x})"
    ...   def __del__(self):
    ...     print('delete: ', self)
    ...
    >>> myc1 = MyC(1)
    >>> myc2 = MyC(2)
    >>> myc2 = None
    delete:  <MyC id=140005502283392(x=2)
    >>> myc1
    <__main__.MyC object at 0x7f55923aef70>
    >>> myc1 = 13
    #
    # ??? why no delete message here ???
    #
    >>> del myc1
    >>> exit()
    delete:  <MyC id=140005502283632(x=1)

```

</details>

So, this is actually unnecessary, and I reckon we would be better off using the CFReader [as a context manager, so it can automatically (always) close after operation](https://github.com/SciTools/iris/blob/4c430bb732f9224963e23d40cd9dc0b693175ca1/lib/iris/fileformats/cf.py#L1094-L1106).  So I am now doing that [here](https://github.com/SciTools/iris/blob/4c430bb732f9224963e23d40cd9dc0b693175ca1/lib/iris/fileformats/netcdf/loader.py#L556).
BTW this is exactly how the [netcdf.saver.Saver is already used](https://github.com/SciTools/iris/blob/5c9d9b30015003dc76f65c90db8e780c80a33062/lib/iris/fileformats/netcdf/saver.py#L2751).

#### nature of errors.
So far, various operations which succeed nevertheless produce a stream of HDF error messages.
Not clear how that can be (if it suceeds), and being multi-threased, naturally it's hard to debug.
Maybe a fault in my code design, or for some reason something else needs guarding 

#### lazy streaming netcdf-to-netcdf in save operations
There were some very obvious things going wrong with save-load-save roundtrips, like [this one](https://github.com/SciTools/iris/blob/5c9d9b30015003dc76f65c90db8e780c80a33062/lib/iris/tests/integration/test_netcdf.py#L98-L111).

This highlights that when streaming _from_ lazy netcdf data to save results, it is necessary that the saving process releases the global netcdf-access lock, since otherwise the lazy input compute will lock up.  
Which means both the lazy-write and the lazy-read operations must independently ensure that the lock is used when reading/writing.
N.B. with synthetic data, like in these tests, you _only_ see this if you save-load-save : The problem is unique to cases where the saved source data is lazy and in a netcdf file :  It's a particularly specific problem.

> it is necessary that the saving process releases the global netcdf-access lock

Which makes changes in the Saver code much more tricky + extensive than those in CFReader
The load_cubes/CFReader action is quite simple, since it involves no "compute" at all, only straight-line code to snapshot the file + create data-proxies.  At the end of which, the file is simply closed.

### Scope
Not considering process-parallel execution yet (dask config "scheduler='processes'", or distributed), which still does not work.
However, I have provided [per-file locking on read](https://github.com/SciTools/iris/blob/4c430bb732f9224963e23d40cd9dc0b693175ca1/lib/iris/fileformats/cf.py#L1067-L1070), and an eventual solution **_may_** need to adopt many of the ideas from #5130 



Discussion / hints:
@bjlittle has shared how xarray handles NetCDF4 not being thread-safe:

https://github.com/pydata/xarray/blob/main/xarray/backends/file_manager.py#L52
Ping @valeriupredoi 👍 
awesome, cheers @bjlittle - also, let me know if I can help :beer: 
# Current thinking

Fix will be a lot of work. We should get started as soon as we can. But we don't need to coincide with the planned release schedule.

## Rationale

### `single-threaded`? ❌

Have confirmed the suggestion that the `single-threaded` Dask scheduler avoids the problem. However much of our user-base would not cope with the slowdowns this would bring, especially if they are depending on large compute clusters to get fast enough execution times.

### File locking? ⏳

File locking (like Xarray) is not going to be a quick fix. Iris' load and save code uses file paths all over the place, which makes it difficult to protect access to a file. @bjlittle and I made a few unsuccessful attempts.

We probably need the common currency in Iris' internals to be a specialised file-handling object, instead of file paths. This new object could control file locking and also pickling, which is important for Dask. This will be a large piece of work, and difficult to estimate since it touches so much of Iris.

### Pinning? 😕

Pinning is awkward. It makes environment resolution less flexible, and can cause all sorts of other problems (see conversation on conda-forge/iris-feedstock#86 for some examples). **So we want to pin for the minimum amount of time**, but missing `v3.4` doesn't mean waiting for `v3.5` (disaster) - the fix will probably be too big for a bugfix but we can make an impromptu minor release whenever is appropriate.
@valeriupredoi will Iris `v3.4` still relying on a pin be a big problem for ESMValTool?
@trexfeathers this was recently mentioned https://github.com/ESMValGroup/ESMValCore/issues/1776#issuecomment-1305549061
Note : IIUC experiment is showing that the problem still occurs with a process-based scheduler, 
so it quite probably **isn't a true "thread-safety" issue**, but some more general problem with parallelisation ?

See also #5031 (lazy saving trial code) : We see that the mutex problems there are solved by xarray with file-specific locks, which is probably what we need here too.  
> Fix will be a lot of work. We should get started as soon as we can. 

**IMHO** code changes for this should probably be well localised -- just the netcdf load+save modules.  So, not a huge rewrite (and not much to document for users, either).  
**But** it isn't so clear whether it will be simple to get to a robust solution -- there could be a lot of research and experimentation hiding in here.

Aside : note also the approach of Xarray, as mentioned above :  Its locking solution also has a (netcdf) file **_caching_** scheme to improve efficiency :  It's probably not a requirement, but it would be interesting to see if this approach would be worth copying too.  After all, the current NetcdfDataProxy objects are opening+closing the file for [_every chunk of data read_](https://github.com/SciTools/iris/blob/421e193d9a73d767c64d989550c816cebdebecab/lib/iris/fileformats/netcdf/loader.py#L74-L82), which at some point has to be pretty inefficient (!)
### On the urgency,
@valeriupredoi has recently said ["pinning the package is not a long term solution"](https://github.com/ESMValGroup/ESMValCore/issues/1727#issuecomment-1285709710)

> @valeriupredoi will Iris `v3.4` still relying on a pin be a big problem for ESMValTool?

sorry I forgot to get meself here, guys! No, not a problem at all - as @zklaus pointed out, moving to Python=3.11 gonna be a bit lengthy anyway since we'd have to wait for all our dependency ducks to get in line with the new Python. But in the long run (ie a few months) it'd be good to have it. Let me know if I can help you in any way BTW :beer: 
Tests have been running for 30mins and counting...
> Tests have been running for 30mins and counting...

It is still deadlocking in some operations, to be investigated + fixed.
That's known about but I [removed the offending one in integration/test_netcdf.py](https://github.com/SciTools/iris/blob/4c430bb732f9224963e23d40cd9dc0b693175ca1/lib/iris/tests/integration/test_netcdf.py#L170-L184).
Other than that, at least _those_ tests now succeed -- but with the unexplained warning messages, as previously described.

I told you it was WIP !

### Notes on existing errors
I have been focussing on the set of tests in `iris.tests.integration.test_netcdf.py` as that has proved a useful testbed.
(where errors emerged in the original CI problems).

 * I've [skipped out](https://github.com/SciTools/iris/blob/4c430bb732f9224963e23d40cd9dc0b693175ca1/lib/iris/tests/integration/test_netcdf.py#L40-L41) all the tests there which now "just work".

 * The remaining (unskipped) tests are producing "interesting" behaviour -- in those cases the test succeeds but a lot of HDF errors appear  Like this...
```
HDF5-DIAG: Error detected in HDF5 (1.12.2) thread 1:
  #000: H5A.c line 528 in H5Aopen_by_name(): can't open attribute
    major: Attribute
    minor: Can't open object
  #001: H5VLcallback.c line 1091 in H5VL_attr_open(): attribute open failed
    major: Virtual Object Layer
    minor: Can't open object
  #002: H5VLcallback.c line 1058 in H5VL__attr_open(): attribute open failed
    major: Virtual Object Layer
    minor: Can't open object
 . . .
```
N.B. although these are _described_ as errors,  something nevertheless "appears to work"? Could be failing to read data, but not checking actual values ?

 * [This test](https://github.com/SciTools/iris/blob/4c430bb732f9224963e23d40cd9dc0b693175ca1/lib/iris/tests/integration/test_netcdf.py#L170-L184) was commented out because this one (still) actually **_hangs_**  -- a lot of roundtrip tests used to do that, before the locks were added, but this shows that something there has still been missed.

#### Update Mon 2022-11-14
Discussed with @bjlittle .  
Particularly [this test](https://github.com/SciTools/iris/blob/4c430bb732f9224963e23d40cd9dc0b693175ca1/lib/iris/tests/integration/test_netcdf.py#L518-L529), which triggers a stream of HDF errors even though it is "only" loading cubes
(i.e. not fetching or streaming lazy data).
Confirmed that data is fetched from NetcdfDataProxy-s, even when "just" loading like this.
@bjlittle points out that the post-load merge operation is a likely culprit.
@pp-mo also noted that the CFVariable holds "back-door" accesses to netCDF4 Variable objects, potentially even after the file is closed.  Which should be addressed (but might not be "the" problem here).  Then tried constructing a list of cubes, then yielding after the CFReader context exits, but that did _not_ give any more errors as we might have expected from this?
We're currently hoping this can be a bugfix for `v3.3` and `v3.4`.

Fix patch:
diff --git a/requirements/ci/py310.yml b/requirements/ci/py310.yml
index 87a16ba18a..6815c7fe6d 100644
--- a/requirements/ci/py310.yml
+++ b/requirements/ci/py310.yml
@@ -16,7 +16,7 @@ dependencies:
   - cftime >=1.5
   - dask-core >=2.26
   - matplotlib >=3.5
-  - netcdf4 !=1.6.1
+  - netcdf4 <1.6.1
   - numpy >=1.19
   - python-xxhash
   - pyproj
diff --git a/requirements/ci/py38.yml b/requirements/ci/py38.yml
index e2d696d30c..316e0868ac 100644
--- a/requirements/ci/py38.yml
+++ b/requirements/ci/py38.yml
@@ -16,7 +16,7 @@ dependencies:
   - cftime >=1.5
   - dask-core >=2.26
   - matplotlib >=3.5
-  - netcdf4 !=1.6.1
+  - netcdf4 <1.6.1
   - numpy >=1.19
   - python-xxhash
   - pyproj
diff --git a/requirements/ci/py39.yml b/requirements/ci/py39.yml
index 500ec2f80f..66e22c230f 100644
--- a/requirements/ci/py39.yml
+++ b/requirements/ci/py39.yml
@@ -16,7 +16,7 @@ dependencies:
   - cftime >=1.5
   - dask-core >=2.26
   - matplotlib >=3.5
-  - netcdf4 !=1.6.1
+  - netcdf4 <1.6.1
   - numpy >=1.19
   - python-xxhash
   - pyproj
diff --git a/setup.cfg b/setup.cfg
index ca35a8eb4e..f6276cb173 100644
--- a/setup.cfg
+++ b/setup.cfg
@@ -52,7 +52,7 @@ install_requires =
     cftime>=1.5.0
     dask[array]>=2.26
     matplotlib>=3.5
-    netcdf4!=1.6.1
+    netcdf4<1.6.1
     numpy>=1.19
     scipy
     shapely!=1.8.3


Test patch:
diff --git a/docs/src/whatsnew/latest.rst b/docs/src/whatsnew/latest.rst
index 57540b9ffe..657b196a94 100644
--- a/docs/src/whatsnew/latest.rst
+++ b/docs/src/whatsnew/latest.rst
@@ -121,7 +121,7 @@ This document explains the changes made to Iris for this release
 #. `@bjlittle`_ and `@trexfeathers`_ (reviewer) fixed an issue which prevented
    uncompressed PP fields with additional trailing padded words in the field
    data to be loaded and saved. (:pull:`5058`)
-   
+
 #. `@lbdreyer`_ and `@trexfeathers`_ (reviewer) fixed the handling of data when
    regridding with :class:`~iris.analysis.UnstructuredNearest` or calling
    :func:`~iris.analysis.trajectory.interpolate` such that the data type and mask is
@@ -184,8 +184,8 @@ This document explains the changes made to Iris for this release
    details.
    (:pull:`4968`)
 
-#. `@trexfeathers`_ introduced the ``netcdf4!=1.6.1`` pin to avoid a problem
-   with segfaults. (:pull:`4968`)
+#. `@trexfeathers`_ introduced the ``netcdf4<1.6.1`` pin to avoid a problem
+   with segfaults. (:pull:`4968`, :pull:`5075`, :issue:`5016`)
 
 #. `@trexfeathers`_ updated the Matplotlib colormap registration in
    :mod:`iris.palette` in response to a deprecation warning. Using the new