Parallelize Flopy GWF 3D model with MF6Splitter

[MVergnd]

Hi there,

I am trying to use Flopy to run a MODFLOW6 model in parallel. I managed to split my model into 4 different models with the function MF6Splitter using a previous discussion on GitHub : #2316

But in my case I have many layers, so I have done this :

mask = np.ones((grid_design.nlay, grid_design.nrow, grid_design.ncol ), dtype=int)
split_row = int(grid_design.nrow/2)
split_col = int(grid_design.ncol/2)

mask[:, 0:split_row, 0:split_col] = 1
mask[:, 0:split_row, split_col:grid_design.ncol] = 2
mask[:, split_row:grid_design.nrow, 0:split_col] = 3
mask[:, split_row:grid_design.nrow, split_col:grid_design.ncol] = 4

mfsplit = Mf6Splitter(sim)
new_sim = mfsplit.split_model(mask)
new_sim.write_simulation()

But I have this error :
--> [217]new_sim = mfsplit.split_model(mask) # split the model, rename sim
[218] new_sim.write_simulation()
...
--> [650]check = arr[mask]
[651] count += np.count_nonzero(check)
[652] if count == 0:

IndexError: index 28167 is out of bounds for axis 0 with size 28167

Which I do not understand because the mask shape is equal to the grid shape. Any help would be greatly appreciated!

Thanks a lot, Maëlle

[jdhughes-usgs]

The splitter works off a 2D mask array. You just have to create a 2D mask array instead of a 3D mask array.

[MVergnd]

Well, thank you for your prompt reply, that's what I did in the first place :

mask = np.ones((grid_design.nrow, grid_design.ncol), dtype=int)
split_row = int(grid_design.nrow/2)
split_col = int(grid_design.ncol/2)

mask[0:split_row, 0:split_col] = 1
mask[0:split_row, split_col:grid_design.ncol] = 2
mask[split_row:grid_design.nrow, 0:split_col] = 3
mask[split_row:grid_design.nrow, split_col:grid_design.ncol] = 4
mfsplit = Mf6Splitter(sim)
new_sim = mfsplit.split_model(mask)
new_sim.write_simulation()

but I have got this error :

--> [217](new_sim = mfsplit.split_model(mask)
...
-> [2290] obsid.append(int(i) - 1)
[2291] except ValueError:
[2292] obsid.append(i)

TypeError: int() argument must be a string, a bytes-like object or a real number, not 'tuple'

[jdhughes-usgs]

@MVergnd Are you using boundnames in your observations? If so do the boundaries associated with boundname exist in more than one domain?

[MVergnd]

@jdhughes-usgs , I'm not using boundnames. And yes the boundaries exist in more than one domain.