bug: export to netCDF file only filled with NaNs

[ThibHlln]

Describe the bug
I have a MODFLOW6 model that I load with flopy and I am trying to export the data contained in each package in netCDF files.

I am taking inspiration from this tutorial https://flopy.readthedocs.io/en/3.10.0/Notebooks/export_tutorial.html. But as I am dealing with a MODFLOW6 and not a MODFLOW2005 model, I am not using:
ml = flopy.modflow.Modflow.load(...) but sim = flopy.mf6.MFSimulation.load(...) + gwf = sim.get_model(sim.model_names[0]) to access the packages.

Then I try to use the .export() methods on the whole model (i.e. gwf.export('model.nc')) or on a specific package (i.e. gwf.dis.export('dis.nc')) but the variables the netCDF file contains seem only initialised and never populated with data (i.e. only containing NaNs) even though my package variables would contain actual data (no NaN).

To Reproduce
Steps to reproduce the behavior:

1. Load the MODFLOW model

import flopy
sim = flopy.mf6.MFSimulation.load(
    sim_name="my-model",
    version="mf6",
    exe_name="mf6.exe",
    sim_ws="mf6",
    verbosity_level=2,
    load_only=['dis'],
    strict=False
)
gwf = sim.get_model(sim.model_names[0])

2. Export the model in a netCDF file

gwf.export('model.nc', fmt='netcdf')

3. Check the topography values in the DIS package

print(gwf.dis.top.array)

[[ 92.98463  92.98469  92.99091 ...  94.98059  94.53332  94.08605]
 [ 92.99647  92.96922  92.9785  ...  95.09064  94.64922  94.21827]
 [ 93.0248   92.99885  93.00288 ...  95.23484  94.80389  94.37293]
 ...
 [102.2147  102.0641  101.9308  ... 101.8908  101.916   101.9429 ]
 [102.2224  102.1369  102.0565  ... 101.861   101.8861  101.9125 ]
 [102.3515  102.2665  102.1895  ... 101.8313  101.8563  101.882  ]]

4. Check the topography values in the netCDF file generated

import xarray as xr
ds = xr.load_dataset('model.nc')
print(ds)

<xarray.Dataset> Size: 17MB
Dimensions:             (layer: 24, y: 161, x: 167, time: 72)
Coordinates:
  * layer               (layer) float32 96B 0.0 1.0 2.0 3.0 ... 21.0 22.0 23.0
    longitude           (y, x) float64 215kB nan nan nan nan ... nan nan nan nan
    latitude            (y, x) float64 215kB nan nan nan nan ... nan nan nan nan
  * time                (time) datetime64[ns] 576B 1970-01-15 ... 1972-10-05
Dimensions without coordinates: y, x
Data variables:
    crs                 int32 4B -2147483647
    elevation           (layer, y, x) float64 5MB 92.63 92.63 ... 56.48 56.49
    x_proj              (y, x) float64 215kB 1.978e+05 1.978e+05 ... 1.991e+05
    y_proj              (y, x) float64 215kB 7.035e+06 7.035e+06 ... 7.034e+06
    latitude_longitude  (time, layer) float64 14kB nan nan nan ... nan nan nan
    delc                (y) float32 644B 5.977 5.977 5.977 ... 5.977 5.977 5.977
    delr                (x) float32 668B 5.987 5.987 5.987 ... 5.987 5.987 5.987
    VerticalTransform   |S1 1B b''
    top                 (y, x) float64 215kB nan nan nan nan ... nan nan nan nan
    botm                (layer, y, x) float64 5MB nan nan nan ... nan nan nan
    idomain             (layer, y, x) float64 5MB nan nan nan ... nan nan nan
Attributes: (12/20)
    Conventions:                     CF-1.6, ACDD-1.3, flopy 3.10.0
    date_created:                    2026-10-21T14:47:51Z
    geospatial_vertical_positive:    up
    geospatial_vertical_min:         55.37115
    geospatial_vertical_max:         107.75425
    geospatial_vertical_resolution:  variable
    ...                              ...
    solver_head_tolerance:           -999
    solver_flux_tolerance:           -999
    flopy_sr_xll:                    198235.0
    flopy_sr_yll:                    7033770.0
    flopy_sr_rotation:               25.0
    start_datetime:                  1970-01-01T00:00:00

print(ds.top.to_numpy())

[[nan nan nan ... nan nan nan]
 [nan nan nan ... nan nan nan]
 [nan nan nan ... nan nan nan]
 ...
 [nan nan nan ... nan nan nan]
 [nan nan nan ... nan nan nan]
 [nan nan nan ... nan nan nan]]

Expected behavior
I would expect actual values in the netCDF files. Is this a bug?

Desktop (please complete the following information):

• OS: Windows 11
• Version 3.10.0

[wpbonelli]

Hi @ThibHlln try this

nc = gwf.export('model.nc', fmt='netcdf')
nc.write()

I think there's a few things going on here. I'll confirm and then follow up.