Allowing for multidimensional doublet_calculations (a2fcab6f) · Commits · AGS_public / pythermogis

src/statistics/calculate_pvalues.py

+0 −1

Original line number	Diff line number	Diff line
		@@ -16,5 +16,4 @@ def generate_thickness_permeability_transmissivity_for_pvalue(thickness_mean, th
		transmissivity_samples = np.sort(np.exp(ln_permeability_dist.rvs(nSamples) + np.log(thickness_dist.rvs(nSamples))))
		transmissivity_pvalue_sampled = transmissivity_samples[int((1 - Pvalue) * nSamples)]


		return thickness_pvalue, permeability_pvalue, transmissivity_pvalue_sampled

src/thermogis_classes/doublet.py

+125 −5

Original line number	Diff line number	Diff line
		import numpy as np
		import xarray as xr

		from src.physics.temperature_grid_calculation import calculate_temperature_from_gradient
		from src.statistics.calculate_pvalues import generate_thickness_permeability_transmissivity_for_pvalue

		def calculate_performance(hydrocarbons: float, depth: float, thickness: float, porosity: float, ntg: float, temperature: float, transmissivity: float, transmissivity_with_ntg: float, doublet=None,

		def calculate_performance_of_single_location(hydrocarbons: float, depth: float, thickness: float, porosity: float, ntg: float, temperature: float, transmissivity: float, transmissivity_with_ntg: float, doublet=None,
		input_params: dict = None):
		# Returns the values from a doublet simulation in order of;
		# power, heat_pump_power, capex, opex, utc, npv, hprod, cop, cophp, pres, flow_rate, welld
		@@ -46,7 +50,124 @@ def calculate_performance(hydrocarbons: float, depth: float, thickness: float, p

		# Reset doublet variables for next calculation
		doublet.setProjectVariables(False, 0.0)
		return output_values["power"], output_values["heat_pump_power"], output_values["capex"], output_values["opex"], output_values["utc"], output_values["npv"], output_values["hprod"], output_values["cop"], output_values["cophp"], output_values["pres"], output_values["flow_rate"], output_values["welld"]
		return output_values["power"], output_values["heat_pump_power"], output_values["capex"], output_values["opex"], output_values["utc"], output_values["npv"], output_values["hprod"], output_values["cop"], output_values[
		"cophp"], output_values["pres"], output_values["flow_rate"], output_values["welld"]


		def calculate_performance_across_dimensions(input_data: xr.Dataset,
		doublet,
		input_params=None,
		p_values=None):
		"""
		Given a set of input_parameters, and a dataset which contains the variables: thickness, thickness_sd, porosity,
		:param input_params:
		:param input_data:
		:return:
		"""
		if p_values is None:
		p_values = [50.0]
		if input_params is None:
		input_params = {"hp_minimum_injection_temperature": 15,
		"return_temperature": 30,
		"surface_temperature": 10,
		"degrees_per_km": 31,
		"max_cooling_temperature_range": 100,
		"stimKhMax": 20,
		"use_stimulation": False,
		"use_heat_pump": False,
		"utc_cutoff_shallow": 5.1,
		"utc_cutoff_deep": 6.5,
		"utc_cutoff_depth": 4000.0}

		if "temperature" not in input_data:
		input_data["temperature"] = calculate_temperature_from_gradient(input_data.depth, input_data.thickness_mean, input_params["degrees_per_km"], input_params["surface_temperature"])

		# Setup output_data
		output_data = input_data.thickness_mean.copy().to_dataset(name="thickness")
		output_data.thickness.data[:] = np.nan
		output_data = output_data.expand_dims({"p_value": p_values})

		# Calculate Thickness, Permeability and Transmissivity for each P-value
		thickness_data = []
		permeability_data = []
		transmissivity_data = []
		for i, p in enumerate(p_values):
		thickness, permeability, transmissivity = xr.apply_ufunc(generate_thickness_permeability_transmissivity_for_pvalue,
		input_data.thickness_mean,
		input_data.thickness_sd,
		input_data.ln_permeability_mean,
		input_data.ln_permeability_sd,
		p,
		input_core_dims=[[], [], [], [], []],
		output_core_dims=[[], [], []],
		vectorize=True,
		)
		thickness_data.append(thickness.data)
		permeability_data.append(permeability.data)
		transmissivity_data.append(transmissivity.data)

		output_data["thickness"] = (output_data.coords, thickness_data)
		output_data["permeability"] = (output_data.coords, permeability_data)
		output_data["transmissivity"] = (output_data.coords, transmissivity_data)
		output_data[f"transmissivity_with_ntg"] = (output_data[f"transmissivity"] * input_data.ntg) / 1e3

		# Calculate performance for each P-value
		power_data = []
		heat_pump_power_data = []
		capex_data = []
		opex_data = []
		utc_data = []
		npv_data = []
		hprod_data = []
		cop_data = []
		cophp_data = []
		pres_data = []
		flow_rate_data = []
		welld_data = []
		for i, p in enumerate(p_values):
		output_data_arrays = xr.apply_ufunc(calculate_performance_of_single_location,
		input_data.hc_accum,
		input_data.depth,
		output_data.thickness.isel(p_value=i),
		input_data.porosity,
		input_data.ntg,
		input_data.temperature,
		output_data.transmissivity.isel(p_value=i),
		output_data.transmissivity_with_ntg.isel(p_value=i),
		kwargs={"doublet": doublet, "input_params": input_params},
		input_core_dims=[[], [], [], [], [], [], [], []],
		output_core_dims=[[], [], [], [], [], [], [], [], [], [], [], []],
		vectorize=True,
		)

		# Assign values from calculate performance to their grids for each p-value
		power_data.append(output_data_arrays[0])
		heat_pump_power_data.append(output_data_arrays[1])
		capex_data.append(output_data_arrays[2])
		opex_data.append(output_data_arrays[3])
		utc_data.append(output_data_arrays[4])
		npv_data.append(output_data_arrays[5])
		hprod_data.append(output_data_arrays[6])
		cop_data.append(output_data_arrays[7])
		cophp_data.append(output_data_arrays[8])
		pres_data.append(output_data_arrays[9])
		flow_rate_data.append(output_data_arrays[10])
		welld_data.append(output_data_arrays[11])

		output_data["power"] = (output_data.coords, power_data)
		output_data["heat_pump_power"] = (output_data.coords, heat_pump_power_data)
		output_data["capex"] = (output_data.coords, capex_data)
		output_data["opex"] = (output_data.coords, opex_data)
		output_data["utc"] = (output_data.coords, utc_data)
		output_data["npv"] = (output_data.coords, npv_data)
		output_data["hprod"] = (output_data.coords, hprod_data)
		output_data["cop"] = (output_data.coords, cop_data)
		output_data["cophp"] = (output_data.coords, cophp_data)
		output_data["pres"] = (output_data.coords, pres_data)
		output_data["flow_rate"] = (output_data.coords, flow_rate_data)
		output_data["welld"] = (output_data.coords, welld_data)

		return output_data


		def set_doublet_parameters(doublet, transmissivity_with_ntg, depth, porosity, ntg, temperature, useStimulation, stimKhMax, surface_temperature, return_temperature, use_heat_pump, max_cooling_temperature_range,
		@@ -68,4 +189,3 @@ def set_doublet_parameters(doublet, transmissivity_with_ntg, depth, porosity, nt

		doublet.doubletCalc1DData.setInjectionTemp(injectionTemp)
		doublet.doubletCalc1DData.setDhReturnTemp(return_temperature)

tests/tests/test_doublet.py

+28 −2

Original line number	Diff line number	Diff line
		@@ -20,7 +20,7 @@ from pygridsio.grid_operations import resample_xarray_grid

		from src.physics.temperature_grid_calculation import calculate_temperature_from_gradient
		from src.statistics.calculate_pvalues import generate_thickness_permeability_transmissivity_for_pvalue
		from src.thermogis_classes.doublet import set_doublet_parameters, calculate_performance
		from src.thermogis_classes.doublet import set_doublet_parameters, calculate_performance_of_single_location, calculate_performance_across_dimensions
		from src.utils.get_repo_root import get_repo_root


		@@ -141,7 +141,7 @@ class PyThermoGIS(TestCase):
		output_grids[f"transmissivity_with_ntg_P{P_value}"] = (output_grids[f"transmissivity_P{P_value}"] * input_grids.ntg) / 1e3

		# Calculate performance for the provided P-value
		output_data_arrays = xr.apply_ufunc(calculate_performance,
		output_data_arrays = xr.apply_ufunc(calculate_performance_of_single_location,
		input_grids.hc_accum,
		input_grids.depth,
		output_grids[f"thickness_P{P_value}"],
		@@ -184,3 +184,29 @@ class PyThermoGIS(TestCase):
		xr.testing.assert_allclose(output_grids.pres, read_grid(self.test_benchmark_output_simplified / f"simplified__pres_P{P_value}.nc"), atol=1.0)
		xr.testing.assert_allclose(output_grids.flow_rate, read_grid(self.test_benchmark_output_simplified / f"simplified__flowr_P{P_value}.nc"), atol=3.0)
		xr.testing.assert_allclose(output_grids.welld, read_grid(self.test_benchmark_output_simplified / f"simplified__welld_P{P_value}.nc"), atol=20)

		def test_doublet_grid_multi_dimension(self):
		# Arrange
		# Import Java Classes
		Logger = JClass("logging.Logger")
		Mockito = JClass("org.mockito.Mockito")
		RNG = JClass("tno.geoenergy.stochastic.RandomNumberGenerator")
		ThermoGISDoublet = JClass("thermogis.calc.doublet.ThermoGisDoublet")
		UTCPropertiesBuilder = JClass("thermogis.properties.builders.UTCPropertiesBuilder")
		propsBuilder = UTCPropertiesBuilder()
		utc_properties = propsBuilder.build()
		# Create an instance of a ThermoGISDoublet
		doublet = ThermoGISDoublet(Mockito.mock(Logger), RNG(123), utc_properties)

		# Act
		# Read Input grids
		input_grids = read_grid(self.model_path / "InputData" / "simplified__thick.zmap").to_dataset(name="thickness_mean")
		input_grids["thickness_sd"] = read_grid(self.model_path / "InputData" / "simplified__thick_sd.zmap")
		input_grids["ntg"] = read_grid(self.model_path / "InputData" / "simplified__ntg.zmap")
		input_grids["porosity"] = read_grid(self.model_path / "InputData" / "simplified__phi.zmap") / 100
		input_grids["depth"] = read_grid(self.model_path / "InputData" / "simplified__top.zmap")
		input_grids["hc_accum"] = read_grid(self.model_path / "InputData" / "simplified__hc_accum.zmap")
		input_grids["ln_permeability_mean"] = np.log(read_grid(self.model_path / "InputData" / "simplified__k.zmap"))
		input_grids["ln_permeability_sd"] = read_grid(self.model_path / "InputData" / "simplified__k_lnsd.zmap")

		output_grids = calculate_performance_across_dimensions(input_grids, doublet, p_values=[10, 50, 90])

tests/tests/test_show.py

0 → 100644

+5 −0

Original line number	Diff line number	Diff line
		from pythermogis import doublet_performace



		output_values = doublet_performace(poro, perm, ntg, ... , input_params)
		No newline at end of file