How to use the flopy.modflow.ModflowBas function in flopy

'''

    ws = os.path.join(model_ws, 'mf2005')
    nm = 'ex01_mf2005'
    exe_name = exe_names['mf2005']
    iu_cbc = 130
    m = flopy.modflow.Modflow(nm, model_ws=ws,
                              exe_name=exe_name)
    flopy.modflow.ModflowDis(m, nlay=nlay, nrow=nrow, ncol=ncol,
                             nper=nper, itmuni=4, lenuni=1,
                             perlen=perlen, nstp=nstp,
                             tsmult=tsmult, steady=True,
                             delr=delr, delc=delc,
                             top=top, botm=botm)
    flopy.modflow.ModflowLpf(m, ipakcb=iu_cbc, laytyp=laytyp, hk=kh, vka=kv)
    flopy.modflow.ModflowBas(m, ibound=1, strt=top)
    # recharge
    flopy.modflow.ModflowRch(m, ipakcb=iu_cbc, rech=rch, nrchop=1)
    # wel
    wd = [i for i in wel_loc] + [wel_q]
    flopy.modflow.ModflowWel(m, ipakcb=iu_cbc, stress_period_data={0: wd})
    # river
    rd = []
    for i in range(nrow):
        rd.append([0, i, ncol - 1, riv_h, riv_c, riv_z])
    flopy.modflow.ModflowRiv(m, ipakcb=iu_cbc, stress_period_data={0: rd})
    # output control
    flopy.modflow.ModflowOc(m, stress_period_data={(0, 0): ['save head',
                                                            'save budget',
                                                            'print head']})
    flopy.modflow.ModflowPcg(m, hclose=1e-6, rclose=1e-3, iter1=100, mxiter=50)

ibnd3 = np.array(ibnd3)
    ibnd1 = ibnd1.reshape(top1.shape)
    ibnd2 = ibnd2.reshape(top1.shape)
    ibnd3 = ibnd3.reshape(top1.shape)
    ibnd = [ibnd1, ibnd2, ibnd3]
    ibnd = np.array(ibnd)
    
    sthd1 = np.array(sthd1)
    sthd1 = sthd1.reshape(top1.shape)
    sthd2 = np.ones(top1.shape) * 210.
    sthd3 = sthd2.copy()
    sthd = [sthd1, sthd2, sthd3]
    sthd = np.array(sthd)
    hdry = -9999
    
    bas = flopy.modflow.ModflowBas(mf, ibound=ibnd, hnoflo=hdry, strt=sthd)
    
    ## Instantiate Lake (LAK) package for MODFLOW-NWT
    
    lkarr1 = [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,

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]
    wel_data["k"][:] = 0
    flopy.modflow.ModflowWel(m,stress_period_data={0:wel_data})
    flopy.modflow.ModflowSfr2(m,nstrm=mo.sfr.nstrm,nss=mo.sfr.nss,istcb2=90,segment_data=mo.sfr.segment_data,reach_data=mo.sfr.reach_data)
    m.write_input()
    pyemu.os_utils.run("mfnwt {0}.nam".format(m.name),cwd=m.model_ws)
    return m

botm = np.arange(-1,-4,-1)
    top = 0.
    
    # creating MODFLOW model
    
    model_ws = os.path.join('.', 'temp', 't068a')
    modelname = 'model_mf'
    
    mf = flopy.modflow.Modflow(modelname, model_ws=model_ws,
                               exe_name=mf_exe_name)
    dis = flopy.modflow.ModflowDis(mf, nlay=nlay, nrow=nrow, ncol=ncol, 
                                   perlen=perlen, nper=nper, botm=botm, top=top, 
                                   steady=False)
    
    bas = flopy.modflow.ModflowBas(mf, ibound=ibound, strt=top)
    lpf = flopy.modflow.ModflowLpf(mf, hk=100, vka=100, ss=0.00001, sy=0.1)
    oc = flopy.modflow.ModflowOc(mf)
    pcg = flopy.modflow.ModflowPcg(mf)
    lmt = flopy.modflow.ModflowLmt(mf)
    
    # recharge
    rchrate = {}
    rchrate[0] = 0.0
    rchrate[5] = 0.001
    rchrate[6] = 0.0
    
    rch = flopy.modflow.ModflowRch(mf, rech=rchrate, nrchop=3)
    
    # define itype
    itype = flopy.mt3d.Mt3dSsm.itype_dict()
    ssm_data = {}

'rclose': 3e-3, 'relax': 1.0,
                     'nbpol': 2, 'damp': 1.0, 'dampt': 1.0}

    # --create model file and run model
    modelname = 'swi2ex5'
    mf_name = 'mf2005'
    if not skipRuns:
        ml = flopy.modflow.Modflow(modelname, version='mf2005',
                                   exe_name=mf_name,
                                   model_ws=dirs[0])
        discret = flopy.modflow.ModflowDis(ml, nrow=nrow, ncol=ncol, nlay=nlay,
                                           delr=delr, delc=delc, top=0,
                                           botm=bot,
                                           laycbd=0, nper=nper, perlen=perlen,
                                           nstp=nstp, steady=steady)
        bas = flopy.modflow.ModflowBas(ml, ibound=ibound, strt=ihead)
        lpf = flopy.modflow.ModflowLpf(ml, hk=kh, vka=kv, ss=ss, sy=ssz,
                                       vkcb=0,
                                       laytyp=0, layavg=1)
        wel = flopy.modflow.ModflowWel(ml, stress_period_data=well_data)
        swi = flopy.modflow.ModflowSwi2(ml, iswizt=55, npln=1, istrat=1,
                                        toeslope=0.025, tipslope=0.025,
                                        nu=[0, 0.025], zeta=z, ssz=ssz,
                                        isource=isource, nsolver=2,
                                        solver2params=solver2params)
        oc = flopy.modflow.ModflowOc(ml, stress_period_data=ocspd)
        pcg = flopy.modflow.ModflowPcg(ml, hclose=1.0e-6, rclose=3.0e-3,
                                       mxiter=100, iter1=50)
        # --write the modflow files
        ml.write_input()
        m = ml.run_model(silent=False)

dmcoef = 0.57024  # m2/day  Could also try 1.62925 as another case of the Henry problem
    hk = 864.  # m/day

    # Create the basic MODFLOW model data
    modelname = 'henry'
    m = flopy.seawat.Seawat(modelname, exe_name="swtv4", model_ws=workspace)

    # Add DIS package to the MODFLOW model
    dis = flopy.modflow.ModflowDis(m, nlay, nrow, ncol, nper=1, delr=delr,
                                   delc=delc, laycbd=0, top=henry_top,
                                   botm=henry_botm, perlen=1.5, nstp=15)

    # Variables for the BAS package
    ibound = np.ones((nlay, nrow, ncol), dtype=np.int32)
    ibound[:, :, -1] = -1
    bas = flopy.modflow.ModflowBas(m, ibound, 0)

    # Add LPF package to the MODFLOW model
    lpf = flopy.modflow.ModflowLpf(m, hk=hk, vka=hk, ipakcb=53)

    # Add PCG Package to the MODFLOW model
    pcg = flopy.modflow.ModflowPcg(m, hclose=1.e-8)

    # Add OC package to the MODFLOW model
    oc = flopy.modflow.ModflowOc(m,
                                 stress_period_data={
                                     (0, 0): ['save head', 'save budget']},
                                 compact=True)

    # Create WEL and SSM data
    itype = flopy.mt3d.Mt3dSsm.itype_dict()
    wel_data = {}

z.append(z1)
    ssz = 0.2
    isource = np.ones((nrow, ncol), 'int')
    isource[0, 0] = 2

    # stratified model
    modelname = 'swiex2_strat'
    print('creating...', modelname)
    ml = flopy.modflow.Modflow(modelname, version='mf2005', exe_name=mf_name,
                               model_ws=dirs[0])
    discret = flopy.modflow.ModflowDis(ml, nlay=1, ncol=ncol, nrow=nrow,
                                       delr=delr,
                                       delc=1,
                                       top=0, botm=[-40.0],
                                       nper=nper, perlen=perlen, nstp=nstp)
    bas = flopy.modflow.ModflowBas(ml, ibound=ibound, strt=0.05)
    bcf = flopy.modflow.ModflowBcf(ml, laycon=0, tran=2 * 40)
    swi = flopy.modflow.ModflowSwi2(ml, iswizt=55, nsrf=nsurf, istrat=1,
                                    toeslope=0.2, tipslope=0.2,
                                    nu=[0, 0.0125, 0.025],
                                    zeta=z, ssz=ssz, isource=isource,
                                    nsolver=1)
    oc = flopy.modflow.ModflowOc(ml,
                                 stress_period_data={(0, 999): ['save head']})
    pcg = flopy.modflow.ModflowPcg(ml)
    ml.write_input()
    # run stratified model
    if not skipRuns:
        m = ml.run_model(silent=False)
    # read stratified results
    zetafile = os.path.join(dirs[0], '{}.zta'.format(modelname))
    zobj = flopy.utils.CellBudgetFile(zetafile)

ssm_data = {0: [0, 0, 0, 35., itype['BAS6']]}

    # print sconcp
    # mt3d print times
    timprs = (np.arange(5) + 1) * 2000.
    nprs = len(timprs)
    # create the MODFLOW files
    m = flopy.seawat.Seawat(modelname, exe_name=swtexe_name, model_ws=dirs[1])
    discret = flopy.modflow.ModflowDis(m, nrow=swt_nrow, ncol=swt_ncol,
                                       nlay=swt_nlay,
                                       delr=swt_delr, delc=swt_delc, laycbd=0,
                                       top=swt_top,
                                       botm=swt_botm,
                                       nper=nper, perlen=perlen, nstp=1,
                                       steady=False)
    bas = flopy.modflow.ModflowBas(m, ibound=swt_ibound, strt=0.05)
    lpf = flopy.modflow.ModflowLpf(m, hk=2.0, vka=2.0, ss=0.0, sy=0.0,
                                   laytyp=0,
                                   layavg=0)
    oc = flopy.modflow.ModflowOc(m, save_every=1, save_types=['save head'])
    pcg = flopy.modflow.ModflowPcg(m)
    # Create the MT3DMS model files
    adv = flopy.mt3d.Mt3dAdv(m, mixelm=-1,  # -1 is TVD
                             percel=0.05,
                             nadvfd=0,
                             # 0 or 1 is upstream; 2 is central in space
                             # particle based methods
                             nplane=4,
                             mxpart=1e7,
                             itrack=2,
                             dceps=1e-4,
                             npl=16,

zini = np.hstack(
        (-9 * np.ones(24), np.arange(-9, -50, -0.5), -50 * np.ones(94)))[
           np.newaxis, :]
    iso = np.zeros((1, 200), dtype=np.int)
    iso[:, :30] = -2
    # model objects
    ml = flopy.modflow.Modflow(modelname, version='mf2005', exe_name=exe_name,
                               model_ws=workspace)
    discret = flopy.modflow.ModflowDis(ml, nrow=nrow, ncol=ncol, nlay=3,
                                       delr=delr,
                                       delc=delc,
                                       laycbd=[0, 0, 0], top=-9.0,
                                       botm=[-29, -30, -50],
                                       nper=2, perlen=[365 * 1000, 1000 * 365],
                                       nstp=[500, 500])
    bas = flopy.modflow.ModflowBas(ml, ibound=1, strt=1.0)
    bcf = flopy.modflow.ModflowBcf(ml, laycon=[0, 0, 0], tran=[40.0, 1, 80.0],
                                   vcont=[0.005, 0.005])
    wel = flopy.modflow.ModflowWel(ml, stress_period_data={0: lrcQ1, 1: lrcQ2})
    ghb = flopy.modflow.ModflowGhb(ml, stress_period_data={0: lrchc})
    swi = flopy.modflow.ModflowSwi2(ml, iswizt=55, nsrf=1, istrat=1,
                                    toeslope=0.01, tipslope=0.04,
                                    nu=[0, 0.025],
                                    zeta=[zini, zini, zini], ssz=0.2,
                                    isource=iso,
                                    nsolver=1)
    oc = flopy.modflow.ModflowOc(ml, save_every=100, save_types=['save head'])
    pcg = flopy.modflow.ModflowPcg(ml)
    # write the model files
    ml.write_input()
    # run the model
    m = ml.run_model(silent=True)

delc = Ly/nrow
delv = (ztop - zbot) / nlay
botm = np.linspace(ztop, zbot, nlay + 1)

# Create the discretization object
dis = flopy.modflow.ModflowDis(mf, nlay, nrow, ncol, delr=delr, delc=delc,
                               top=ztop, botm=botm[1:])

# Variables for the BAS package
ibound = np.ones((nlay, nrow, ncol), dtype=np.int32)
ibound[:, :, 0] = -1
ibound[:, :, -1] = -1
strt = np.ones((nlay, nrow, ncol), dtype=np.float32)
strt[:, :, 0] = 10.
strt[:, :, -1] = 0.
bas = flopy.modflow.ModflowBas(mf, ibound=ibound, strt=strt)

# Add LPF package to the MODFLOW model
lpf = flopy.modflow.ModflowLpf(mf, hk=10., vka=10., ipakcb=53)

# Add OC package to the MODFLOW model
spd = {(0, 0): ['print head', 'print budget', 'save head', 'save budget']}
oc = flopy.modflow.ModflowOc(mf, stress_period_data=spd, compact=True)

# Add PCG package to the MODFLOW model
pcg = flopy.modflow.ModflowPcg(mf)

# Write the MODFLOW model input files
mf.write_input()

# Run the MODFLOW model
success, buff = mf.run_model()