Skip to content

Commit c010735

Browse files
authored
Compute at compuation time (#1042)
# Description New PR for the evaluating the computation time for custom grids. Test now includes model computation for pure model, after adding a topography and using the `compute_model_at` function with 1000 random points. I removed the test that deep created to clean up. Just fyi - here are the resulting times on my computer: ![image](https://github.com/user-attachments/assets/c4efea44-bdb2-4738-b936-bcdd65ff5d01) Added a warning message to `compute_model_at` function to alert users about potential side effects when setting a custom grid and computing the model. Modified the `set_custom_grid` function to include a `reset` parameter, which is now set to `True` in `compute_model_at` to ensure proper grid reset. Removed code in the `solutions` setter that was setting solutions per group, as this was causing issues with the custom grid functionality. Added comprehensive tests to measure and compare computation times across different grid configurations: - Dense grid computation (with consistency checks) - Topography grid vs dense grid computation - Custom grid computation with random coordinates # Checklist - [x] My code uses type hinting for function and method arguments and return values. - [x] I have created tests which cover my code. - [x] The test code either 1. demonstrates at least one valuable use case (e.g. integration tests) or 2. verifies that outputs are as expected for given inputs (e.g. unit tests). - [x] New tests pass locally with my changes. Relates to: gempy-project/gempy_engine#16 (comment)
2 parents 0ed5afd + 65743a9 commit c010735

File tree

5 files changed

+127
-13
lines changed

5 files changed

+127
-13
lines changed

gempy/API/compute_API.py

Lines changed: 5 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -82,11 +82,14 @@ def compute_model_at(gempy_model: GeoModel, at: np.ndarray,
8282
Returns:
8383
np.ndarray: The computed geological model at the specified coordinates.
8484
"""
85+
86+
print("WARNING: This function sets a custom grid and computes the model so be wary of side effects.")
8587
set_custom_grid(
8688
grid=gempy_model.grid,
87-
xyz_coord=at
89+
xyz_coord=at,
90+
reset=True
8891
)
89-
92+
9093
sol = compute_model(gempy_model, engine_config, validate_serialization=True)
9194
return sol.raw_arrays.custom
9295

gempy/API/grid_API.py

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -88,11 +88,11 @@ def set_topography_from_file(grid: Grid, filepath: str, crop_to_extent: Union[Se
8888
return set_topography_from_subsurface_structured_grid(grid, struct)
8989

9090

91-
def set_custom_grid(grid: Grid, xyz_coord: np.ndarray):
91+
def set_custom_grid(grid: Grid, xyz_coord: np.ndarray, reset: bool = False):
9292
custom_grid = CustomGrid(values=xyz_coord)
9393
grid.custom_grid = custom_grid
9494

95-
set_active_grid(grid, [Grid.GridTypes.CUSTOM])
95+
set_active_grid(grid, grid_type=[Grid.GridTypes.CUSTOM], reset=reset)
9696
return grid.custom_grid
9797

9898

gempy/core/data/geo_model.py

Lines changed: 0 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -101,14 +101,6 @@ def solutions(self, value: Solutions):
101101

102102
self._solutions = value
103103

104-
# * Set solutions per group
105-
if self._solutions.raw_arrays is not None:
106-
for e, group in enumerate(self.structural_frame.structural_groups):
107-
group.kriging_solution = RawArraysSolution( # ? Maybe I need to add more fields, but I am not sure yet
108-
scalar_field_matrix=self._solutions.raw_arrays.scalar_field_matrix[e],
109-
block_matrix=self._solutions.raw_arrays.block_matrix[e],
110-
)
111-
112104
# * Set solutions per element
113105
for e, element in enumerate(self.structural_frame.structural_elements[:-1]): # * Ignore basement
114106
element.scalar_field_at_interface = value.scalar_field_at_surface_points[e]

test/test_modules/_geophysics_TO_UPDATE/test_gravity.test_gravity.verify/2-layers.approved.txt

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -162,7 +162,7 @@
162162
"mesh_extraction_masking_options": 3,
163163
"mesh_extraction_fancy": true,
164164
"evaluation_chunk_size": 500000,
165-
"compute_scalar_gradient": true,
165+
"compute_scalar_gradient": false,
166166
"verbose": false
167167
},
168168
"debug": true,
Lines changed: 119 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,119 @@
1+
import numpy as np
2+
3+
import gempy as gp
4+
import time
5+
6+
PLOT = True
7+
DENSE_RESOLUTION = [20, 20, 20]
8+
9+
10+
def test_compute_time_dense_dense():
11+
geo_model: gp.data.GeoModel = _setup_model()
12+
geo_model.interpolation_options.evaluation_options.mesh_extraction = True
13+
14+
geo_model.grid.active_grids = gp.data.Grid.GridTypes.DENSE
15+
start_time = time.perf_counter()
16+
gp.compute_model(geo_model)
17+
print(f"Computing model on grids: {geo_model.grid.active_grids}")
18+
end_time = time.perf_counter()
19+
computation_time_dense_I = end_time - start_time
20+
21+
start_time = time.perf_counter()
22+
gp.compute_model(geo_model)
23+
print(f"Computing model on grids: {geo_model.grid.active_grids}")
24+
end_time = time.perf_counter()
25+
computation_time_dense_II = end_time - start_time
26+
27+
print(f"Computation only model dense grid 125*50*50: {computation_time_dense_I:.2f} seconds")
28+
print(f"Computation only model dense grid 125*50*50: {computation_time_dense_II:.2f} seconds")
29+
30+
# Assert that it is not too different
31+
assert abs(computation_time_dense_I - computation_time_dense_II) < 0.2
32+
33+
34+
def test_compute_time_topo_dense_grid():
35+
geo_model: gp.data.GeoModel = _setup_model()
36+
37+
geo_model.grid.active_grids = gp.data.Grid.GridTypes.DENSE
38+
start_time = time.perf_counter()
39+
gp.compute_model(geo_model)
40+
print(f"Computing model on grids: {geo_model.grid.active_grids}")
41+
end_time = time.perf_counter()
42+
computation_time_dense = end_time - start_time
43+
44+
# Compute a solution for the model
45+
geo_model.grid.active_grids = gp.data.Grid.GridTypes.TOPOGRAPHY
46+
start_time = time.perf_counter()
47+
gp.compute_model(geo_model)
48+
print(f"Computing model on grids: {geo_model.grid.active_grids}")
49+
end_time = time.perf_counter()
50+
computation_time_topo = end_time - start_time
51+
52+
# Recompute model as a new grid was added
53+
geo_model.grid.active_grids = gp.data.Grid.GridTypes.TOPOGRAPHY | gp.data.Grid.GridTypes.DENSE
54+
start_time = time.perf_counter()
55+
gp.compute_model(geo_model)
56+
print(f"Computing model on grids: {geo_model.grid.active_grids}")
57+
end_time = time.perf_counter()
58+
computation_time_topo_dense = end_time - start_time
59+
60+
print(f"Computation only model dense grid 125*50*50: {computation_time_dense:.2f} seconds")
61+
print(f"Computation time with topography 125*50: {computation_time_topo:.2f} seconds")
62+
print(f"Computation time with topography and dense grid 125*50*50: {computation_time_topo_dense:.2f} seconds")
63+
64+
# Assert that dense takes longer than topo and that the sum of both is close
65+
assert computation_time_dense > computation_time_topo
66+
assert computation_time_topo_dense > computation_time_dense
67+
68+
69+
def test_compute_time_custom_dense_grid():
70+
geo_model: gp.data.GeoModel = _setup_model()
71+
72+
# numpy array with random coordinates within the extent of the model
73+
custom_coordinates = np.random.uniform(
74+
low=geo_model.grid.extent[:3],
75+
high=geo_model.grid.extent[3:],
76+
size=(1000, 3)
77+
)
78+
79+
start_time = time.perf_counter()
80+
gp.compute_model_at(geo_model, custom_coordinates)
81+
end_time = time.perf_counter()
82+
computation_time_at = end_time - start_time
83+
84+
print(f"Computation compute_at with 1000 custom points: {computation_time_at:.2f} seconds")
85+
86+
87+
def _setup_model():
88+
# Define the path to data
89+
data_path = 'https://raw.githubusercontent.com/cgre-aachen/gempy_data/master/'
90+
path_to_data = data_path + "/data/input_data/jan_models/"
91+
# Create a GeoModel instance
92+
geo_model = gp.create_geomodel(
93+
project_name='EGU_example',
94+
extent=[0, 2500, 0, 1000, 0, 1000],
95+
resolution=DENSE_RESOLUTION,
96+
importer_helper=gp.data.ImporterHelper(
97+
path_to_orientations=path_to_data + "model7_orientations.csv",
98+
path_to_surface_points=path_to_data + "model7_surface_points.csv"
99+
)
100+
)
101+
# Map geological series to surfaces
102+
gp.map_stack_to_surfaces(
103+
gempy_model=geo_model,
104+
mapping_object={
105+
"Fault_Series" : ('fault'),
106+
"Strat_Series1": ('rock3'),
107+
"Strat_Series2": ('rock2', 'rock1'),
108+
}
109+
)
110+
# Define youngest structural group as fault
111+
gp.set_is_fault(geo_model, ["Fault_Series"])
112+
# Setting a randomly generated topography
113+
gp.set_topography_from_random(
114+
grid=geo_model.grid,
115+
fractal_dimension=2,
116+
d_z=np.array([700, 950]),
117+
topography_resolution=np.array([125, 50])
118+
)
119+
return geo_model

0 commit comments

Comments
 (0)