Skip to content

Commit

Permalink
Merge pull request #32 from daavid00/dataReport
Browse files Browse the repository at this point in the history 
Horizontal well (sources) for the spe11c cornerp grid
  • Loading branch information
@daavid00
daavid00 authored Feb 4, 2024
2 parents 837d240 + bccbf4f commit e7f25c0
Show file tree
Hide file tree
Showing 3 changed files with 108 additions and 82 deletions.
38 changes: 26 additions & 12 deletions src/pyopmspe11/utils/mapproperties.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -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()
Expand Down Expand Up @@ -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(
Expand All @@ -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()
Expand Down Expand Up @@ -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]
Expand Down Expand Up @@ -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
)
Expand All @@ -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()
Expand All @@ -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

Expand Down
38 changes: 12 additions & 26 deletions src/pyopmspe11/visualization/data.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -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

Expand All @@ -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
Expand Down Expand Up @@ -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:
Expand Down Expand Up @@ -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}, "
Expand All @@ -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}"
Expand Down Expand Up @@ -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
Expand Down Expand Up @@ -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}, "
Expand Down
114 changes: 70 additions & 44 deletions src/pyopmspe11/visualization/plotting.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -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",
Expand All @@ -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":
Expand Down Expand Up @@ -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]")
Expand Down

0 comments on commit e7f25c0

Please sign in to comment.