How to use flopy - 10 common examples

def freyberg_prior_build_test():
    import numpy as np
    import pandas as pd
    pd.set_option('display.max_rows', 500)
    pd.set_option('display.max_columns', 500)
    pd.set_option('display.width', 1000)
    try:
        import flopy
    except:
        return

    org_model_ws = os.path.join("..", "examples", "freyberg_sfr_update")
    nam_file = "freyberg.nam"
    m = flopy.modflow.Modflow.load(nam_file, model_ws=org_model_ws,
                                   check=False, forgive=False,
                                   exe_name=mf_exe_path)
    flopy.modflow.ModflowRiv(m, stress_period_data={
        0: [[0, 0, 0, m.dis.top.array[0, 0], 1.0, m.dis.botm.array[0, 0, 0]],
            [0, 0, 1, m.dis.top.array[0, 1], 1.0, m.dis.botm.array[0, 0, 1]],
            [0, 0, 1, m.dis.top.array[0, 1], 1.0, m.dis.botm.array[0, 0, 1]]]})

    welsp = m.wel.stress_period_data.data.copy()
    addwell = welsp[0].copy()
    addwell['k'] = 1
    welsp[0] = np.rec.array(np.concatenate([welsp[0], addwell]))
    samewell = welsp[1].copy()
    samewell['flux'] *= 10
    welsp[1] = np.rec.array(np.concatenate([welsp[1], samewell]))
    m.wel.stress_period_data = welsp

def freyberg_test():
    import numpy as np
    import pandas as pd
    pd.set_option('display.max_rows', 500)
    pd.set_option('display.max_columns', 500)
    pd.set_option('display.width', 1000)
    try:
        import flopy
    except:
        return

    org_model_ws = os.path.join("..", "examples", "freyberg_sfr_update")
    nam_file = "freyberg.nam"
    m = flopy.modflow.Modflow.load(nam_file, model_ws=org_model_ws,
                                   check=False, forgive=False,
                                   exe_name=mf_exe_path)
    flopy.modflow.ModflowRiv(m, stress_period_data={
        0: [[0, 0, 0, m.dis.top.array[0, 0], 1.0, m.dis.botm.array[0, 0, 0]],
            [0, 0, 1, m.dis.top.array[0, 1], 1.0, m.dis.botm.array[0, 0, 1]],
            [0, 0, 1, m.dis.top.array[0, 1], 1.0, m.dis.botm.array[0, 0, 1]]]})

    org_model_ws = "temp_pst_from"
    if os.path.exists(org_model_ws):
        shutil.rmtree(org_model_ws)
    m.external_path = "."
    m.change_model_ws(org_model_ws)
    m.write_input()
    print("{0} {1}".format(mf_exe_path, m.name + ".nam"), org_model_ws)
    os_utils.run("{0} {1}".format(mf_exe_path, m.name + ".nam"),
                 cwd=org_model_ws)

def pst_from_flopy_geo_draw_test():
    import shutil
    import numpy as np
    import pandas as pd
    try:
        import flopy
    except:
        return
    import pyemu
    org_model_ws = os.path.join("..", "examples", "freyberg_sfr_update")
    nam_file = "freyberg.nam"
    m = flopy.modflow.Modflow.load(nam_file, model_ws=org_model_ws, check=False)
    flopy.modflow.ModflowRiv(m, stress_period_data={0: [[0, 0, 0, 30.0, 1.0, 25.0],
                                                        [0, 0, 1, 31.0, 1.0, 25.0],
                                                        [0, 0, 1, 31.0, 1.0, 25.0]]})
    org_model_ws = "temp"
    m.change_model_ws(org_model_ws)
    m.write_input()

    new_model_ws = "temp_pst_from_flopy"

    hds_kperk = []
    for k in range(m.nlay):
        for kper in range(m.nper):
            hds_kperk.append([kper, k])
    temp_list_props = [["wel.flux", None]]
    spat_list_props = [["riv.cond", 0], ["riv.stage", 0]]
    ph = pyemu.helpers.PstFromFlopyModel(nam_file, new_model_ws=new_model_ws,

def mflist_budget_test():
    import pyemu
    import os
    try:
        import flopy
    except:
        print("no flopy...")
        return
    model_ws = os.path.join("..","examples","Freyberg_transient")
    ml = flopy.modflow.Modflow.load("freyberg.nam",model_ws=model_ws,check=False)
    list_filename = os.path.join(model_ws,"freyberg.list")
    assert os.path.exists(list_filename)
    pyemu.gw_utils.setup_mflist_budget_obs(list_filename,start_datetime=ml.start_datetime)

def kl_test():
    import os
    import numpy as np
    import pandas as pd
    import pyemu
    try:
        import flopy
    except:
        print("flopy not imported...")
        return
    model_ws = os.path.join("..","verification","Freyberg","extra_crispy")
    ml = flopy.modflow.Modflow.load("freyberg.nam",model_ws=model_ws)
    str_file = os.path.join("..","verification","Freyberg","structure.dat")
    arr_dict = {"test":np.ones((ml.nrow,ml.ncol))}
    arr_dict["hk_tru"] = np.loadtxt(os.path.join("..","verification",
                                                 "Freyberg","extra_crispy",
                                                 "hk.truth.ref"))
    basis_file = os.path.join("utils","basis.dat")
    tpl_file = os.path.join("utils","test.tpl")
    back_dict = pyemu.utils.helpers.kl_setup(num_eig=800,sr=ml.sr,
                                             struct_file =str_file,
                                             array_dict=arr_dict,
                                             basis_file=basis_file,
                                             tpl_file=tpl_file)
    for name in back_dict.keys():
        diff = np.abs((arr_dict[name] - back_dict[name])).sum()
        print(diff)
        assert np.abs(diff) < 1.0e-2

def pcg_fmt_test():
    # mf2k container - this will pass
    m2k = fp.modflow.Modflow(version='mf2k')
    m2k.pcg = fp.modflow.ModflowPcg.load(model=m2k, 
                                     f=pcg_fname)
    
    # mf2005 container
    m05 = fp.modflow.Modflow(version='mf2005')
    m05.pcg = fp.modflow.ModflowPcg.load(model=m05, 
                                 f=pcg_fname)
    # this will exit with ValueError without the except block added in pull req
    
    assert m2k.pcg.rclose == m05.pcg.rclose
    assert m2k.pcg.damp == m05.pcg.damp
    
    return

def setup_freyberg_transient_model():


    org_model_ws = os.path.join("..", "examples", "freyberg_sfr_update")
    nam_file = "freyberg.nam"
    mo = flopy.modflow.Modflow.load(nam_file, model_ws=org_model_ws, check=False,forgive=False)

    perlen = np.ones((365))
    m = flopy.modflow.Modflow("freyberg_transient",model_ws=os.path.join("da","freyberg","truth"),version="mfnwt",
                              external_path=".")

    flopy.modflow.ModflowDis(m,nrow=mo.nrow,ncol=mo.ncol,nlay=1,delr=mo.dis.delr,
                             delc=mo.dis.delc,top=mo.dis.top,botm=mo.dis.botm[-1],nper=len(perlen),perlen=perlen)
    flopy.modflow.ModflowBas(m,ibound=mo.bas6.ibound[0],strt=mo.bas6.strt[0])
    flopy.modflow.ModflowUpw(m,laytyp=mo.upw.laytyp,hk=mo.upw.hk[0],vka=mo.upw.vka[0],ss=0.00001,sy=0.01)
    flopy.modflow.ModflowNwt(m)
    oc_data = {}
    for iper in range(m.nper):
        oc_data[iper,0] = ["save head","save budget"]
    flopy.modflow.ModflowOc(m,stress_period_data=oc_data)
    flopy.modflow.ModflowRch(m,rech=mo.rch.rech.array[0])
    wel_data = mo.wel.stress_period_data[0]

return
    import pyemu

    ext = ''
    bin_path = os.path.join("..", "..", "bin")
    if "linux" in platform.platform().lower():
        bin_path = os.path.join(bin_path, "linux")
    elif "darwin" in platform.platform().lower():
        bin_path = os.path.join(bin_path, "mac")
    else:
        bin_path = os.path.join(bin_path, "win")
        ext = '.exe'

    org_model_ws = os.path.join("..", "examples", "freyberg_sfr_update")
    nam_file = "freyberg.nam"
    m = flopy.modflow.Modflow.load(nam_file, model_ws=org_model_ws, check=False, forgive=False,
                                   exe_name=mf_exe_name)
    org_model_ws = "temp"

    m.change_model_ws(org_model_ws)
    m.write_input()
    print("{0} {1}".format(mf_exe_name, m.name + ".nam"), org_model_ws)
    pyemu.os_utils.run("{0} {1}".format(mf_exe_name, m.name + ".nam"), cwd=org_model_ws)
    hds_file = "freyberg.hds"
    list_file = "freyberg.list"
    for f in [hds_file, list_file]:
        assert os.path.exists(os.path.join(org_model_ws, f))

    new_model_ws = "template1"

    props = [["upw.hk", None], ["upw.vka", None], ["upw.ss", None], ["rch.rech", None]]

# except:
    #     pass
    os.chdir(bd)

    org_model_ws = os.path.join("..", "examples", "freyberg_sfr_update")
    nam_file = "freyberg.nam"
    m = flopy.modflow.Modflow.load(nam_file, model_ws=org_model_ws, check=False)

    helper = pyemu.helpers.PstFromFlopyModel(nam_file, new_model_ws, org_model_ws,
                                             hds_kperk=[0, 0], remove_existing=True,
                                             model_exe_name="mfnwt", sfr_pars=True, sfr_obs=True)
    pe = helper.draw(100)

    # go again testing passing list to sfr_pars
    m = flopy.modflow.Modflow.load(nam_file, model_ws=org_model_ws, check=False)

    helper = pyemu.helpers.PstFromFlopyModel(nam_file, new_model_ws, org_model_ws,
                                             hds_kperk=[0, 0], remove_existing=True,
                                             model_exe_name="mfnwt", sfr_pars=['flow', 'not_a_par'], sfr_obs=True)
    try:
        pe = helper.draw(100)
    except:
        pass
    else:
        raise Exception()



    # go again passing bumph to sfr_par
    m = flopy.modflow.Modflow.load(nam_file, model_ws=org_model_ws, check=False)

def test_rasters():
    from flopy.utils import Raster
    import os
    import flopy as fp

    ws = os.path.join("..", "examples", "data", "options")
    raster_name = "dem.img"

    try:
        rio = Raster.load(os.path.join(ws, "dem", raster_name))
    except ImportError:
        return

    ml = fp.modflow.Modflow.load("sagehen.nam", version="mfnwt",
                                 model_ws=os.path.join(ws, 'sagehen'))
    xoff = 214110
    yoff = 4366620
    ml.modelgrid.set_coord_info(xoff, yoff)

    # test sampling points and polygons
    val = rio.sample_point(xoff + 2000, yoff + 2000, band=1)
    print(val - 2336.3965)
    if abs(val - 2336.3965) > 1e-4:
        raise AssertionError

    x0, x1, y0, y1 = rio.bounds

    x0 += 1000
    y0 += 1000
    x1 -= 1000