Repository: LSSTDESC/NaMaster (language: python)

Issue:
NmtField lmax or lmax_mask negtive don't work
This code fails:
```
f = nmt.NmtField(m, None, spin=0, lmax_mask=-1)
```
(also for `lmax`), with error. If it is intentional, the documentation needs to be updated.  

```
File ~/mambaforge/envs/tjpcov/lib/python3.12/site-packages/pymaster/field.py:185, in NmtField.__init__(self, mask, maps, spin, templates, beam, purify_e, purify_b, n_iter, n_iter_mask, tol_pinv, wcs, lmax, lmax_mask, masked_on_input, lite, mask_22, mask_12)
    183 if lmax_mask is None:
    184     lmax_mask = self.minfo.get_lmax()
--> 185 self.ainfo = ut.NmtAlmInfo(lmax)
    186 self.ainfo_mask = ut.NmtAlmInfo(lmax_mask)
    188 # Beam

File ~/mambaforge/envs/tjpcov/lib/python3.12/site-packages/pymaster/utils.py:452, in NmtAlmInfo.__init__(self, lmax)
    450 m = np.arange(self.mmax+1)
    451 self.mstart = (m*(2*self.lmax+1-m)//2).astype(np.uint64, copy=False)
--> 452 self.nelem = int(np.max(self.mstart) + (self.lmax+1))

File ~/mambaforge/envs/tjpcov/lib/python3.12/site-packages/numpy/core/fromnumeric.py:2810, in max(a, axis, out, keepdims, initial, where)
   2692 @array_function_dispatch(_max_dispatcher)
   2693 @set_module('numpy')
   2694 def max(a, axis=None, out=None, keepdims=np._NoValue, initial=np._NoValue,
   2695          where=np._NoValue):
   2696     """
   2697     Return the maximum of an array or maximum along an axis.
   2698 
   (...)
   2808     5
   2809     """
-> 2810     return _wrapreduction(a, np.maximum, 'max', axis, None, out,
   2811                           keepdims=keepdims, initial=initial, where=where)

File ~/mambaforge/envs/tjpcov/lib/python3.12/site-packages/numpy/core/fromnumeric.py:88, in _wrapreduction(obj, ufunc, method, axis, dtype, out, **kwargs)
     85         else:
     86             return reduction(axis=axis, out=out, **passkwargs)
---> 88 return ufunc.reduce(obj, axis, dtype, out, **passkwargs)

ValueError: zero-size array to reduction operation maximum which has no identity
```


Discussion / hints:
It might be related with this error. I believe that the check  `lmax <= 0` is not working at all. It seems that for `lmax=0`it actually computes it up to `l=0`

I got this error when I set `lmax = 0` in the field and reuse a precomputed covariance workspace
```
    def comp_gaussian_covariance(cw, spin_a1, spin_a2, spin_b1, spin_b2, wa, wb, ncl11, ncl12, ncl21, ncl22, dout):
>       return _nmtlib.comp_gaussian_covariance(cw, spin_a1, spin_a2, spin_b1, spin_b2, wa, wb, ncl11, ncl12, ncl21, ncl22, dout)
E       RuntimeError: Coupling coefficients only computed up to l=0, but you requirelmax=95. Recompute this workspace with a larger lmax
```
@carlosggarcia I think this is a problem with the documentation. If you want to use the default lmaxs, you should set them to `None`, which is the default. Does that make sense? I'll open a PR to fix the docs.

Fix patch:
diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml
index c22857a..dd9feb5 100644
--- a/.github/workflows/ci.yml
+++ b/.github/workflows/ci.yml
@@ -48,7 +48,7 @@ jobs:
           channel-priority: strict
           show-channel-urls: true
           miniforge-version: latest
-          miniforge-variant: Mambaforge
+          miniforge-variant: Miniforge3
 
       - name: Set current date as environment variable
         run: echo "TODAY=$(date +'%Y.%m.%d')" >> $GITHUB_ENV
diff --git a/pymaster/field.py b/pymaster/field.py
index 9bc6398..2028122 100644
--- a/pymaster/field.py
+++ b/pymaster/field.py
@@ -70,11 +70,11 @@ class NmtField(object):
             :meth:`~pymaster.utils.get_default_params`,
             and modified via :meth:`~pymaster.utils.set_n_iter_default`.
         lmax (:obj:`int`): Maximum multipole up to which map power spectra
-            will be computed. If negative or zero, the maximum multipole given
+            will be computed. If ``None``, the maximum multipole given
             the map resolution will be used (e.g. :math:`3N_{\\rm side}-1`
             for HEALPix maps).
         lmax_mask (:obj:`int`): Maximum multipole up to which the power
-            spectrum of the mask will be computed. If negative or zero, the
+            spectrum of the mask will be computed. If ``None``, the
             maximum multipole given the map resolution will be used (e.g.
             :math:`3N_{\\rm side}-1` for HEALPix maps).
         tol_pinv (:obj:`float`): When computing the pseudo-inverse of the
@@ -182,6 +182,9 @@ def __init__(self, mask, maps, *, spin=None, templates=None, beam=None,
             lmax = self.minfo.get_lmax()
         if lmax_mask is None:
             lmax_mask = self.minfo.get_lmax()
+        if (lmax <= 0) or (lmax_mask <= 0):
+            raise ValueError("`lmax` and `lmax_mask` must be positive.")
+
         self.ainfo = ut.NmtAlmInfo(lmax)
         self.ainfo_mask = ut.NmtAlmInfo(lmax_mask)
 
diff --git a/pyproject.toml b/pyproject.toml
index 347efae..0677000 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -7,7 +7,7 @@ build-backend = "setuptools.build_meta"
 
 [project]
 name = "pymaster"
-version = "2.3.2"
+version = "2.3.3"
 authors = [
   {name="David Alonso", email="david.alonso@physics.ox.ac.uk"}
 ]
diff --git a/src/nmt_covar.c b/src/nmt_covar.c
index 0fcecae..75a74ee 100644
--- a/src/nmt_covar.c
+++ b/src/nmt_covar.c
@@ -96,7 +96,7 @@ void  nmt_compute_gaussian_covariance_coupled(nmt_covar_workspace *cw,
                                               flouble *covar_out)
 {
   if((cw->lmax<wa->bin->ell_max) || (cw->lmax<wb->bin->ell_max))
-    report_error(NMT_ERROR_COVAR,"Coupling coefficients only computed up to l=%d, but you require"
+    report_error(NMT_ERROR_COVAR,"Coupling coefficients only computed up to l=%d, but you require "
 		 "lmax=%d. Recompute this workspace with a larger lmax\n",cw->lmax,wa->bin->ell_max);
   if((cw->spin0_only) && (spin_a+spin_b+spin_c+spin_d != 0))
     report_error(NMT_ERROR_COVAR,"Coupling coefficients only computed for spin=0\n");
@@ -191,7 +191,7 @@ void  nmt_compute_gaussian_covariance(nmt_covar_workspace *cw,
 				      flouble *covar_out)
 {
   if((cw->lmax<wa->bin->ell_max) || (cw->lmax<wb->bin->ell_max))
-    report_error(NMT_ERROR_COVAR,"Coupling coefficients only computed up to l=%d, but you require"
+    report_error(NMT_ERROR_COVAR,"Coupling coefficients only computed up to l=%d, but you require "
 		 "lmax=%d. Recompute this workspace with a larger lmax\n",cw->lmax,wa->bin->ell_max);
   if((cw->spin0_only) && (spin_a+spin_b+spin_c+spin_d != 0))
     report_error(NMT_ERROR_COVAR,"Coupling coefficients only computed for spin=0\n");


Test patch:
diff --git a/pymaster/tests/test_field.py b/pymaster/tests/test_field.py
index 9659217..f8bc12e 100644
--- a/pymaster/tests/test_field.py
+++ b/pymaster/tests/test_field.py
@@ -308,3 +308,9 @@ def test_field_error():
     with pytest.raises(ValueError):
         f = nmt.NmtField(FT.msk, [FT.mps[1], FT.mps[2]], spin=1,
                          purify_b=True, n_iter=0)
+
+    # lmax must be zero
+    with pytest.raises(ValueError):
+        nmt.NmtField(FT.msk, [FT.mps[0]], lmax=0)
+    with pytest.raises(ValueError):
+        nmt.NmtField(FT.msk, [FT.mps[0]], lmax_mask=0)