Writing a .gsf file from a MODFLOW-USG quadtree model built using FloPy

Question

I've built a MODFLOW USG-Transport model in FloPY (quadtree grid using Gridgen) and I want to be able to open it in a GUI such as Vistas so that others who don't use FloPy can use it.

I need a .gsf file to be able to do that. I tried using the PEST utility gridgen2gsf but not all nodal connections are included in that file.

Has anyone written a gsf file for a USG model and be willing to give me hand?

cheers, Emma

Answer 1

If it helps, the Groundwater Data Utility (or PEST utility) "gridgen2gsf" written by John Doherty is made specifically for this purpose - to create *.GSF files from *.DISU files with Quadtree meshes.

It can be downloaded for free from https://pesthomepage.org/

Note that this utility does not work with Voronoi meshes.

Answer 2

I ran into the same problem. I needed a file with the extension .gsf for importing the flopy model to gms. I wrote the following code:

g=Gridgen(dis,exe_name='gridgen_x64', model_ws=gridgen_dir)
ncells = g.get_nodes()
verts, iverts = g.get_verts_iverts(ncells)
# vertices=g.get_vertices(ncells)
vertices = g.get_gridprops_vertexgrid()
lenverts=len(verts)
f = open("test.gsf", "a")
f.write('UNSTRUCTURED\n' + str(g.get_nodes()) + ' ' + str(g.get_nlay()) + ' ' + 
        str(1) + ' ' + str(1) + '\n' + str(
        lenverts * (nlay+1)) + '\n')

tuples = []
count = 0
for i in g.get_top():
    if count != ncells:
       for j in g._vertdict[count]:
           # np.savetxt(f, (j+(i,),), delimiter=' ', fmt='%g')
           tuples.append(j + (i,))
    count += 1
count=ncells-ncells/nlay
for i in list(g.get_bot()[int(ncells-ncells/nlay):int(ncells)]):
    for j in g._vertdict[count]:
        tuples.append(j + (i,))
    count += 1

nvertexes = list(dict.fromkeys(tuples))

for i in range(0, ncells):
    iverts[i] = [x + 1 for x in iverts[i]]
    del iverts[i][-1]
for i in nvertexes:
    np.savetxt(f, (i,), delimiter=' ', fmt='%g')
nodesgfs = []
nodeiter=0
for j in range(0,nlay):
    for i in g.get_nod_recarray()[nodeiter:]:
        ivertex=[x + lenverts*(j) for x in iverts[nodeiter]]
        nodesgfs.append((i[0] + 1,) + (i[2],) + (i[3],) + (i[4],) + (i[1] + 1,) 
                       + (len(ivertex)*2,) + tuple([x + lenverts*(j+1) for x in 
                        iverts[nodeiter]])+tuple(ivertex))
        nodeiter += 1
        if nodeiter in [numcell for numcell in range(int(ncells/nlay),ncells,int(ncells/nlay))]:
            break

nodesgfs = list(dict.fromkeys(nodesgfs))
for i in nodesgfs:
    np.savetxt(f, (i,), delimiter=' ', fmt='%g')