Merge branch 'master' of github.com:Loop3D/LoopStructural

.github/workflows/release-please.yml

@@ -92,9 +92,11 @@ jobs:
           folder: docs/build/html # The folder the action should deploy.
 
   conda-deploy:
-    name: Building conda package
-    # needs: ["documentation-test", "continuous-integration"]
-    runs-on: ubuntu-latest
+    name: Building conda package for python ${{ matrix.os }})
+    needs: ["documentation-test", "continuous-integration", "release-please"]
+    runs-on: ${{matrix.os}}
+    if: github.ref == 'refs/heads/master'
+
     strategy:
       fail-fast: false
       matrix:
@@ -114,8 +116,8 @@ jobs:
           ANACONDA_API_TOKEN: ${{ secrets.ANACONDA_TOKEN  }}
         shell: bash -l {0}
         run: |
-          conda install -c conda-forge conda-build scikit-build numpy cython anaconda-client  -y
-          conda build -c anaconda -c conda-forge -c loop3d --output-folder conda conda --python ${{ matrix.python-version }}
+          conda install -c conda-forge conda-build scikit-build numpy cython anaconda-client -y
+          conda build -c anaconda -c conda-forge -c loop3d --output-folder conda conda  --python ${{matrix.python-version}}
           conda convert -p all conda/linux-64/*.tar.bz2 -f -o conda
       - name: upload artifacts
         uses: actions/upload-artifact@v4
@@ -124,7 +126,7 @@ jobs:
           path: conda
 
   make_sdist:
-    needs: ["documentation-test", "continuous-integration"]
+    needs: ["documentation-test", "continuous-integration", "release-please"]
     name: Make SDist
     runs-on: ubuntu-latest
 

CHANGELOG.md

@@ -1,5 +1,14 @@
 # Changelog
 
+## [1.6.2](https://github.com/Loop3D/LoopStructural/compare/v1.6.1...v1.6.2) (2024-08-06)
+
+
+### Bug Fixes
+
+* extra import ([7d10434](https://github.com/Loop3D/LoopStructural/commit/7d10434eb11631fa501275c14d617ed014f092a7))
+* tuple to Tuple ([e567810](https://github.com/Loop3D/LoopStructural/commit/e567810e4fafd36da7ccf5696dd9245a904d4462))
+* updating type hints ([a064224](https://github.com/Loop3D/LoopStructural/commit/a0642243fac0bd7e90f28957b95d68e31bac0af7))
+
 ## [1.6.1](https://github.com/Loop3D/LoopStructural/compare/v1.6.0...v1.6.1) (2024-08-05)
 
 

LoopStructural/datatypes/_bounding_box.py

@@ -4,6 +4,7 @@
 from LoopStructural.utils import rng
 from LoopStructural.datatypes._structured_grid import StructuredGrid
 import numpy as np
+import copy
 
 from LoopStructural.utils.logging import getLogger
 
@@ -412,11 +413,16 @@ def vtk(self):
     def structured_grid(
         self, cell_data: Dict[str, np.ndarray] = {}, vertex_data={}, name: str = "bounding_box"
     ):
+        # python is passing a reference to the cell_data, vertex_data dicts so we need to
+        # copy them to make sure that different instances of StructuredGrid are not sharing the same
+        # underlying objects
+        _cell_data = copy.deepcopy(cell_data)
+        _vertex_data = copy.deepcopy(vertex_data)
         return StructuredGrid(
             origin=self.global_origin,
             step_vector=self.step_vector,
             nsteps=self.nsteps,
-            cell_properties=cell_data,
-            properties=vertex_data,
+            cell_properties=_cell_data,
+            properties=_vertex_data,
             name=name,
         )

LoopStructural/interpolators/_discrete_interpolator.py

@@ -62,7 +62,6 @@ def __init__(self, support, data={}, c=None, up_to_date=False):
         self.interpolation_weights = {}
         logger.info("Creating discrete interpolator with {} degrees of freedom".format(self.nx))
         self.type = InterpolatorType.BASE_DISCRETE
-        self.c = np.zeros(self.support.n_nodes)
 
     @property
     def nx(self) -> int:
@@ -133,14 +132,36 @@ def set_interpolation_weights(self, weights):
             self.up_to_date = False
             self.interpolation_weights[key] = weights[key]
 
+    def _pre_solve(self):
+        """
+        Pre solve function to be run before solving the interpolation
+        """
+        self.c = np.zeros(self.support.n_nodes)
+        self.c[:] = np.nan
+        return True
+
+    def _post_solve(self):
+        """Post solve function(s) to be run after the solver has been called"""
+        self.clear_constraints()
+        return True
+
+    def clear_constraints(self):
+        """
+        Clear the constraints from the interpolator, this makes sure we are not storing
+        the constraints after the solver has been run
+        """
+        self.constraints = {}
+        self.ineq_constraints = {}
+        self.equal_constraints = {}
+
     def reset(self):
         """
         Reset the interpolation constraints
 
         """
         self.constraints = {}
         self.c_ = 0
-        logger.debug("Resetting interpolation constraints")
+        logger.info("Resetting interpolation constraints")
 
     def add_constraints_to_least_squares(self, A, B, idc, w=1.0, name="undefined"):
         """
@@ -511,14 +532,15 @@ def build_inequality_matrix(self):
             mats.append(c['matrix'])
             bounds.append(c['bounds'])
         if len(mats) == 0:
-            return None, None
+            return sparse.csr_matrix((0, self.nx), dtype=float), np.zeros((0, 3))
         Q = sparse.vstack(mats)
         bounds = np.vstack(bounds)
         return Q, bounds
 
     def solve_system(
         self,
         solver: Optional[Union[Callable[[sparse.csr_matrix, np.ndarray], np.ndarray], str]] = None,
+        tol: Optional[float] = None,
         solver_kwargs: dict = {},
     ) -> bool:
         """
@@ -540,37 +562,46 @@ def solve_system(
 
         """
         starttime = time()
-        self.c = np.zeros(self.support.n_nodes)
-        self.c[:] = np.nan
+        if not self._pre_solve():
+            raise ValueError("Pre solve failed")
+
         A, b = self.build_matrix()
         Q, bounds = self.build_inequality_matrix()
         if callable(solver):
             logger.warning('Using custom solver')
             self.c = solver(A.tocsr(), b)
             self.up_to_date = True
-
-            return True
-        ## solve with lsmr
-        if isinstance(solver, str):
+        elif isinstance(solver, str) or solver is None:
             if solver not in ['cg', 'lsmr', 'admm']:
                 logger.warning(
                     f'Unknown solver {solver} using cg. \n Available solvers are cg and lsmr or a custom solver as a callable function'
                 )
                 solver = 'cg'
         if solver == 'cg':
             logger.info("Solving using cg")
-            ATA = A.T.dot(A)
-            ATB = A.T.dot(b)
-            res = sparse.linalg.cg(ATA, ATB, **solver_kwargs)
+            if 'atol' not in solver_kwargs or 'rtol' not in solver_kwargs:
+                if tol is not None:
+                    solver_kwargs['atol'] = tol
+
+            logger.info(f"Solver kwargs: {solver_kwargs}")
+
+            res = sparse.linalg.cg(A.T @ A, A.T @ b, **solver_kwargs)
             if res[1] > 0:
                 logger.warning(
                     f'CG reached iteration limit ({res[1]})and did not converge, check input data. Setting solution to last iteration'
                 )
             self.c = res[0]
             self.up_to_date = True
-            return True
+
         elif solver == 'lsmr':
             logger.info("Solving using lsmr")
+            if 'atol' not in solver_kwargs:
+                if tol is not None:
+                    solver_kwargs['atol'] = tol
+            if 'btol' not in solver_kwargs:
+                if tol is not None:
+                    solver_kwargs['btol'] = tol
+            logger.info(f"Solver kwargs: {solver_kwargs}")
             res = sparse.linalg.lsmr(A, b, **solver_kwargs)
             if res[1] == 1 or res[1] == 4 or res[1] == 2 or res[1] == 5:
                 self.c = res[0]
@@ -585,8 +616,7 @@ def solve_system(
                 )
                 self.c = res[0]
             self.up_to_date = True
-            logger.info("Interpolation took %f seconds" % (time() - starttime))
-            return True
+
         elif solver == 'admm':
             logger.info("Solving using admm")
 
@@ -601,28 +631,35 @@ def solve_system(
 
             try:
                 from loopsolver import admm_solve
+
+                try:
+                    linsys_solver = solver_kwargs.pop('linsys_solver', 'lsmr')
+                    res = admm_solve(
+                        A,
+                        b,
+                        Q,
+                        bounds,
+                        x0=x0,
+                        admm_weight=solver_kwargs.pop('admm_weight', 0.01),
+                        nmajor=solver_kwargs.pop('nmajor', 200),
+                        linsys_solver_kwargs=solver_kwargs,
+                        linsys_solver=linsys_solver,
+                    )
+                    self.c = res
+                    self.up_to_date = True
+                except ValueError as e:
+                    logger.error(f"ADMM solver failed: {e}")
+                    self.up_to_date = False
             except ImportError:
                 logger.warning(
                     "Cannot import admm solver. Please install loopsolver or use lsmr or cg"
                 )
-                return False
-            try:
-                res = admm_solve(
-                    A,
-                    b,
-                    Q,
-                    bounds,
-                    x0=x0,
-                    admm_weight=solver_kwargs.pop('admm_weight', 0.01),
-                    nmajor=solver_kwargs.pop('nmajor', 200),
-                    linsys_solver_kwargs=solver_kwargs,
-                )
-                self.c = res
-                self.up_to_date = True
-            except ValueError as e:
-                logger.error(f"ADMM solver failed: {e}")
-                return False
-        return False
+                self.up_to_date = False
+        else:
+            logger.error(f"Unknown solver {solver}")
+            self.up_to_date = False
+        # self._post_solve()
+        return self.up_to_date
 
     def update(self) -> bool:
         """
@@ -641,7 +678,8 @@ def update(self) -> bool:
             return False
         if not self.up_to_date:
             self.setup_interpolator()
-            return self.solve_system(self.solver)
+            self.up_to_date = self.solve_system(self.solver)
+            return self.up_to_date
 
     def evaluate_value(self, locations: np.ndarray) -> np.ndarray:
         """Evaluate the value of the interpolator at location

LoopStructural/modelling/features/builders/_base_builder.py

@@ -46,6 +46,7 @@ def build_arguments(self, build_arguments):
         # self._build_arguments = {}
         for k, i in build_arguments.items():
             if i != self._build_arguments.get(k, None):
+                logger.info(f"Setting {k} to {i} for {self.name}")
                 self._build_arguments[k] = i
                 ## if build_arguments change then flag to reinterpolate
                 self._up_to_date = False
@@ -107,5 +108,6 @@ def add_fault(self, fault):
         -------
 
         """
+        logger.info(f'Adding fault {fault.name} to {self.name}')
         self._up_to_date = False
         self.faults.append(fault)

LoopStructural/modelling/features/builders/_geological_feature_builder.py

@@ -109,6 +109,7 @@ def interpolation_region(self, interpolation_region):
         else:
             self._interpolation_region = RegionEverywhere()
             self._interpolator.set_region(region=self._interpolation_region)
+        logger.info(f'Setting interpolation region {self.name}')
         self._up_to_date = False
 
     def add_data_from_data_frame(self, data_frame, overwrite=False):
@@ -156,6 +157,8 @@ def add_orthogonal_feature(self, feature, w=1.0, region=None, step=1, B=0):
             logger.error("Cannot cast {} as integer, setting step to 1".format(step))
             step = 1
         self._orthogonal_features[feature.name] = [feature, w, region, step, B]
+
+        logger.info(f"Adding orthogonal constraint {feature.name} to {self.name}")
         self._up_to_date = False
 
     def add_data_to_interpolator(self, constrained=False, force_constrained=False, **kwargs):
@@ -321,6 +324,7 @@ def install_gradient_constraint(self):
 
     def add_equality_constraints(self, feature, region, scalefactor=1.0):
         self._equality_constraints[feature.name] = [feature, region, scalefactor]
+        logger.info(f'Adding equality constraints to {self.name}')
         self._up_to_date = False
 
     def install_equality_constraints(self):
@@ -536,7 +540,9 @@ def build(self, fold=None, fold_weights={}, data_region=None, **kwargs):
         logger.info(f'running interpolation for {self.name}')
 
         self.interpolator.solve_system(
-            solver=kwargs.get('solver', None), solver_kwargs=kwargs.get('solver_kwargs', {})
+            solver=kwargs.get('solver', None),
+            tol=kwargs.get('tol', None),
+            solver_kwargs=kwargs.get('solver_kwargs', {}),
         )
         logger.info(f'Finished building  {self.name}')
         self._up_to_date = True

LoopStructural/modelling/features/fault/_fault_segment.py

@@ -282,6 +282,7 @@ def apply_to_points(self, points, reverse=False):
         -------
 
         """
+        logger.info(f'Applying fault {self.name} to points {points.shape}')
         steps = self.steps
         newp = np.copy(points).astype(float)
         # evaluate fault function for all points

LoopStructural/modelling/input/process_data.py

@@ -441,6 +441,12 @@ def _stratigraphic_value(self):
             for g in reversed(sg):
                 if g not in self.thicknesses:
                     logger.warning(f"No thicknesses for {g}")
+                    stratigraphic_value[g] = np.nan
+                if self.thicknesses[g] <= 0:
+                    logger.error(
+                        f"Thickness for {g} is less than or equal to 0\n Update the thickness value for {g} before continuing"
+                    )
+
                     stratigraphic_value[g] = np.nan
                 else:
                     stratigraphic_value[g] = value

LoopStructural/version.py

@@ -1 +1 @@
-__version__ = "1.6.1"
+__version__ = "1.6.2"