|
# 导入Flask类
|
import numpy as np
|
import flopy
|
import flopy.utils.binaryfile as bf
|
import os
|
import json
|
import Base as base
|
import CalHead
|
import ModelPeriod
|
import OpenExcel
|
|
|
base_init_year=["2020","2021","2022"]
|
river_start_index = 454
|
river_end_index =562
|
|
#预测周期数
|
predict_per = 12
|
|
#降水量
|
|
base_water = base.prefix + 'base_water.ini'
|
def predict_water_chart(base_year,start_time ,end_time,value):
|
water_array = np.loadtxt(base_water, dtype=str,encoding='utf-8')
|
y_data=[]
|
x_data= ModelPeriod.get_months_in_year()
|
water= water_array[int(base_year)]
|
for data in water:
|
float_data = float(data)/9
|
float_data= round(float_data*value,2)
|
y_data.append(float_data)
|
|
result = {"y_data": y_data, "x_data": x_data}
|
return result
|
|
#河流的折线图
|
|
base_river = base.prefix + 'base_river.ini'
|
def predict_river_chart(base_year,start_time ,end_time,value):
|
river_array = np.loadtxt(base_river, dtype=str,encoding='utf-8')
|
y_data=[]
|
x_data= ModelPeriod.get_months_in_year()
|
for data in river_array:
|
#万立方米
|
float_data = float(data)/4/10000
|
float_data= round(float_data*value,2)
|
y_data.append(float_data)
|
|
result = {"y_data": y_data, "x_data": x_data}
|
return result
|
|
|
#开采量折线图
|
base_mining = base.prefix + 'base_mining.ini'
|
def predict_well_chart(base_year,start_time ,end_time,data):
|
# 2022/海淀/昌平/门头沟/石景山
|
area_names=["全部区域","海淀区","昌平区","门头沟区","石景山区"]
|
mining_array = np.loadtxt(base_mining, dtype=str,encoding='utf-8')
|
x_data= ModelPeriod.get_months_in_year()
|
result = {"x_data": x_data}
|
y_data=[]
|
for item in data:
|
name = item["area_name"]
|
value = float(item["value"])
|
if name =='朝阳区':
|
continue
|
index = area_names.index(name)
|
row_data = mining_array[index]
|
float_data=[]
|
for item in row_data:
|
x = round(float(item)/10000*value,2)
|
float_data.append(x)
|
dicts={"name":name,"data":float_data}
|
y_data.append(dicts)
|
result["y_data"] = y_data
|
return result
|
|
|
def run_model(model_name):
|
if model_name == base.not_allowed_model or model_name in base.archive_models:
|
return "本模型为验证模型,不允许修改！"
|
|
dicts= run_model_predict(model_name)
|
if dicts["code"] == 400:
|
return dicts
|
#导出csv文件
|
CalHead.exportCsV(model_name)
|
|
# #更新模型三维网格配置
|
base.updateModelConfig(model_name)
|
|
# #创建模型的三维网格
|
filedir = base.model3d_path + model_name
|
|
if not os.path.exists(filedir):
|
os.makedirs(filedir, exist_ok=True)
|
|
base.callModelexe()
|
#计算水资源量和水均衡
|
CalHead.run_zonebudget_bal(model_name)
|
CalHead.run_zonebudget_res(model_name)
|
dicts= {"code":200,"msg":"预测模型运行成功！" }
|
return dicts
|
|
|
def run_model_predict(model_name):
|
|
predictiondata=""
|
prediction_path = base.model_dir + model_name +"\\prediction.json"
|
if os.path.exists(prediction_path):
|
with open(prediction_path,encoding='utf-8') as f:
|
predictiondata = json.load(f)
|
|
periods =CalHead.get_model_period(model_name)
|
periods_len= len(periods)
|
|
if predictiondata:
|
updateDisFile(model_name,periods_len)
|
try:
|
flag = updateBase6File(model_name,predictiondata)
|
if bool(flag)==False:
|
dicts= {"code":400,"msg":"请检查模型的初始水是否设置！" }
|
return dicts
|
except:
|
dicts= {"code":400,"msg":"请检查模型的初始水是否设置！" }
|
return dicts
|
|
try:
|
updateRchFile(model_name,predictiondata)
|
except:
|
print("RchFile无预测参数，无需修改！")
|
|
try:
|
updateRiverFile(model_name,predictiondata)
|
except:
|
print("RiverFile无预测参数，无需修改！")
|
else:
|
print("prediction.json 预测场景文件为空，无需更改相应文件")
|
|
|
model_ws = base.model_dir + model_name
|
|
ml = flopy.modflow.Modflow.load("modflow.nam", model_ws=model_ws,
|
exe_name="mf2005", verbose=True, version="mf2005", check=False)
|
ml.run_model(report = True)
|
dicts= {"code":200,"msg":"预测模型运行成功！" }
|
return dicts
|
|
|
|
#获取 area的 name--> ratio 的结构
|
def get_area_dict(area):
|
result ={}
|
|
for i in range(len(area)):
|
name = area[i]["name"]
|
rt = area[i]["ratio"]
|
result[name]= rt
|
return result
|
|
|
#获取区县的 row+column --> name结构
|
def get_distric_dict():
|
data = base.district
|
result = {}
|
for row ,column ,id ,name in data:
|
key = str(row)+","+str(column)
|
result[key]= name
|
return result
|
|
|
#根据 row clomn 获取 ratio
|
def get_row_column_ratio(row, column ,district_dict, area_dict ):
|
key = str(row) +"," + str(column)
|
if area_dict.__contains__("全部区域"):
|
return area_dict["全部区域"]
|
|
if district_dict.__contains__(key):
|
name = district_dict[key]
|
ratio = area_dict[name]
|
return float(ratio)
|
|
return float(1.0)
|
|
|
|
|
def updateRiverFile(model_name,predictiondata):
|
|
flag = check_rain_param(predictiondata)
|
|
if flag == "true":
|
|
rain_ratio = float(predictiondata["rain"]["ratio"])
|
rain_base_year = predictiondata["rain"]["base_year"]
|
if rain_base_year=='4':
|
rain_base_year="1"
|
|
river_ratio= float(predictiondata["river"]["ratio"])
|
area= predictiondata["mine"]["area"]
|
|
# ws = base.predictParamModel + rain_base_year
|
ws=""
|
start_time = predictiondata["start_time"]
|
end_time = predictiondata["end_time"]
|
count = ModelPeriod.get_months_in_range_count(start_time, end_time)
|
if count==12:
|
ws= base.predictParamModel + rain_base_year
|
else:
|
ws= base.predictModel60
|
|
baseMdoel = flopy.modflow.Modflow.load("modflow.nam", model_ws= ws,
|
exe_name="mf2005", verbose=True, version="mf2005", check=False)
|
|
update_model_ws = base.model_dir + model_name
|
updateMdoel = flopy.modflow.Modflow.load("modflow.nam", model_ws= update_model_ws,
|
exe_name="mf2005", verbose=True, version="mf2005", check=False)
|
|
district_dict = get_distric_dict()
|
|
area_dict = get_area_dict(area)
|
|
lrcq = {}
|
|
periods =CalHead.get_model_period(model_name)
|
periods_len= len(periods)
|
|
for per in range(periods_len):
|
wel = []
|
array2d = []
|
|
wel = baseMdoel.wel.stress_period_data.__getitem__(kper = per)
|
wel_len = len(wel)
|
|
#侧向边界
|
for i in range (0,453):
|
wel[i][3] = wel[i][3] * rain_ratio
|
|
#河流
|
for i in range(453, 562):
|
wel[i][3] = wel[i][3] * river_ratio
|
|
#抽水井
|
for i in range(562,wel_len):
|
|
r = (float) (get_row_column_ratio(wel[i][1], wel[i][2], district_dict, area_dict))
|
wel[i][3] = wel[i][3] * r
|
|
|
#重置数组
|
for Layer, Row, Column, Q in wel:
|
array = [Layer, Row, Column, Q]
|
array2d.append(array)
|
|
#追加井
|
flex_data= getFlexdata(model_name)
|
print("==============")
|
print(flex_data)
|
|
for i in range(len(flex_data)):
|
array2d.append(flex_data[i])
|
|
lrcq[per] = array2d
|
|
flopy.modflow.ModflowWel(updateMdoel,
|
ipakcb= baseMdoel.wel.ipakcb,
|
dtype=baseMdoel.wel.dtype,
|
options=baseMdoel.wel.options,
|
stress_period_data=lrcq)
|
updateMdoel.write_input()
|
|
else:
|
print("Well--River文件无需修改！")
|
|
#追加复杂原汇项信息
|
def getFlexdata(model_name):
|
welldata=""
|
well_path = base.model_dir + model_name +"\\pump_well.json"
|
data=[]
|
if os.path.exists(well_path):
|
with open(well_path,encoding='utf-8') as f:
|
welldata = json.load(f)
|
wel= welldata["well"]
|
|
for i in range (len(wel)):
|
layer = int (wel[i]['layer'])-1
|
row= int(wel[i]['row'])-1
|
column = int(wel[i]['column'])-1
|
v = float(wel[i]['value'])
|
arr = [layer,row, column, v]
|
data.append(arr)
|
|
return data
|
|
|
def updateRchFile(model_name,predictiondata):
|
flag = check_rain_param(predictiondata)
|
if flag == "true":
|
#丰水年 枯水年
|
base_year = predictiondata["rain"]["base_year"]
|
if base_year =="1" or base_year =="2" or base_year =="3":
|
updateRchBaseYear(model_name,predictiondata)
|
elif base_year =="4":
|
#接入传感器数据
|
updateRchRealData(model_name,predictiondata)
|
|
else:
|
print("Rch文件无需修改！")
|
|
#更新平水年 丰水年 枯水年数据
|
def updateRchBaseYear(model_name,predictiondata):
|
#丰水年 枯水年
|
base_year = predictiondata["rain"]["base_year"]
|
ratio= float(predictiondata["rain"]["ratio"])
|
#数据来源的模型文件夹
|
base_ws=""
|
start_time = predictiondata["start_time"]
|
end_time = predictiondata["end_time"]
|
count = ModelPeriod.get_months_in_range_count(start_time, end_time)
|
if count==12:
|
base_ws= base.predictParamModel + base_year
|
else:
|
base_ws= base.predictModel60
|
|
baseMdoel = flopy.modflow.Modflow.load("modflow.nam", model_ws= base_ws,
|
exe_name="mf2005", verbose=True, version="mf2005", check=False)
|
|
update_model_ws = base.model_dir + model_name
|
updateMdoel = flopy.modflow.Modflow.load("modflow.nam", model_ws= update_model_ws,
|
exe_name="mf2005", verbose=True, version="mf2005", check=False)
|
|
|
periods =CalHead.get_model_period(model_name)
|
periods_len= len(periods)
|
for per in range(periods_len):
|
item = baseMdoel.rch.rech.__getitem__(kper = per)
|
array2d = item.get_value()
|
array2d_len = len(array2d)
|
|
for i in range(array2d_len):
|
array_len = len(array2d[i])
|
for j in range(array_len):
|
if str(base.area_array[i][j]) != '-9999':
|
array2d[i][j] = array2d[i][j] * ratio
|
|
updateMdoel.rch.rech.__setitem__(key = per, value=array2d)
|
rch = flopy.modflow.ModflowRch(updateMdoel,
|
nrchop=baseMdoel.rch.nrchop,
|
ipakcb=baseMdoel.rch.ipakcb,
|
rech=updateMdoel.rch.rech,
|
irch =baseMdoel.rch.irch)
|
|
rch.write_file(check=False)
|
|
|
#接入传感器数据
|
def updateRchRealData(model_name,predictiondata):
|
|
# 取平水年的rch文件
|
base_year = "1"
|
ratio= float(predictiondata["rain"]["ratio"])
|
if not ratio:
|
ratio = 1
|
excel_data = OpenExcel.read_jyl_excel(model_name)
|
array_data =[]
|
|
if not excel_data:
|
print("降水模板数据未获取到！")
|
return "降水模板数据未获取到！"
|
for i in range(1,len(excel_data)):
|
temp =[]
|
for j in range(1,len(excel_data[i])):
|
data = round(float(excel_data[i][j]),8)
|
temp.append(data)
|
array_data.append(temp)
|
|
base_ws=""
|
start_time = predictiondata["start_time"]
|
end_time = predictiondata["end_time"]
|
count = ModelPeriod.get_months_in_range_count(start_time, end_time)
|
if count==12:
|
base_ws= base.predictParamModel + base_year
|
else:
|
base_ws= base.predictModel60
|
|
baseMdoel = flopy.modflow.Modflow.load("modflow.nam", model_ws= base_ws,
|
exe_name="mf2005", verbose=True, version="mf2005", check=False)
|
|
update_model_ws = base.model_dir + model_name
|
|
updateMdoel = flopy.modflow.Modflow.load("modflow.nam", model_ws= update_model_ws,
|
exe_name="mf2005", verbose=True, version="mf2005", check=False)
|
|
#降雨分区
|
area_dicts = base.getAreas()
|
|
|
periods =CalHead.get_model_period(model_name)
|
periods_len= len(periods)
|
|
for per in range(periods_len):
|
#16个分区某个周期的 数据 下标(0,15) 和area_dicts字典降雨分区的key对应
|
water_per_data = get_Array2_column(array_data,per)
|
|
item = baseMdoel.rch.rech.__getitem__(kper = per)
|
array2d = item.get_value()
|
for key in area_dicts:
|
tuples= area_dicts[key]
|
values = water_per_data[int(key)]
|
for i in range(len(tuples)):
|
x = tuples[i][0]
|
y = tuples[i][1]
|
array2d[x][y]= values*ratio
|
|
updateMdoel.rch.rech.__setitem__(key = per, value=array2d)
|
|
rch = flopy.modflow.ModflowRch(updateMdoel,
|
nrchop=baseMdoel.rch.nrchop,
|
ipakcb=baseMdoel.rch.ipakcb,
|
rech=updateMdoel.rch.rech,
|
irch =baseMdoel.rch.irch)
|
|
rch.write_file(check=False)
|
print("降水分区数据更新完毕!")
|
return "降水分区数据更新完毕！"
|
|
#获取列数据
|
def get_Array2_column(array_data,column):
|
arr = np.array(array_data)
|
column_data = arr[:, column]
|
return column_data
|
|
|
def check_rain_param(predictiondata):
|
|
rain = predictiondata["rain"]
|
if not rain:
|
print("Rch预测参数为空，无需要修改")
|
return "false"
|
|
base_year = predictiondata["rain"]["base_year"]
|
if not base_year :
|
print(" Rch : base_year预测参数为空，无需要修改")
|
return "false"
|
|
ratio= predictiondata["rain"]["ratio"]
|
if not ratio :
|
print(" Rch : ratio预测参数为空，无需要修改")
|
return "false"
|
|
return "true"
|
|
|
#更新bas6文件 初始水头信息
|
def updateBase6File(model_name,predictdata):
|
model_ws = base.model_dir + model_name
|
ml = flopy.modflow.Modflow.load("modflow.nam", model_ws=model_ws,
|
exe_name="mf2005", verbose=True, version="mf2005", check=False)
|
#初始水头
|
|
if "initHeader" not in predictdata:
|
print("=============has no initHeader=============")
|
return False
|
init_header = predictdata["initHeader"]
|
|
dir = base.model_dir + init_header + "\\modflow.head"
|
head = bf.HeadFile(dir)
|
alldata = head.get_alldata()
|
|
last_index = len(alldata)-1
|
last_array3= alldata[last_index]
|
|
strt = ml.bas6.strt
|
# strs = ml.bas6.strt.__getitem__(2)
|
# print(strs.get_value())
|
strt.__setitem__(0,last_array3[0])
|
strt.__setitem__(1,last_array3[1])
|
strt.__setitem__(2,last_array3[2])
|
|
|
mfBase6 = flopy.modflow.ModflowBas(
|
ml,
|
strt= strt,
|
ibound=ml.bas6.ibound,
|
hnoflo=ml.bas6.hnoflo,
|
extension="bas6",)
|
|
mfBase6.write_file(check=False)
|
return True
|
|
|
#修改dis 文件
|
def updateDisFile(model_name, per):
|
|
model_ws = base.model_dir + model_name
|
ml = flopy.modflow.Modflow.load("modflow.nam", model_ws=model_ws,
|
exe_name="mf2005", verbose=True, version="mf2005", check=False)
|
|
mfDis = flopy.modflow.ModflowDis(
|
ml,
|
nlay=ml.dis.nlay,
|
nrow=ml.dis.nrow,
|
ncol=ml.dis.ncol,
|
nper=per,
|
delr=ml.dis.delr,
|
delc=ml.dis.delc,
|
top=ml.dis.top,
|
botm=ml.dis.botm,
|
perlen=ml.dis.perlen,
|
nstp=ml.dis.nstp,
|
tsmult=ml.dis.tsmult,
|
steady=ml.dis.steady,
|
itmuni=ml.dis.itmuni,
|
lenuni=ml.dis.lenuni,
|
extension="dis")
|
|
mfDis.write_file(check=False)
|
|
