TNO Intern

Commit 65d43a41 authored by Hen Brett's avatar Hen Brett 🐔
Browse files

Merge branch '24-add-direct-heat-and-chiller-tests' into 'main' 

Resolve "Add direct heat and chiller tests"

Closes #24

See merge request AGS/pythermogis!24
parents 548f7562 3bb9f9df

pyproject.toml

+1 −1
Original line number Diff line number Diff line
@@ -4,7 +4,7 @@ build-backend = "setuptools.build_meta" 


[project]

name = "pythermogis"

version = "0.1.15"

version = "0.1.16"

description = "This is a Python API to the ThermoGIS Java code"

authors = [{ name = "Brett", email = "hen.brett@tno.nl" }]

readme = "README.md"

tests/test_ThermoGISDoublet_Benchmark.py

0 → 100644
+251 −0
Original line number Diff line number Diff line 
from unittest import TestCase



import numpy as np

from jpype.types import *



from pythermogis.thermogis_classes.jvm_start import start_jvm

from pythermogis.thermogis_classes.utc_properties import instantiate_utc_properties_builder





class ThermoGISDoubletBenchmark(TestCase):

    """

    This is a series of tests which copy exactly tests found in the Java code; this ensures that the Jpype implementation of Java still produces the same results as the native java implementation

    of the ThermoGISDoublet class.

    """



    def test_calculateDoubletPerformanceTest(self):

        """

        This is a copy of a test from the Java ThermoGisDoubletTest.java script; to validate that this python implementation of the ThermoGIS Doublet

        returns the same values.

        """

        # Arrange

        # Import Java Classes

        start_jvm()

        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")

        ViscosityMode = JClass("tno.geoenergy.ViscosityMode")



        # Instantiate the UTC properties class

        utc_properties = UTCPropertiesBuilder().build()



        # Create an instance of a ThermoGISDoublet

        doublet = ThermoGISDoublet(Mockito.mock(Logger), RNG(0), utc_properties)



        # Setting input values

        thickness = 100

        permeability = 175

        transmissivity = thickness * permeability

        doublet.doubletCalc1DData.setReservoirTemp(76)

        doublet.doubletCalc1DData.setDepth(2000)

        doublet.doubletCalc1DData.setInjectionTemp(40)

        doublet.doubletCalc1DData.setNtg(1.0)

        doublet.doubletCalc1DData.setViscosityMode(ViscosityMode.KESTIN)

        doublet.setNoStimulation()



        # Act

        doublet.calculateDoubletPerformance(-999, thickness, transmissivity)



        # Assert

        self.assertTrue(np.isclose(17500, transmissivity, 0.001))

        self.assertTrue(np.isclose(227.2757568359375, doublet.doubletCalc1DData.getFlowrate(), 1))

        self.assertTrue(np.isclose(60, doublet.doubletCalc1DData.getDrawdownPressure() / 1e5, 0.001))

        self.assertTrue(np.isclose(6.616096470753937, doublet.getUtcPeurctkWh(), 0.001))

        self.assertTrue(np.isclose(1159.17968, doublet.doubletCalc1DData.getWellDistance(), 0.001))

        self.assertTrue(np.isclose(8.636903762817383, doublet.doubletCalc1DData.getHeatPowerPerDoublet(), 0.001))

        self.assertTrue(np.isclose(13.627557754516602, doublet.doubletCalc1DData.getCop(), 0.001))



    def test_calculateDoubletPerformance_directHeat(self):

        """

        This is a copy of a test from the Java ThermoGisDoubletTest.java script; to validate that this python implementation of the ThermoGIS Doublet

        returns the same values.

        """

        # Arrange

        start_jvm()

        Logger = JClass("logging.Logger")

        Mockito = JClass("org.mockito.Mockito")

        RNG = JClass("tno.geoenergy.stochastic.RandomNumberGenerator")

        ThermoGISDoublet = JClass("thermogis.calc.doublet.ThermoGisDoublet")

        ViscosityMode = JClass("tno.geoenergy.ViscosityMode")



        # Instantiate the UTC properties class

        utc_properties = self.setup_template_utc_properties_builder().setOpexBase(0).build()



        # Create an instance of a ThermoGISDoublet

        doublet = ThermoGISDoublet(Mockito.mock(Logger), RNG(0), utc_properties)



        # Setting input values

        thickness = 100

        permeability = 175

        transmissivity = thickness * permeability

        doublet.doubletCalc1DData.setReservoirTemp(76)

        doublet.doubletCalc1DData.setDepth(2000)

        doublet.doubletCalc1DData.setInjectionTemp(40)

        doublet.doubletCalc1DData.setNtg(1.0)

        doublet.doubletCalc1DData.setViscosityMode(ViscosityMode.KESTIN)

        doublet.setNoStimulation()



        # Act

        doublet.calculateDoubletPerformance(-999, thickness, transmissivity)



        # Assert

        self.assertTrue(np.isclose(17500, transmissivity, 0.001))

        self.assertTrue(np.isclose(227.2757568359375, doublet.doubletCalc1DData.getFlowrate(), 1))

        self.assertTrue(np.isclose(60, doublet.doubletCalc1DData.getDrawdownPressure() / 1e5, 0.001))

        self.assertTrue(np.isclose(8.624696731567383, doublet.doubletCalc1DData.getHeatPowerPerDoublet(), 0.001))

        self.assertTrue(np.isclose(1159.1796, doublet.doubletCalc1DData.getWellDistance(), 0.001))

        self.assertTrue(np.isclose(13.623167037963867, doublet.doubletCalc1DData.getCop(), 0.001))

        self.assertTrue(np.isclose(5.229816400909403, doublet.getUtcPeurctkWh(), 0.001))

        self.assertTrue(np.isclose(16.439682499211536, doublet.economicalData.getSumcapex(), 0.001))



    def test_calculateDoubletPerformance_chiller(self):

        """

        This is a copy of a test from the Java ThermoGisDoubletTest.java script; to validate that this python implementation of the ThermoGIS Doublet

        returns the same values.

        """

        # Arrange

        start_jvm()

        Logger = JClass("logging.Logger")

        Mockito = JClass("org.mockito.Mockito")

        RNG = JClass("tno.geoenergy.stochastic.RandomNumberGenerator")

        ThermoGISDoublet = JClass("thermogis.calc.doublet.ThermoGisDoublet")

        ViscosityMode = JClass("tno.geoenergy.ViscosityMode")



        # Instantiate the UTC properties class

        utc_properties = (self.setup_template_utc_properties_builder().setCapexConst(0.5)

                          .setCapexVariable(1100)

                          .setHeatExchangerEfficiency(0.4)

                          .setHeatExchangerParasitic(0.1).build())



        # Create an instance of a ThermoGISDoublet

        doublet = ThermoGISDoublet(Mockito.mock(Logger), RNG(0), utc_properties)



        # Setting input values

        thickness = 100

        permeability = 175

        transmissivity = thickness * permeability

        doublet.doubletCalc1DData.setReservoirTemp(76)

        doublet.doubletCalc1DData.setDepth(2000)

        doublet.doubletCalc1DData.setInjectionTemp(60)

        doublet.doubletCalc1DData.setNtg(1.0)

        doublet.doubletCalc1DData.setViscosityMode(ViscosityMode.KESTIN)

        doublet.setNoStimulation()



        # Act

        doublet.calculateDoubletPerformance(-999, thickness, transmissivity)



        # Assert

        self.assertTrue(np.isclose(17500, transmissivity, 0.001))

        self.assertTrue(np.isclose(256.59625244140625, doublet.doubletCalc1DData.getFlowrate(), 1))

        self.assertTrue(np.isclose(60, doublet.doubletCalc1DData.getDrawdownPressure() / 1e5, 0.001))

        self.assertTrue(np.isclose(1.6594701766967774, doublet.doubletCalc1DData.getHeatPowerPerDoublet(), 0.001))

        self.assertTrue(np.isclose(1227.1484375, doublet.doubletCalc1DData.getWellDistance(), 0.001))

        self.assertTrue(np.isclose(1.8887300754346803, doublet.doubletCalc1DData.getCop(), 0.001))

        self.assertTrue(np.isclose(20.470115103822685, doublet.getUtcPeurctkWh(), 0.001))

        self.assertTrue(np.isclose(12.748051248109613, doublet.economicalData.getSumcapex(), 0.001))



    def test_calculateDoubletPerformance_directheatHP(self):

        """

        This is a copy of a test from the Java ThermoGisDoubletTest.java script; to validate that this python implementation of the ThermoGIS Doublet

        returns the same values.

        """

        # Arrange

        start_jvm()

        Logger = JClass("logging.Logger")

        Mockito = JClass("org.mockito.Mockito")

        RNG = JClass("tno.geoenergy.stochastic.RandomNumberGenerator")

        ThermoGISDoublet = JClass("thermogis.calc.doublet.ThermoGisDoublet")

        ViscosityMode = JClass("tno.geoenergy.ViscosityMode")



        # Instantiate the UTC properties class

        utc_properties = self.setup_template_utc_properties_builder().setOpexPerPower(100).setOpexBase(0).build()

        # Create an instance of a ThermoGISDoublet

        doublet = ThermoGISDoublet(Mockito.mock(Logger), RNG(0), utc_properties)



        # Setting input values

        thickness = 100

        permeability = 175

        transmissivity = thickness * permeability

        doublet.doubletCalc1DData.setReservoirTemp(50)

        doublet.doubletCalc1DData.setDepth(2000)

        doublet.doubletCalc1DData.setDhReturnTemp(35)

        doublet.doubletCalc1DData.setInjectionTemp(doublet.calculateInjectionTempWithHeatPump(50, 70))

        doublet.doubletCalc1DData.setNtg(1.0)

        doublet.doubletCalc1DData.setViscosityMode(ViscosityMode.KESTIN)

        doublet.setNoStimulation()



        # Act

        doublet.calculateDoubletPerformance(-999, thickness, transmissivity)



        # Assert

        self.assertTrue(np.isclose(17499.99940, transmissivity, 0.001))

        self.assertTrue(np.isclose(163.99771118164062, doublet.doubletCalc1DData.getFlowrate(), 1))

        self.assertTrue(np.isclose(60, doublet.doubletCalc1DData.getDrawdownPressure() / 1e5, 0.001))

        self.assertTrue(np.isclose(4.97566556930542, doublet.doubletCalc1DData.getHeatPowerPerDoublet(), 0.001))

        self.assertTrue(np.isclose(989.2578125, doublet.doubletCalc1DData.getWellDistance(), 0.001))

        self.assertTrue(np.isclose(10.922356605529785, doublet.doubletCalc1DData.getCop(), 0.001))

        self.assertTrue(np.isclose(8.45577723033848, doublet.getUtcPeurctkWh(), 0.001))

        self.assertTrue(np.isclose(15.046041336220961, doublet.economicalData.getSumcapex(), 0.001))



    def test_calculateDoubletPerformance_ORC(self):

        """

        This is a copy of a test from the Java ThermoGisDoubletTest.java script; to validate that this python implementation of the ThermoGIS Doublet

        returns the same values.

        """

        # Arrange

        start_jvm()

        Logger = JClass("logging.Logger")

        Mockito = JClass("org.mockito.Mockito")

        RNG = JClass("tno.geoenergy.stochastic.RandomNumberGenerator")

        ThermoGISDoublet = JClass("thermogis.calc.doublet.ThermoGisDoublet")

        ViscosityMode = JClass("tno.geoenergy.ViscosityMode")



        # Instantiate the UTC properties class

        utc_properties = (self.setup_template_utc_properties_builder().setCapexConst(0)

                          .setCapexVariable(2300)

                          .setHeatExchangerEfficiency(0.6)

                          .setUseORC(True)

                          .setHeatExchangerBasetemp(20).build())



        # Create an instance of a ThermoGISDoublet

        doublet = ThermoGISDoublet(Mockito.mock(Logger), RNG(0), utc_properties)



        # Setting input values

        thickness = 100

        permeability = 175

        transmissivity = thickness * permeability

        doublet.doubletCalc1DData.setReservoirTemp(100)

        doublet.doubletCalc1DData.setDepth(2000)

        doublet.doubletCalc1DData.setInjectionTemp(60)

        doublet.doubletCalc1DData.setNtg(1.0)

        doublet.doubletCalc1DData.setViscosityMode(ViscosityMode.KESTIN)

        doublet.setNoStimulation()



        # Act

        doublet.calculateDoubletPerformance(-999, thickness, transmissivity)



        # Assert

        self.assertTrue(np.isclose(17499.99940, transmissivity, 0.001))

        self.assertTrue(np.isclose(293.6246643066406, doublet.doubletCalc1DData.getFlowrate(), 1))

        self.assertTrue(np.isclose(60, doublet.doubletCalc1DData.getDrawdownPressure() / 1e5, 0.001))

        self.assertTrue(np.isclose(0.05274631495788107, doublet.doubletCalc1DData.getHeatPowerPerDoublet(), 0.01))

        self.assertTrue(np.isclose(1306.4453125, doublet.doubletCalc1DData.getWellDistance(), 0.001))

        self.assertTrue(np.isclose(0.06459403120103477, doublet.doubletCalc1DData.getCop(), 0.001))

        self.assertTrue(np.isclose(36.98296076530068, doublet.getUtcPeurctkWh(), 0.001))

        self.assertTrue(np.isclose(12.44409167482118, doublet.economicalData.getSumcapex(), 0.001))



    def setup_template_utc_properties_builder(self):

        return (instantiate_utc_properties_builder()

         .setSalinityGradient(47.0)

         .setDrillingTime(1)

         .setTaxRate(20)

         .setEquityReturn(15)

         .setOpexBase(10000)

         .setOpexPerPower(50)

         .setWellCostScaling(1)

         .setOpexPerEnergy(0)

         .setWellCostConst(0.375)

         .setWellCostZ(1050)

         .setWellCostZ2(0.3))
 
 No newline at end of file

tests/test_doublet.py

+26 −70
Original line number Diff line number Diff line
@@ -36,50 +36,6 @@ class PyThermoGIS(TestCase): 
    def tearDown(self):

        shutil.rmtree(self.test_files_out_path)



    def test_calculateDoubletPerformanceTest(self):

        """

        This is a copy of a test from the Java ThermoGisDoubletTest.java script; to validate that this python implementation of the ThermoGIS Doublet

        returns the same values. This does not use any of the extra interfaces incorporated here.

        """

        # Arrange

        # Import Java Classes

        start_jvm()

        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")

        ViscosityMode = JClass("tno.geoenergy.ViscosityMode")



        # Instantiate the UTC properties class

        utc_properties = UTCPropertiesBuilder().build()



        # Create an instance of a ThermoGISDoublet

        doublet = ThermoGISDoublet(Mockito.mock(Logger), RNG(0), utc_properties)



        # Setting input values

        thickness = 100

        permeability = 175

        transmissivity = thickness * permeability

        doublet.doubletCalc1DData.setReservoirTemp(76)

        doublet.doubletCalc1DData.setDepth(2000)

        doublet.doubletCalc1DData.setInjectionTemp(40)

        doublet.doubletCalc1DData.setNtg(1.0)

        doublet.doubletCalc1DData.setViscosityMode(ViscosityMode.KESTIN)

        doublet.setNoStimulation()



        # Act

        doublet.calculateDoubletPerformance(-999, thickness, transmissivity)



        # Assert

        self.assertTrue(np.isclose(17500, transmissivity, 0.001))

        self.assertTrue(np.isclose(227.2757568359375, doublet.doubletCalc1DData.getFlowrate(), 1))

        self.assertTrue(np.isclose(60, doublet.doubletCalc1DData.getDrawdownPressure() / 1e5, 0.001))

        self.assertTrue(np.isclose(6.616096470753937, doublet.getUtcPeurctkWh(), 0.001))

        self.assertTrue(np.isclose(1159.17968, doublet.doubletCalc1DData.getWellDistance(), 0.001))

        self.assertTrue(np.isclose(8.636903762817383, doublet.doubletCalc1DData.getHeatPowerPerDoublet(), 0.001))

        self.assertTrue(np.isclose(13.627557754516602, doublet.doubletCalc1DData.getCop(), 0.001))



    def test_doublet_grid_basecase(self):

        # This tests that the python API produces (approximately) the same output as the Java code, when running on the same set of input data for the base case scenario

        # Run calculation across all dimensions of input_grids, and all provided P_values
@@ -137,32 +93,6 @@ class PyThermoGIS(TestCase): 
            output_grids_p_value = output_grids_p_value.drop_vars("p_value")

            self.assert_output_and_benchmark_is_close(self.test_benchmark_output_HP_STIM, output_grids_p_value, p_value, "_HP_STIM")



    def assert_output_and_benchmark_is_close(self, benchmark_path, output_grids, p_value, scenario):

        xr.testing.assert_allclose(output_grids[f"permeability"], read_grid(benchmark_path / f"simplified__k_P{p_value}{scenario}.nc"), atol=0.00000001)

        xr.testing.assert_allclose(output_grids[f"thickness"], read_grid(benchmark_path / f"simplified__h_P{p_value}{scenario}.nc"), atol=0.00000001)

        # The large differences in the output grids come from the fact that the way we calculate the Transmissivity P-values in both this python implementation and the java code is not deterministic;

        # we use a sampling method which is imprecise, and without the same random number generator + seed, will not produce exactly the same values between the python and java code

        xr.testing.assert_allclose(output_grids[f"transmissivity_with_ntg"], read_grid(benchmark_path / f"simplified__kh_P{p_value}{scenario}.nc"), atol=4)

        xr.testing.assert_allclose(output_grids.power, read_grid(benchmark_path / f"simplified__power_P{p_value}{scenario}.nc"), atol=1.0)

        xr.testing.assert_allclose(output_grids.utc, read_grid(benchmark_path / f"simplified__utc_P{p_value}{scenario}.nc"), atol=0.1)

        xr.testing.assert_allclose(output_grids.npv, read_grid(benchmark_path / f"simplified__npv_P{p_value}{scenario}.nc"), atol=0.5)

        xr.testing.assert_allclose(output_grids.hprod, read_grid(benchmark_path / f"simplified__hprod_P{p_value}{scenario}.nc"), atol=30000)

        xr.testing.assert_allclose(output_grids.cop, read_grid(benchmark_path / f"simplified__cop_P{p_value}{scenario}.nc"), atol=1.0)

        xr.testing.assert_allclose(output_grids.pres, read_grid(benchmark_path / f"simplified__pres_P{p_value}{scenario}.nc"), atol=2.0)

        xr.testing.assert_allclose(output_grids.flow_rate, read_grid(benchmark_path / f"simplified__flowr_P{p_value}{scenario}.nc"), atol=5)

        xr.testing.assert_allclose(output_grids.welld, read_grid(benchmark_path / f"simplified__welld_P{p_value}{scenario}.nc"), atol=30)



    def read_input_grids(self):

        input_grids = read_grid(self.input_data / "InputData" / "simplified__thick.zmap").to_dataset(name="thickness_mean")

        input_grids["thickness_sd"] = read_grid(self.input_data / "InputData" / "simplified__thick_sd.zmap")

        input_grids["ntg"] = read_grid(self.input_data / "InputData" / "simplified__ntg.zmap")

        input_grids["porosity"] = read_grid(self.input_data / "InputData" / "simplified__phi.zmap") / 100

        input_grids["depth"] = read_grid(self.input_data / "InputData" / "simplified__top.zmap")

        input_grids["mask"] = read_grid(self.input_data / "InputData" / "simplified__hc_accum.zmap")

        input_grids["ln_permeability_mean"] = np.log(read_grid(self.input_data / "InputData" / "simplified__k.zmap"))

        input_grids["ln_permeability_sd"] = read_grid(self.input_data / "InputData" / "simplified__k_lnsd.zmap")

        return input_grids



    def test_doublet_single_values(self):

        # This tests that the python API runs on a simple set of input with single values

        input_data = xr.Dataset({
@@ -193,3 +123,29 @@ class PyThermoGIS(TestCase): 
        # Act & Assert

        with self.assertRaises(ValueError):

            calculate_doublet_performance(input_data)



    def assert_output_and_benchmark_is_close(self, benchmark_path, output_grids, p_value, scenario):

        xr.testing.assert_allclose(output_grids[f"permeability"], read_grid(benchmark_path / f"simplified__k_P{p_value}{scenario}.nc"), atol=0.00000001)

        xr.testing.assert_allclose(output_grids[f"thickness"], read_grid(benchmark_path / f"simplified__h_P{p_value}{scenario}.nc"), atol=0.00000001)

        # The large differences in the output grids come from the fact that the way we calculate the Transmissivity P-values in both this python implementation and the java code is not deterministic;

        # we use a sampling method which is imprecise, and without the same random number generator + seed, will not produce exactly the same values between the python and java code

        xr.testing.assert_allclose(output_grids[f"transmissivity_with_ntg"], read_grid(benchmark_path / f"simplified__kh_P{p_value}{scenario}.nc"), atol=5)

        xr.testing.assert_allclose(output_grids.power, read_grid(benchmark_path / f"simplified__power_P{p_value}{scenario}.nc"), atol=1.0)

        xr.testing.assert_allclose(output_grids.utc, read_grid(benchmark_path / f"simplified__utc_P{p_value}{scenario}.nc"), atol=0.1)

        xr.testing.assert_allclose(output_grids.npv, read_grid(benchmark_path / f"simplified__npv_P{p_value}{scenario}.nc"), atol=0.5)

        xr.testing.assert_allclose(output_grids.hprod, read_grid(benchmark_path / f"simplified__hprod_P{p_value}{scenario}.nc"), atol=30000)

        xr.testing.assert_allclose(output_grids.cop, read_grid(benchmark_path / f"simplified__cop_P{p_value}{scenario}.nc"), atol=1.0)

        xr.testing.assert_allclose(output_grids.pres, read_grid(benchmark_path / f"simplified__pres_P{p_value}{scenario}.nc"), atol=2.0)

        xr.testing.assert_allclose(output_grids.flow_rate, read_grid(benchmark_path / f"simplified__flowr_P{p_value}{scenario}.nc"), atol=5)

        xr.testing.assert_allclose(output_grids.welld, read_grid(benchmark_path / f"simplified__welld_P{p_value}{scenario}.nc"), atol=30)



    def read_input_grids(self):

        input_grids = read_grid(self.input_data / "InputData" / "simplified__thick.zmap").to_dataset(name="thickness_mean")

        input_grids["thickness_sd"] = read_grid(self.input_data / "InputData" / "simplified__thick_sd.zmap")

        input_grids["ntg"] = read_grid(self.input_data / "InputData" / "simplified__ntg.zmap")

        input_grids["porosity"] = read_grid(self.input_data / "InputData" / "simplified__phi.zmap") / 100

        input_grids["depth"] = read_grid(self.input_data / "InputData" / "simplified__top.zmap")

        input_grids["mask"] = read_grid(self.input_data / "InputData" / "simplified__hc_accum.zmap")

        input_grids["ln_permeability_mean"] = np.log(read_grid(self.input_data / "InputData" / "simplified__k.zmap"))

        input_grids["ln_permeability_sd"] = read_grid(self.input_data / "InputData" / "simplified__k_lnsd.zmap")

        return input_grids
 
 No newline at end of file