devbox
Cpp五条
这个屌丝很懒，什么也没留下！
关注作者
热门标签
热门文章
当前位置:   article > 正文

集总式新安江代码（二、三水源)python_新安江模型python代码

作者：Cpp五条 | 2024-04-12 18:01:29

赞

踩

新安江模型python代码

首先设置模型参数

  1. import numpy as np
  2. import pandas as pd
  3. import math
  4. df = pd.read_excel(r"C:\Users\Admin\Desktop\pet.xlsx")
  5. df = df.fillna(0)
  6. print(df)
  7. WM,WUM,WLM,WDM,b,K,C,FC,F,T,CG,L,CR = 120,15,85,20,0.3,0.95,0.14,15,2856,24,0.9,5,0.5    # F单位为km^2 T单位为h
  8. Sm,KI,KG,EX,FR1,CI,Q = 20,0.2,0.5,1.5,0.4175,0.9,5.6
  9. WMM = (1+b)*WM
  10. U = F / (3.6 * T)
  11. SMM = Sm * (1+EX)

我的数据是excel数据其中有流域面平均雨量和蒸发资料

首先是设置一个for循坏我一共5000多行降雨数据，设置一个for循坏

  1. for i in range(0,5479):
  2. # 计算土湿
  3.     if i == 0:
  4.         df.loc[i, "WU"] = 0
  5.         df.loc[i,"WL"] = 2.2
  6.         df.loc[i, "WD"] = 20
  7.         df.loc[i, "EP"] = K * df.loc[i, "E0"]
  8.         if df.loc[i, "WU"] + df.loc[i, "P"] >= df.loc[i, "EP"]:
  9.             df.loc[i, "EU"] = df.loc[i, "EP"]
  10.             df.loc[i, "EL"] = 0
  11.             df.loc[i, "ED"] = 0
  12.         if df.loc[i, "WU"] + df.loc[i, "P"] < df.loc[i, "EP"] and df.loc[i, "WL"] >= C * WLM:
  13.             df.loc[i, "EU"] = df.loc[i, "WU"] + df.loc[i, "P"]
  14.             df.loc[i, "EL"] = (df.loc[i, "EP"] - df.loc[i, "EU"]) * df.loc[i, "WL"] / WLM
  15.             df.loc[i, "ED"] = 0
  16.         if df.loc[i, "WU"] + df.loc[i, "P"] < df.loc[i, "EP"] and df.loc[i, "WL"] < C * WLM and df.loc[i, "WL"] >= C * (df.loc[i, "EP"] - df.loc[i, "EU"]):
  17.             df.loc[i, "EU"] = df.loc[i, "WU"] + df.loc[i, "P"]
  18.             df.loc[i, "EL"] = C * (df.loc[i, "EP"] - df.loc[i, "EU"])
  19.             df.loc[i, "ED"] = 0
  20.         if df.loc[i, "WU"] + df.loc[i, "P"] < df.loc[i, "EP"] and df.loc[i, "WL"] < C * (
  21.                 df.loc[i, "EP"] - df.loc[i, "EU"]):
  22.             df.loc[i, "EU"] = df.loc[i, "WU"] + df.loc[i, "P"]
  23.             df.loc[i, "EL"] = df.loc[i, "WL"]
  24.             df.loc[i, "ED"] = C * (df.loc[i, "EP"] - df.loc[i, "EU"]) - df.loc[i, "EL"]
  25.         df.loc[i, "E"] = df.loc[i, "EU"] + df.loc[i, "EL"] + df.loc[i, "ED"]
  26.         df.loc[i, "PE"] = df.loc[i, "P"] - df.loc[i, "E"]
  27.         df.loc[i, "W"] = df.loc[i, "WU"] + df.loc[i, "WL"] + df.loc[i, "WD"]
  28.     else:
  29.         dw = df.loc[i-1,"PE"] - df.loc[i-1,"R"] #这边需要连接下面R，否则会报错
  30.         df.loc[i,"WU"] = df.loc[i-1,"WU"] + dw
  31.         df.loc[i,"WL"] = df.loc[i-1,"WL"]
  32.         df.loc[i,"WD"] = df.loc[i-1,"WD"]
  33.         if df.loc[i,"WU"] < 0 :
  34.             df.loc[i,"WL"] = df.loc[i-1,"WL"] + dw + df.loc[i-1,"WU"]
  35.             df.loc[i,"WD"] = df.loc[i-1,"WD"]
  36.             df.loc[i,"WU"] = 0
  37.             if df.loc[i,"WL"] < 0:
  38.                 df.loc[i,"WD"] = df.loc[i-1,"WD"] + dw + df.loc[i-1,"WU"] + df.loc[i-1,"WL"]
  39.                 df.loc[i,"WU"] = 0
  40.                 df.loc[i,"WL"] = 0
  41.                 if df.loc[i,"WD"] < 0:
  42.                     df.loc[i,"WU"] = 0
  43.                     df.loc[i,"WD"] = 0
  44.                     df.loc[i,"WL"] = 0
  45.         if df.loc[i,"WU"] > WUM:
  46.             df.loc[i,"WU"] = WUM
  47.             df.loc[i,"WL"] = df.loc[i-1,"WL"] + dw - WUM + df.loc[i-1,"WU"]
  48.             df.loc[i,"WD"] = df.loc[i-1,"WD"]
  49.             if df.loc[i,"WL"] > WLM:
  50.                 df.loc[i,"WU"] = WUM
  51.                 df.loc[i,"WL"] = WLM
  52.                 df.loc[i,"WD"] = df.loc[i-1,"WD"] + dw - WLM - WUM + df.loc[i-1,"WU"] + df.loc[i-1,"WL"]
  53.         df.loc[i,"EP"] = K * df.loc[i,"E0"]
  54.         if df.loc[i,"WU"] + df.loc[i,"P"] >= df.loc[i,"EP"]:
  55.             df.loc[i,"EU"] = df.loc[i,"EP"]
  56.             df.loc[i,"EL"] = 0
  57.             df.loc[i,"ED"] = 0
  58.         if df.loc[i,"WU"] + df.loc[i,"P"] < df.loc[i,"EP"] and df.loc[i,"WL"] >= C * WLM:
  59.             df.loc[i,"EU"] = df.loc[i,"WU"] + df.loc[i,"P"]
  60.             df.loc[i,"EL"] = (df.loc[i,"EP"] - df.loc[i,"EU"]) * df.loc[i,"WL"] / WLM
  61.             df.loc[i,"ED"] = 0
  62.         if df.loc[i,"WU"] + df.loc[i,"P"] < df.loc[i,"EP"] and df.loc[i,"WL"] < C * WLM and df.loc[i,"WL"] >= C * (df.loc[i,"EP"]-df.loc[i,"EU"]):
  63.             df.loc[i,"EU"] = df.loc[i,"WU"] + df.loc[i,"P"]
  64.             df.loc[i,"EL"] = C * (df.loc[i,"EP"] - df.loc[i,"EU"])
  65.             df.loc[i,"ED"] = 0
  66.         if df.loc[i,"WU"] + df.loc[i,"P"] < df.loc[i,"EP"] and df.loc[i,"WL"] < C * (df.loc[i,"EP"]-df.loc[i,"EU"]):
  67.             df.loc[i,"EU"] = df.loc[i,"WU"] + df.loc[i,"P"]
  68.             df.loc[i,"EL"] = df.loc[i,"WL"]
  69.             df.loc[i,"ED"] = C * (df.loc[i,"EP"] - df.loc[i,"EU"]) - df.loc[i,"EL"]
  70.         df.loc[i,"E"] = df.loc[i,"EU"] + df.loc[i,"EL"] + df.loc[i,"ED"]
  71.         df.loc[i,"PE"] = df.loc[i,"P"] - df.loc[i,"E"]
  72.         df.loc[i,"W"] = df.loc[i,"WU"] + df.loc[i,"WL"] + df.loc[i,"WD"]

再根据自己的要求计算R

  1. # 计算R
  2.     a = WMM * (1- math.pow((1-df.loc[i,"W"]/WM),1/(1+b)))
  3.     if i ==0:
  4.         if df.loc[i,"PE"] > 0:
  5.             if a + df.loc[i,"PE"] <= WMM:
  6.                 df.loc[i,"R"] = df.loc[i,"PE"] + df.loc[i,"W"] - WM + WM * math.pow(1-((df.loc[i,"PE"]+a)/WMM),b+1)  # df["R"][i] = df["PE"][i] + df["W"][i] - WM + WM * math.pow(1 - ((df["PE"][i] + a) / WMM), b + 1)
  7.             if a + df.loc[i,"PE"] > WMM:
  8.                 df.loc[i,"R"] = df.loc[i,"PE"] - (WM - df.loc[i,"W"])
  9.         if df.loc[i,"R"] < 0:
  10.             df.loc[i,"R"] = 0
  11.     else:
  12.         if df.loc[i,"PE"] > 0:
  13.             if a + df.loc[i, "PE"] <= WMM:
  14.                  df.loc[i, "R"] = df.loc[i, "PE"] + df.loc[i, "W"] - WM + WM * math.pow(1 - ((df.loc[i, "PE"] + a) / WMM),b + 1)  # df["R"][i] = df["PE"][i] + df["W"][i] - WM + WM * math.pow(1 - ((df["PE"][i] + a) / WMM), b + 1)
  15.             if a + df.loc[i, "PE"] > WMM:
  16.                 df.loc[i, "R"] = df.loc[i, "PE"] - (WM - df.loc[i, "W"])
  17.         if df.loc[i, "R"] < 0:
  18.             df.loc[i, "R"] = 0

在计算产流面积

  1. # 计算流域产流面积
  2.     if i == 0:
  3.         if df.loc[i,"R"] > 0 :
  4.             df.loc[i,"FR"] = df.loc[i,"R"] / df.loc[i,"PE"]
  5.             if df.loc[i,"FR"] >1:
  6.                 df.loc[i,"FR"] =1
  7.  
  8.         else:
  9.             df.loc[i,"FR"] = FR1
  10.     else:
  11.         if df.loc[i, "R"] > 0:
  12.             df.loc[i, "FR"] = df.loc[i, "R"] / df.loc[i, "PE"]
  13.             if df.loc[i, "FR"] > 1:
  14.                 df.loc[i, "FR"] = 1
  15.         else:
  16.             df.loc[i, "FR"] = df.loc[i-1,"FR"]

这边是二水源代码

  1. # 划分二水源
  2. #     if i == 0:
  3. #         if df.loc[i, "R"] > 0:
  4. #             df.loc[i, "FR"] = df.loc[i, "R"] / df.loc[i, "PE"]
  5. #             if df.loc[i,"FR"] > 1:
  6. #                 df.loc[i,"FR"] = 1
  7. #             if df.loc[i, "FR"] > 0:
  8. #                 if df.loc[i, "PE"] > FC:
  9. #                     df.loc[i, "RS"] = (df.loc[i, "PE"] - FC) * df.loc[i, "FR"]
  10. #                     df.loc[i, "RG"] = df.loc[i, "R"] - df.loc[i, "RS"]
  11. #                 if df.loc[i, "PE"] <= FC:
  12. #                     df.loc[i, "RS"] = 0
  13. #                     df.loc[i, "RG"] = df.loc[i, "R"]
  14. #             if df.loc[i, "FR"] < 0:
  15. #                 df.loc[i, "FR"] = 0
  16. #                 df.loc[i, "RS"] = 0
  17. #                 df.loc[i, "RG"] = 0
  18. #         if df.loc[i, "R"] == 0:
  19. #             df.loc[i, "RS"] = 0
  20. #             df.loc[i, "RG"] = 0
  21. #     else:
  22. #         if df.loc[i,"R"] > 0 :
  23. #             df.loc[i,"FR"] = df.loc[i,"R"] / df.loc[i,"PE"]
  24. #             if df.loc[i,"FR"] > 0:
  25. #                 if df.loc[i,"PE"] > FC:
  26. #                     df.loc[i,"RS"] = (df.loc[i,"PE"] - FC) * df.loc[i,"FR"]
  27. #                     df.loc[i,"RG"] = df.loc[i,"R"] - df.loc[i,"RS"]
  28. #                 if df.loc[i,"PE"] <= FC:
  29. #                     df.loc[i,"RS"] = 0
  30. #                     df.loc[i,"RG"] = df.loc[i,"R"]
  31. #             if df.loc[i,"FR"] < 0:
  32. #                 df.loc[i,"FR"] = 0
  33. #                 df.loc[i,"RS"] = 0
  34. #                 df.loc[i,"RG"] = 0
  35. #         if df.loc[i,"R"] == 0:
  36. #             df.loc[i,"RS"] = 0
  37. #             df.loc[i,"RG"] = 0

三水源

  1. # 划分三水源
  2.     if i == 0:     # 第一时段计算
  3.         df.loc[i, "S1"] = 3.4525
  4.         if df.loc[i,"PE"] > 0:
  5.             AU = SMM * (1 - math.pow((1-(((df.loc[i,"S1"]*FR1)/df.loc[i,"FR"])/Sm)),1/(1+EX)))
  6.             if df.loc[i,"PE"] + AU < SMM:
  7.                 df.loc[i,"RS"] = df.loc[i,"FR"] * (df.loc[i,"PE"] +(df.loc[i,"S1"] * FR1)/df.loc[i,"FR"] - Sm +
  8.                                                  Sm *math.pow((1-(df.loc[i,"PE"] +AU)/SMM),1+EX))
  9.             if df.loc[i,"PE"] + AU >= SMM:
  10.                 df.loc[i,"RS"] = df.loc[i,"FR"] * (df.loc[i,"PE"] + (df.loc[i,"S1"] * FR1)/df.loc[i,"FR"] - Sm)
  11.             S = (df.loc[i,"S1"] * FR1)/df.loc[i,"FR"] + (df.loc[i,"R"] - df.loc[i,"RS"])/df.loc[i,"FR"]
  12.             df.loc[i,"RI"] = KI * S *df.loc[i,"FR"]
  13.             df.loc[i,"RG"] = KG * S *df.loc[i,"FR"]
  14.             df.loc[i+1,"S1"] = S*(1-KI-KG)
  15.         else:
  16.             S = (df.loc[i, "S1"] * FR1) / df.loc[i, "FR"]
  17.             df.loc[i + 1, "S1"] = S * (1 - KG - KI)
  18.             df.loc[i, "RS"] = 0
  19.             df.loc[i, "RG"] = KI*S * df.loc[i,"FR"]
  20.             df.loc[i, "RI"] = KG * S *df.loc[i,"FR"]
  21.     else: #其余时段计算
  22.         if df.loc[i, "PE"] > 0:
  23.             AU = SMM * (1 - math.pow((1-(((df.loc[i,"S1"]*df.loc[i-1,"FR"])/df.loc[i,"FR"])/Sm)),1/(1+EX)))
  24.             if df.loc[i,"PE"] + AU < SMM:
  25.                 df.loc[i,"RS"] = df.loc[i,"FR"] * (df.loc[i,"PE"] +(df.loc[i,"S1"] * df.loc[i-1,"FR"])/df.loc[i,"FR"] - Sm +
  26.                                                  Sm *math.pow((1-(df.loc[i,"PE"] +AU)/SMM),1+EX))
  27.             if df.loc[i,"PE"] + AU >= SMM:
  28.                 df.loc[i,"RS"] = df.loc[i,"FR"] * (df.loc[i,"PE"] + (df.loc[i,"S1"] * df.loc[i-1,"FR"])/df.loc[i,"FR"] - Sm)
  29.             S = (df.loc[i,"S1"] * df.loc[i-1,"FR"])/df.loc[i,"FR"] + (df.loc[i,"R"] - df.loc[i,"RS"])/df.loc[i,"FR"]
  30.             df.loc[i,"RI"] = KI * S *df.loc[i,"FR"]
  31.             df.loc[i,"RG"] = KG * S *df.loc[i,"FR"]
  32.             df.loc[i+1,"S1"] = S*(1-KI-KG)
  33.         else:
  34.             S = (df.loc[i, "S1"] * df.loc[i-1,"FR"]) / df.loc[i, "FR"]
  35.             df.loc[i + 1, "S1"] =S * (1 - KG - KI)
  36.             df.loc[i, "RS"] = 0
  37.             df.loc[i, "RG"] = KG *S *df.loc[i,"FR"]
  38.             df.loc[i, "RI"] = KI * S *df.loc[i,"FR"]

二水源汇流

  1. # # 计算二水源汇流
  2. #     if i == 0:
  3. #         df.loc[i, "Qs"] = df.loc[i, "RS"] * U
  4. #         df.loc[i,"Qg"] = 5.6
  5. #         df.loc[i, "QT"] = df.loc[i, "Qs"] + df.loc[i, "Qg"]
  6. #         df.loc[i,"Qt"] = 5.6
  7. #     else:
  8. #         df.loc[i,"Qs"] = df.loc[i,"RS"] * U
  9. #         df.loc[i,"Qg"] = CG * df.loc[i-1, "Qg"] + (1-CG) * df.loc[i,"RG"] * U
  10. #         df.loc[i,"QT"] = df.loc[i,"Qs"] + df.loc[i,"Qg"]
  11. #         df.loc[i,"Qt"] = CR * df.loc[i-1,"Qt"] + (1 - CR) * df.loc[i-L,"QT"]

三水源汇流

  1. # 计算三水源汇流
  2.     if i == 0:
  3.         df.loc[i,"QS"] = df.loc[i,"RS"] * U
  4.         df.loc[i,"QI"] = 1/3 * Q
  5.         df.loc[i,"QG"] = 1/3 * Q
  6.         df.loc[i, "QT"] = df.loc[i, "QS"] + df.loc[i, "QI"] + df.loc[i, "QG"]
  7.     if i >0:
  8.         df.loc[i, "QS"] = df.loc[i, "RS"] * U
  9.         df.loc[i,"QI"] = CI * df.loc[i-1,"QI"] + (1-CI)*df.loc[i,"RI"]*U
  10.         df.loc[i,"QG"] = CG* df.loc[i-1,"QG"]+(1-CG)*df.loc[i,"RG"] *U
  11.         df.loc[i,"QT"] = df.loc[i,"QS"]+df.loc[i,"QI"]+df.loc[i,"QG"]
  12.     if i >=0 and i <= L:
  13.         df.loc[i,"Qt"] = Q
  14.     if i >= L+1:
  15.         df.loc[i,"Qt"] = df.loc[i-1,"Qt"]*CR+(1-CR) * df.loc[i-L,"QT"]

最后写入excel表，保留4位小数

  1. ds = pd.DataFrame(df)
  2. ds.to_excel(r"C:\Users\Admin\Desktop\新安江1.xlsx",float_format='%.4f')

最后三水源结果

写到这感觉自己写的太麻烦了，希望以后能够改进吧

这边蒸散发还是借鉴了一个哥们的代码，链接在这边：(26条消息) Python数据分析实例，利用Pandas建立流域三层蒸发和蓄满产流模型_三层蒸发模型代码_Joyonlyonly的博客-CSDN博客

声明：本文内容由网友自发贡献，不代表【wpsshop博客】立场，版权归原作者所有，本站不承担相应法律责任。如您发现有侵权的内容，请联系我们。转载请注明出处：https://www.wpsshop.cn/w/Cpp五条/article/detail/411961
推荐阅读
相关标签

Copyright © 2003-2013 www.wpsshop.cn 版权所有，并保留所有权利。

  

闽ICP备14008679号