Horizontal well (sources) for the spe11c cornerp grid

src/pyopmspe11/utils/mapproperties.py

@@ -324,11 +324,12 @@ def corner_point_handling_spe11bc(dic):
         dic (dict): Global dictionary with new added parameters
 
     """
-    well1, well2, sensor1, sensor2, xtemp, centers = [], [], [], [], [], []
+    well1, well2, sensor1, sensor2, xtemp, ztemp, centers = [], [], [], [], [], [], []
     dic["wellijk"] = [[] for _ in range(len(dic["wellCoord"]))]
     for i in range(dic["no_cells"]):
         dic = get_cell_info(dic, i)
         xtemp.append(dic["xyz"][0])
+        ztemp.append(dic["xyz"][2])
         fgl = pd.Series(
             ((dic["cxc1"] - dic["xyz"][0]) ** 2 + (dic["czc1"] - dic["xyz"][2]) ** 2)
         ).argmin()
@@ -371,7 +372,7 @@ def corner_point_handling_spe11bc(dic):
                 for names in ["satnum", "poro", "permx", "thconr"]:
                     dic[f"{names}"].extend(dic[f"{names}"][-dic["noCells"][0] :])
                 for i_i in range(dic["noCells"][0]):
-                    z_c = dic["xyz"][2]
+                    z_c = ztemp[i_i]
                     if dic["spe11"] == "spe11c":
                         z_c -= map_z(dic, j + 1)
                     dic = boxes(
@@ -384,7 +385,7 @@ def corner_point_handling_spe11bc(dic):
                     centers.append(
                         str([xtemp[i_i], dic["ymy_center"][j + 1], z_c])[1:-1]
                     )
-            xtemp = []
+            xtemp, ztemp = [], []
     dic["pop1"] = pd.Series(sensor1).argmin()
     dic["pop2"] = pd.Series(sensor2).argmin()
     dic["well1"] = pd.Series(well1).argmin()
@@ -429,10 +430,23 @@ def locate_wells_sensors(dic):
         dic["wellijk"][1][1] = wji
         dic["wellijkf"][0][1] = wjf
         dic["wellijkf"][1][1] = wjf
-        dic["wellkh"] = [
-            dic["wellijk"][0][2]
-            for _ in range(dic["wellijk"][0][1], dic["wellijkf"][0][1] + 1)
-        ]
+        dic["wellkh"] = []
+        z_centers = []
+        for k in range(dic["noCells"][2]):
+            dic = get_cell_info(dic, well1ijk[0] + k * dic["noCells"][0])
+            z_centers.append(dic["xyz"][2])
+        for j in range(dic["wellijk"][0][1], dic["wellijkf"][0][1] + 1):
+            midpoints = z_centers - map_z(dic, j - 1)
+            dic["wellkh"].append(
+                pd.Series(
+                    abs(
+                        dic["wellCoord"][0][2]
+                        - map_z(dic, dic["wellijk"][0][1] - 1)
+                        - midpoints
+                    )
+                ).argmin()
+                + 1
+            )
     dic["fipnum"][
         dic["sensorijk"][0][0]
         + dic["sensorijk"][0][1] * dic["noCells"][0]
@@ -591,7 +605,7 @@ def wells(dic):
                 )
             if j == 0:
                 well_coord = dic["wellCoord"][j][2]
-                midpoints = dic["zmz_center"] - map_z(dic, dic["wellijk"][j][1])
+                midpoints = dic["zmz_center"] - map_z(dic, dic["wellijk"][j][1] - 1)
                 dic["wellijk"][j].append(
                     pd.Series(abs(well_coord - midpoints)).argmin() + 1
                 )
@@ -604,12 +618,12 @@ def wells(dic):
     if dic["spe11"] == "spe11c":
         dic["wellkh"] = []
         for j in range(dic["wellijk"][0][1], dic["wellijkf"][0][1] + 1):
-            midpoints = dic["zmz_center"] - map_z(dic, j)
+            midpoints = dic["zmz_center"] - map_z(dic, j - 1)
             dic["wellkh"].append(
                 pd.Series(
                     abs(
                         dic["wellCoord"][0][2]
-                        - map_z(dic, dic["wellijk"][0][1])
+                        - map_z(dic, dic["wellijk"][0][1] - 1)
                         - midpoints
                     )
                 ).argmin()
@@ -633,8 +647,8 @@ def map_z(dic, j):
     z_pos = (
         -dic["elevation"]
         + dic["elevation"]
-        * (1.0 - (dic["ymy_center"][j - 1] / (0.5 * dic["dims"][1]) - 1) ** 2.0)
-        - dic["ymy_center"][j - 1] * dic["backElevation"] / dic["dims"][1]
+        * (1.0 - (dic["ymy_center"][j] / (0.5 * dic["dims"][1]) - 1) ** 2.0)
+        - dic["ymy_center"][j] * dic["backElevation"] / dic["dims"][1]
     )
     return z_pos
 

src/pyopmspe11/visualization/data.py

@@ -653,25 +653,14 @@ def get_centers(dig, dil):
     """Centers from the simulation grid"""
     for i in ["x", "y", "z"]:
         dil[f"sim{i}cent"] = [0.0] * dig["nocellst"]
+    with open(f"{dig['path']}/deck/centers.txt", "r", encoding="utf8") as file:
+        for j, row in enumerate(csv.reader(file)):
+            dil["simxcent"][j] = float(row[0])
+            dil["simycent"][j] = float(row[1])
+            dil["simzcent"][j] = dig["dims"][2] - float(row[2])
     if dig["use"] == "opm":
-        # opm.io.ecl.EGrid.xyz_from_ijk is returning nan values, that's why we do this
-        with open(f"{dig['path']}/deck/centers.txt", "r", encoding="utf8") as file:
-            for j, row in enumerate(csv.reader(file)):
-                dil["simxcent"][j] = float(row[0])
-                dil["simycent"][j] = float(row[1])
-                dil["simzcent"][j] = dig["dims"][2] - float(row[2])
         dil["satnum"] = list(dig["init"]["SATNUM"])
     else:
-        for cell in dig["egrid"].cells():
-            (
-                dil["simxcent"][cell.global_index],
-                dil["simycent"][cell.global_index],
-                dil["simzcent"][cell.global_index],
-            ) = (
-                cell.coordinate[0],
-                cell.coordinate[1],
-                dig["dims"][2] - cell.coordinate[2],
-            )
         dil["satnum"] = list(dig["init"].iget_kw("SATNUM")[0])
     return dil
 
@@ -696,7 +685,7 @@ def dense_data(dig):
                 dil["refzgrid"][ind] = k
                 ind += 1
     ind, dil["cell_ind"] = 0, np.zeros(dig["nocellst"], dtype=int)
-    for i, j, k in zip(dil["simxcent"], dil["simycent"], np.flip(dil["simzcent"])):
+    for i, j, k in zip(dil["simxcent"], dil["simycent"], dil["simzcent"]):
         dil["cell_ind"][ind] = pd.Series(
             (dil["refxgrid"] - i) ** 2
             + (dil["refygrid"] - j) ** 2
@@ -747,8 +736,9 @@ def static_map_to_report_grid_performance_spatial(dig, dil):
             if j > 0:
                 infotimes.append(86400.0 * float((row[0].strip()).split()[0]))
                 tsteps.append(86400.0 * float((row[0].strip()).split()[1]))
+    infotimes = np.array(infotimes)
     for time in dig["dense_t"][:-1]:
-        ind = infotimes.index(time + dig["time_initial"])
+        ind = pd.Series(np.abs(infotimes - (time + dig["time_initial"]))).argmin()
         if ind > 0:
             dil["latest_dts"].append(tsteps[ind - 1])
         else:
@@ -846,9 +836,7 @@ def write_dense_data_performance_spatial(dig, dil, i):
         idxy = 0
         for ycord in dil["refycent"]:
             for xcord in dil["refxcent"]:
-                idc = (
-                    dig["nxyz"][0] * dig["nxyz"][1] * (dig["nxyz"][2] - idz - 1) + idxy
-                )
+                idc = -dig["nxyz"][0] * dig["nxyz"][1] * (dig["nxyz"][2] - idz) + idxy
                 if dig["case"] != "spe11c":
                     text.append(
                         f"{xcord:.3e}, {zcord:.3e}, "
@@ -860,7 +848,7 @@ def write_dense_data_performance_spatial(dig, dil, i):
                 else:
                     text.append(
                         f"{xcord:.3e}, {ycord:.3e}, {zcord:.3e}, "
-                        + f"{dil['cvol_refg'][idc] :.3e}, {dil['arrat_refg'][idc] :.3e}, "
+                        + f"{dil['cvol_refg'][idc] :.3e}, {dil['arat_refg'][idc] :.3e}, "
                         + f"{dil['co2mn_refg'][idc] :.3e}, {dil['h2omn_refg'][idc] :.3e}, "
                         + f"{dil['co2mb_refg'][idc] :.3e}, {dil['h2omb_refg'][idc] :.3e}, "
                         + f"{np.nan}"
@@ -947,7 +935,7 @@ def handle_inactive_mapping(dig, dil):
             dil["refygrid"][i] = np.inf
             dil["refzgrid"][i] = np.inf
     ind, dil["cell_ind"] = 0, np.zeros(dig["nocellst"], dtype=int)
-    for i, j, k in zip(dil["simxcent"], dil["simycent"], np.flip(dil["simzcent"])):
+    for i, j, k in zip(dil["simxcent"], dil["simycent"], dil["simzcent"]):
         dil["cell_ind"][ind] = pd.Series(
             (dil["refxgrid"] - i) ** 2
             + (dil["refygrid"] - j) ** 2
@@ -1005,9 +993,7 @@ def write_dense_data(dig, dil, i):
         idxy = 0
         for ycord in dil["refycent"]:
             for xcord in dil["refxcent"]:
-                idc = (
-                    dig["nxyz"][0] * dig["nxyz"][1] * (dig["nxyz"][2] - idz - 1) + idxy
-                )
+                idc = -dig["nxyz"][0] * dig["nxyz"][1] * (dig["nxyz"][2] - idz) + idxy
                 if dig["case"] == "spe11a":
                     text.append(
                         f"{xcord:.3e}, {zcord:.3e}, {dil['pressure_refg'][idc] :.3e}, "

src/pyopmspe11/visualization/plotting.py

@@ -316,48 +316,61 @@ def generate_grid(dic):
     return dic
 
 
+def handle_kind(dic, kind):
+    """Identify between dense and performance-spatial"""
+    if kind == "":
+        dic["quantities"] = [
+            "pressure",
+            "sgas",
+            "xco2",
+            "xh20",
+            "gden",
+            "wden",
+            "tco2",
+        ]
+        dic["units"] = [
+            "[Pa]",
+            "[-]",
+            "[-]",
+            "[-]",
+            r"[kg/m$^3$]",
+            r"[kg/m$^3$]",
+            "[kg]",
+        ]
+        if dic["case"] != "spe11a":
+            dic["quantities"] += ["temp"]
+            dic["units"] += ["C"]
+    else:
+        dic["quantities"] = [
+            "cvol",
+            "arat",
+            "CO2 max_norm_res",
+            "H2O max_norm_res",
+            "CO2 mb_error",
+            "H2O mb_error",
+        ]
+        dic["units"] = [r"[m$^3$]", "[-]", "[-]", "[-]", "[-]", "[-]"]
+    return dic
+
+
+def ini_quantity_plot(dic):
+    """Initialize the figure"""
+    if dic["case"] != "spe11a":
+        dic["fig"] = plt.figure(figsize=(50, 3 * len(dic["times"])))
+    else:
+        dic["fig"] = plt.figure(figsize=(45, 6.5 * len(dic["times"])))
+    for name in ["plot", "min", "max", "sum"]:
+        dic[f"{name}"] = []
+    return dic
+
+
 def dense_data(dic):
     """2D spatial maps"""
     dic = generate_grid(dic)
     for kind in dic["kinds"]:
-        if kind == "":
-            dic["quantities"] = [
-                "pressure",
-                "sgas",
-                "xco2",
-                "xh20",
-                "gden",
-                "wden",
-                "tco2",
-            ]
-            dic["units"] = [
-                "[Pa]",
-                "[-]",
-                "[-]",
-                "[-]",
-                r"[kg/m$^3$]",
-                r"[kg/m$^3$]",
-                "[kg]",
-            ]
-            if dic["case"] != "spe11a":
-                dic["quantities"] += ["temp"]
-                dic["units"] += ["C"]
-        else:
-            dic["quantities"] = [
-                "cvol",
-                "arat",
-                "CO2 max_norm_res",
-                "H2O max_norm_res",
-                "CO2 mb_error",
-                "H2O mb_error",
-            ]
-            dic["units"] = [r"[m$^3$]", "[-]", "[-]", "[-]", "[-]", "[-]"]
+        dic = handle_kind(dic, kind)
         for k, quantity in enumerate(dic["quantities"]):
-            if dic["case"] != "spe11a":
-                dic["fig"] = plt.figure(figsize=(50, 3 * len(dic["times"])))
-            else:
-                dic["fig"] = plt.figure(figsize=(45, 6.5 * len(dic["times"])))
-            dic["plot"] = []
+            dic = ini_quantity_plot(dic)
             csv = np.genfromtxt(
                 f"{dic['exe']}/{dic['folders'][0]}/data/{dic['case']}{kind}_spatial_map_"
                 + f"0{dic['tlabel']}.csv",
@@ -377,8 +390,12 @@ def dense_data(dic):
                     skip_header=1,
                 )
                 quan = np.array([csv[i][dic["dims"] + k] for i in range(csv.shape[0])])
-                dic["minc"] = min(dic["minc"], quan[quan >= 0].min())
-                dic["maxc"] = max(dic["maxc"], quan[quan >= 0].max())
+                dic["min"].append(quan[quan >= 0].min())
+                dic["max"].append(quan[quan >= 0].max())
+                if quantity == "tco2":
+                    dic["sum"].append(quan[quan >= 0].sum())
+                dic["minc"] = min(dic["minc"], dic["min"][-1])
+                dic["maxc"] = max(dic["maxc"], dic["max"][-1])
                 dic["plot"].append(np.zeros([len(dic["zmz"]) - 1, len(dic["xmx"]) - 1]))
                 for i in np.arange(0, len(dic["zmz"]) - 1):
                     if dic["case"] != "spe11c":
@@ -406,11 +423,20 @@ def dense_data(dic):
                     shading="flat",
                     cmap=dic["cmaps"][k],
                 )
-                axis.set_title(
-                    f"{time}{dic['tlabel']}, {quantity} "
-                    + dic["units"][k]
-                    + f", {dic['case']} ({dic['folders'][0]})"
-                )
+                if quantity == "tco2":
+                    axis.set_title(
+                        f"{time}{dic['tlabel']}, {quantity} "
+                        + dic["units"][k]
+                        + f"(sum={dic['sum'][j]:.1E})"
+                        + f", {dic['case']} ({dic['folders'][0]})"
+                    )
+                else:
+                    axis.set_title(
+                        f"{time}{dic['tlabel']}, {quantity} "
+                        + dic["units"][k]
+                        + f"(min={dic['min'][j]:.1E}, max={dic['max'][j]:.1E})"
+                        + f", {dic['case']} ({dic['folders'][0]})"
+                    )
                 axis.axis("scaled")
                 axis.set_xlabel("x [m]")
                 axis.set_ylabel("z [m]")