diff --git a/requirements.txt b/requirements.txt new file mode 100644 index 0000000..bc5aa98 --- /dev/null +++ b/requirements.txt @@ -0,0 +1,5 @@ +pytest>=8.3 +numpy>=1.24 +pftools>=1.3.11 +subsettools>=2.0 +hf_hydrodata>=1.3 \ No newline at end of file diff --git a/workflow/conftest.py b/workflow/conftest.py new file mode 100644 index 0000000..fb8156e --- /dev/null +++ b/workflow/conftest.py @@ -0,0 +1,46 @@ +"""Common pytest fixtures.""" +import os +import pytest +from parflow.tools.fs import mkdir, rm + +@pytest.fixture(scope="session") +def set_parflow_dir(): + """Set ParFlow directory environment variable.""" + old_parflow_dir = os.environ.get("PARFLOW_DIR") + if not old_parflow_dir: + # set PARFLOW_DIR to local installation where the tests will be run + os.environ["PARFLOW_DIR"] = "/home/SHARED/software/parflow/3.10.0" + yield None + if old_parflow_dir: + os.environ["PARFLOW_DIR"] = old_parflow_dir + else: + del os.environ["PARFLOW_DIR"] + + +@pytest.fixture(scope="module") +def setup_dir_structure(set_parflow_dir): + """Create and return paths for directories necessary for ParFlow run.""" + home = os.path.expanduser("~") + base_dir = os.path.join(home, "test_output") + input_dir = os.path.join(base_dir, "inputs") + output_dir = os.path.join(base_dir, "outputs") + + def setup(run_name): + static_write_dir = os.path.join(input_dir, run_name, "static") + mkdir(static_write_dir) + forcing_dir = os.path.join(input_dir, run_name, "forcing") + mkdir(forcing_dir) + pf_out_dir = os.path.join(output_dir, run_name) + mkdir(pf_out_dir) + correct_output_dir = os.path.abspath(os.path.join("correct_output", run_name)) + target_runscript = os.path.join(pf_out_dir, run_name + ".yaml") + return ( + static_write_dir, + forcing_dir, + pf_out_dir, + correct_output_dir, + target_runscript, + ) + + yield setup + rm(base_dir) diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00000.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00000.pfb new file mode 100644 index 0000000..f403947 Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00000.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00001.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00001.pfb new file mode 100644 index 0000000..2fb699d Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00001.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00002.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00002.pfb new file mode 100644 index 0000000..8878fde Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00002.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00003.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00003.pfb new file mode 100644 index 0000000..53e9a59 Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00003.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00004.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00004.pfb new file mode 100644 index 0000000..2a8f0a5 Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00004.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00005.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00005.pfb new file mode 100644 index 0000000..91208f5 Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00005.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00006.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00006.pfb new file mode 100644 index 0000000..bf0ebe8 Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00006.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00007.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00007.pfb new file mode 100644 index 0000000..6e7e2ae Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00007.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00008.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00008.pfb new file mode 100644 index 0000000..199cda1 Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00008.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00009.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00009.pfb new file mode 100644 index 0000000..5bc6ae0 Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00009.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00010.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00010.pfb new file mode 100644 index 0000000..32f12bd Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00010.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00011.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00011.pfb new file mode 100644 index 0000000..d0c27f9 Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00011.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00012.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00012.pfb new file mode 100644 index 0000000..4d0088f Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00012.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00013.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00013.pfb new file mode 100644 index 0000000..59b75b7 Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00013.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00014.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00014.pfb new file mode 100644 index 0000000..403a2b7 Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00014.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00015.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00015.pfb new file mode 100644 index 0000000..5fc5e7f Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00015.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00016.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00016.pfb new file mode 100644 index 0000000..a1de9f8 Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00016.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00017.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00017.pfb new file mode 100644 index 0000000..ac09835 Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00017.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00018.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00018.pfb new file mode 100644 index 0000000..d11164f Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00018.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00019.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00019.pfb new file mode 100644 index 0000000..dc0b664 Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00019.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00020.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00020.pfb new file mode 100644 index 0000000..22cd416 Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00020.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00021.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00021.pfb new file mode 100644 index 0000000..fe9c8dc Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00021.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00022.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00022.pfb new file mode 100644 index 0000000..ac3dd00 Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00022.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00023.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00023.pfb new file mode 100644 index 0000000..f1ec382 Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00023.pfb differ diff --git a/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00024.pfb b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00024.pfb new file mode 100644 index 0000000..c648312 Binary files /dev/null and b/workflow/correct_output/conus1_15060201/conus1_15060201.out.press.00024.pfb differ diff --git a/workflow/test_parflow_conus1_headwater.py b/workflow/test_parflow_conus1_headwater.py new file mode 100644 index 0000000..dd34ce6 --- /dev/null +++ b/workflow/test_parflow_conus1_headwater.py @@ -0,0 +1,204 @@ +""" +Tests to compare a full ParFlow-CONUS1 subsetting run against +established outputs for headwater watershed(s). +""" +import os +import pytest +import numpy as np +from parflow import Run, read_pfb +from parflow.tools.settings import set_working_directory +import subsettools as st +import hf_hydrodata as hf +from testutils import pf_test_msig_diff + +# This is the list of HUCs included in pytest parameterize +HUC_LIST = ["15060201"] + +# Number of timesteps (hours) to run simulation for and compare against +NUM_TIMESTEPS = 24 + +@pytest.fixture(scope="module", autouse=True) +def setup_run(setup_dir_structure, huc_id): + """ + Fixture to set up and run ParFlow. + This pytest fixture will be run automatically for use of all + tests within this module and setting it up in this way + ensures it is only run once per module. So all + tests within this module will use the same ParFlow run, + without having to run ParFlow multiple times. + """ + + # Define run_name based on grid and HUC provided + grid = "conus1" + run_name = f"{grid}_{huc_id}" + + # These could be parameterized if in the future multiple time ranges are desired to be tested + start = "2002-10-01" + end = "2002-10-03" + + + run_ds = "conus1_baseline_mod" + var_ds = "conus1_domain" + forcing_ds = "NLDAS2" + + # Set up directories for ParFlow run + ( + static_write_dir, + forcing_dir, + pf_out_dir, + correct_output_dir, + target_runscript, + ) = setup_dir_structure(run_name) + + # set cluster topology + P = 1 + Q = 1 + + # load template runscript + reference_run = st.get_template_runscript( + grid, "transient", "solid", pf_out_dir + ) + + # get ParFlow i/j bounding box + ij_bounds, mask = st.define_huc_domain(hucs=[huc_id], grid=grid) + + # create mask and solid files + st.write_mask_solid(mask=mask, grid=grid, write_dir=static_write_dir) + + # subset static inputs and pressure + st.subset_static( + ij_bounds, + dataset=var_ds, + write_dir=static_write_dir, + var_list=( + "slope_x", + "slope_y", + "pf_indicator", + "pme", + "ss_pressure_head", + ), + ) + init_press_filepath = st.subset_press_init( + ij_bounds, + dataset=run_ds, + date=start, + write_dir=static_write_dir, + time_zone="UTC", + ) + + # Configure CLM drivers + st.config_clm( + ij_bounds, start=start, end=end, dataset=run_ds, write_dir=static_write_dir + ) + + # Subset climate forcing + st.subset_forcing( + ij_bounds, + grid=grid, + start=start, + end=end, + dataset=forcing_ds, + write_dir=forcing_dir, + ) + + # set up baseline run from template + target_runscript = st.edit_runscript_for_subset( + ij_bounds, + runscript_path=reference_run, + runname=run_name, + forcing_dir=forcing_dir, + ) + + # copy static files into run directory + st.copy_files(read_dir=static_write_dir, write_dir=pf_out_dir) + + # change file names + target_runscript = st.change_filename_values( + runscript_path=target_runscript, + init_press=os.path.basename(init_press_filepath), + ) + + # change processor topology and re-distribute forcing + target_runscript = st.dist_run( + topo_p=P, + topo_q=Q, + runscript_path=target_runscript, + dist_clim_forcing=True, + ) + + # Set ParFlow output dir + set_working_directory(pf_out_dir) + + # Run ParFlow + run = Run.from_definition(target_runscript) + run.TimingInfo.StopTime = NUM_TIMESTEPS + run.run(working_directory=pf_out_dir) + + return (pf_out_dir, run_name, correct_output_dir) + +@pytest.mark.parametrize("huc_id", HUC_LIST, scope="module") +def test_parflow_conus1_results(setup_run, huc_id): + """ + Test to compare run of subsettools workflow against + published ParFlow-CONUS1 results. + """ + + # Get run information + (pf_out_dir, run_name, _) = setup_run + + # Compare to published ParFlow-CONUS simulation results + conus1_results = hf.get_gridded_data(dataset="conus1_baseline_mod", variable="pressure_head", + huc_id=huc_id, + start_time="2002-10-01", end_time="2002-10-03", + temporal_resolution="hourly") + + # Read in mask. We will mask results at each timestep so the area + # outside of the HUC domain is set to NaN + mask = read_pfb(f"{pf_out_dir}/{run_name}.out.mask.pfb") + + # Test each timestep + for i in range(NUM_TIMESTEPS): + timestep = str(i).rjust(5, "0") + parflow_subset_results = read_pfb(f"{pf_out_dir}/{run_name}.out.press.{timestep}.pfb") + masked_subset_results = np.where(mask == 0, np.nan, parflow_subset_results) + + print(f">>>>>> checking timestep {i}") + + # ensure the shapes of the arrays returned are the same + assert masked_subset_results.shape == conus1_results[i, :, :, :].shape + + # The value bounds on these tests are broad and probably specific to the HUC + # currently chosen. There is an RSE-hydrologist buddy task defined in Notion + # to investigate the cause of these differences. Once those are determined, + # this section of the test can be updated to decrease the allowed threshold + # of difference. Ideally, the test could use pf_test_msig_diff or similar, + # to test closeness of arrays in the way that the ParFlow testing suite does. + # The arrays are currently too different to implement that kind of test meaningfully + # at the moment. + for j in range(5): + diff = masked_subset_results[j] - conus1_results[i, j, :, :] + if j == 0: + assert (np.nanmax(diff) > -0.6) & (np.nanmax(diff) < 0.75) + else: + assert (np.nanmax(diff) > -0.2) & (np.nanmax(diff) < 1.55) + + +@pytest.mark.parametrize("huc_id", HUC_LIST, scope="module") +def test_parflow_conus1_reproducible(setup_run, huc_id): + """ + Test to compare run of subsettools workflow against + pre-produced and hydrologist-reviewed outputs of a + workflow run to ensure the workflow is reproducible. + """ + + # Get run information + (pf_out_dir, run_name, correct_output_dir) = setup_run + + # Test each timestep + for i in range(NUM_TIMESTEPS): + timestep = str(i).rjust(5, "0") + run_results = read_pfb(f"{pf_out_dir}/{run_name}.out.press.{timestep}.pfb") + correct_output = read_pfb(f"{correct_output_dir}/{run_name}.out.press.{timestep}.pfb") + + print(f">>>>>> checking timestep {i}") + assert pf_test_msig_diff(run_results, correct_output, 4) diff --git a/workflow/testutils.py b/workflow/testutils.py new file mode 100644 index 0000000..14fd55b --- /dev/null +++ b/workflow/testutils.py @@ -0,0 +1,262 @@ +""" + Utility functions for workflow testing. + Several of these functions come from the ParFlow repo + and are commonly used for testing ParFlow outputs against + each other. +""" +# pylint: disable=W0632,W0718 + +import os +from parflow.tools.io import read_pfb + + +def msig_diff(data1, data2, m, abs_zero=0.0): + """Python version of the C MSigDiff function. + + Two nd arrays are given as the first two arguments. + The grid point at which the number of digits in agreement + (significant digits) is fewest is determined. If m >= 0 + then the coordinate whose two values differ in more than + m significant digits will be computed. If m < 0 then the + coordinate whose values have a minimum number of significant + digits will be computed. The number of the fewest significant + digits is determined, and the maximum absolute difference is + computed. The only coordinates that will be considered will be + those whose differences are greater than absolute zero. + + Args: + data1 (numpy.ndarray): first ndarray + data2 (numpy.ndarray): second ndarray + m (int): number of significant digits + abs_zero (float): threshold below which all numbers are considered to be 0 + + Returns: + A list of the following form is returned upon success: + + [[i j k s] max_adiff] + + where i, j, and k are the coordinates computed, sd is the + minimum number of significant digits computed, and max_adiff + is the maximum absolute difference computed. + + Raises: + ValueError: If data1 and data2 do not have the same shape. + """ + assert isinstance(abs_zero, float) + assert abs_zero >= 0 + + if not data1.shape == data2.shape: + print("data shape:") + print(data1.shape) + print("correct data shape:") + print(data2.shape) + raise ValueError("Error: Data arrays must have the same dimensions.") + + nx, ny, nz = data1.shape + + if m >= 0: + sig_dig_rhs = 0.5 / 10**m + else: + sig_dig_rhs = 0.0 + + m_sig_digs_everywhere = True + max_sdiff = 0.0 + max_adiff = 0.0 + + for k in range(nz): + for j in range(ny): + for i in range(nx): + adiff = abs(data1[i, j, k] - data2[i, j, k]) + amax = max(abs(data1[i, j, k]), abs(data2[i, j, k])) + + max_adiff = max(max_adiff, adiff) + + m_sig_digs = True + if amax > abs_zero: + sdiff = adiff / amax + if sdiff > sig_dig_rhs: + m_sig_digs = False + + if not m_sig_digs: + if sdiff > max_sdiff: + max_sdiff = sdiff + mi, mj, mk = i, j, k + + m_sig_digs_everywhere = False + + result = [] + if not m_sig_digs_everywhere: + sig_digs = 0 + sdiff = max_sdiff + while sdiff <= 0.5e-01: + sdiff *= 10.0 + sig_digs += 1 + + result.append([mi, mj, mk, sig_digs]) + result.append(max_adiff) + + return result + +def pf_test_msig_diff(data, correct_data, message, sig_digits=6): + """Wrapper around msig_diff function to produce boolean response. + + Two file paths are given as the first two arguments. + The function reads them into ndarrays and calls a comparison + function (msig_diff) to check if the files differ by + more than sig_digits significant digits in any coordinate. + If they do, the function prints an error message and the + coordinate in which the files have the greatest difference. + + Args: + file (string): path to test output file + correct_file (string): path to reference output file + message (string): message to print in case of test failure + sig_digits (int): number of significant digits + + Returns: + bool: True if files differ in no more than sig_digits + significant digits in all coordinates. False otherwise. + """ + result = msig_diff(data, correct_data, sig_digits) + if (len(result)) == 0: + return True + + m_sig_digs, max_abs_diff = result + + i, j, k, sig_digs = m_sig_digs + + print(f"FAILED : {message}") + print(f"\tMinimum significant digits at ({i:3d}, {j:3d}, {k:3d}) = {sig_digs:2d}") + print(f"\tCorrect value {correct_data[i, j, k]:e}") + print(f"\tComputed value {data[i, j, k]:e}") + + elt_diff = abs(data[i, j, k] - correct_data[i, j, k]) + print(f"\tDifference {elt_diff:e}") + + print(f"\tMaximum absolute difference = {max_abs_diff:e}") + + return False + +def pf_test_file(file, correct_file, message, sig_digits=6): + """Python version of the tcl pftestFile procedure. + + Two file paths are given as the first two arguments. + The function reads them into ndarrays and calls a comparison + function (msig_diff) to check if the files differ by + more than sig_digits significant digits in any coordinate. + If they do, the function prints an error message and the + coordinate in which the files have the greatest difference. + + Args: + file (string): path to test output file + correct_file (string): path to reference output file + message (string): message to print in case of test failure + sig_digits (int): number of significant digits + + Returns: + bool: True if files differ in no more than sig_digits + significant digits in all coordinates. False otherwise. + + Raises: + FileNotFoundError: If file or correct_file do not exist. + Exception: If the function fails to read either file into an ndarray. + """ + if not os.path.exists(file): + raise FileNotFoundError(f"FAILED : output file <{file}> not created") + + if not os.path.exists(correct_file): + raise FileNotFoundError( + f"FAILED : regression check output file <{correct_file}> does not exist" + ) + + try: + data = read_pfb(file) + correct_data = read_pfb(correct_file) + except Exception as e: + print("Error: Failed to load data from files.", e) + return False + + result = msig_diff(data, correct_data, sig_digits) + if (len(result)) == 0: + return True + + m_sig_digs, max_abs_diff = result + + i, j, k, sig_digs = m_sig_digs + + print(f"FAILED : {message}") + print(f"\tMinimum significant digits at ({i:3d}, {j:3d}, {k:3d}) = {sig_digs:2d}") + print(f"\tCorrect value {correct_data[i, j, k]:e}") + print(f"\tComputed value {data[i, j, k]:e}") + + elt_diff = abs(data[i, j, k] - correct_data[i, j, k]) + print(f"\tDifference {elt_diff:e}") + + print(f"\tMaximum absolute difference = {max_abs_diff:e}") + + return False + + +def pf_test_file_with_abs(file, correct_file, message, abs_value, sig_digits=6): + """Python version of the tcl pftestFileWithAbs procedure. + + Two file paths are given as the first two arguments. + The function reads them into ndarrays and calls a comparison + function (msig_diff) to check if the files differ by + more than sig_digits significant digits in any coordinate. + If they do, the function checks if the difference in that + coordinate is greater than abs_value. If it is, it prints + an error message and the coordinate in which the files have + the greatest difference. + + Args: + file (string): path to test output file + correct_file (string): path to reference output file + message (string): message to print in case of test failure + abs_value (float): threshold to determine if two files differ + enough to fail the test + sig_digits (int): number of significant digits + + Returns: + bool: True if files differ in no more than sig_digits + significant digits in all coordinates. False otherwise. + + Raises: + FileNotFoundError: If file or correct_file do not exist. + Exception: If the function fails to read either file into an ndarray. + """ + if not os.path.exists(file): + raise FileNotFoundError(f"FAILED : output file <{file}> not created") + + if not os.path.exists(correct_file): + raise FileNotFoundError( + f"FAILED : regression check output file <{correct_file}> does not exist" + ) + + try: + data = read_pfb(file) + correct_data = read_pfb(correct_file) + except Exception as e: + print("Error: Failed to load data from files.", e) + return False + + result = msig_diff(data, correct_data, sig_digits) + if len(result) != 0: + m_sig_digs, max_abs_diff = result + + i, j, k, sig_digs = m_sig_digs + + elt_diff = abs(data[i, j, k] - correct_data[i, j, k]) + + if elt_diff > abs_value: + print(f"FAILED : {message}") + print( + f"\tMinimum significant digits at ({i:3d}, {j:3d}, {k:3d}) = {sig_digs:2d}" + ) + print(f"\tCorrect value {correct_data[i, j, k]:e}") + print(f"\tComputed value {data[i, j, k]:e}") + print(f"\tDifference {elt_diff:e}") + print(f"\tMaximum absolute difference = {max_abs_diff:e}") + return False + + return True