When adding a slippling fault perturbation to a dataset it creates some instabilities when evaluating the the model for its stratigraphy on a regular grid. When the fault slip is set to 0, the evaluation runs without problem. When the faultslip is set to 1500, then have a value error as described below.
The data perturbation basically moves all data points on on side of the fault by the faultslip along the slip vector.
data.loc[groupSix,"scalar_field"] = model.evaluate_feature_value('stratigraphy',xyzS,scale=True)#,scale=True)
data.loc[groupSix,"scalar_field"] = model.evaluate_feature_value('stratigraphy',xyzS,scale=True)#,scale=True)
---------------------------------------------------------------------------
ValueError Traceback (most recent call last)
<ipython-input-3-874a38232180> in <module>
----> 1 data.loc[groupSix,"scalar_field"] = model.evaluate_feature_value('stratigraphy',xyzS,scale=True)#,scale=True)
C:\ProgramData\Anaconda3\lib\site-packages\LoopStructural\modelling\core\geological_model.py in evaluate_feature_value(self, feature_name, xyz, scale)
1470 if scale:
1471 scaled_xyz = self.scale(xyz,inplace=False)
-> 1472 return feature.evaluate_value(scaled_xyz)
1473 else:
1474 return np.zeros(xyz.shape[0])
C:\ProgramData\Anaconda3\lib\site-packages\LoopStructural\modelling\features\geological_feature.py in evaluate_value(self, evaluation_points)
124 #if evaluation_points is not a numpy array try and convert
125 #otherwise error
--> 126 self.builder.up_to_date()
127 # check if the points are within the display region
128 v = np.zeros(evaluation_points.shape[0])
C:\ProgramData\Anaconda3\lib\site-packages\LoopStructural\modelling\features\geological_feature_builder.py in up_to_date(self)
88 #has anything changed in the builder since we built the feature? if so update
89 if self._up_to_date == False:
---> 90 self.update()
91 #check if the interpolator is up to date, if not solve
92 if self._interpolator.up_to_date == False:
C:\ProgramData\Anaconda3\lib\site-packages\LoopStructural\modelling\features\geological_feature_builder.py in update(self)
77
78 def update(self):
---> 79 self.build(**self.build_arguments)
80 @property
81 def name(self):
C:\ProgramData\Anaconda3\lib\site-packages\LoopStructural\modelling\features\geological_feature_builder.py in build(self, fold, fold_weights, data_region, **kwargs)
405
406
--> 407 self.add_data_to_interpolator(**kwargs)
408 if data_region is not None:
409 xyz = self.interpolator.get_data_locations()
C:\ProgramData\Anaconda3\lib\site-packages\LoopStructural\modelling\features\geological_feature_builder.py in add_data_to_interpolator(self, constrained, force_constrained, **kwargs)
172 for f in self.faults:
173 data.loc[:,xyz_names()] = f.apply_to_points(
--> 174 self.get_data_locations())
175 # Now check whether there are enough constraints for the
176 # interpolator to be able to solve
C:\ProgramData\Anaconda3\lib\site-packages\LoopStructural\modelling\fault\fault_segment.py in apply_to_points(self, points)
275 gz_future = executor.submit(self.faultframe.features[2].evaluate_value, newp)
276 gx = gx_future.result()
--> 277 gy = gy_future.result()
278 gz = gz_future.result()
279 else:
C:\ProgramData\Anaconda3\lib\concurrent\futures\_base.py in result(self, timeout)
426 raise CancelledError()
427 elif self._state == FINISHED:
--> 428 return self.__get_result()
429
430 self._condition.wait(timeout)
C:\ProgramData\Anaconda3\lib\concurrent\futures\_base.py in __get_result(self)
382 def __get_result(self):
383 if self._exception:
--> 384 raise self._exception
385 else:
386 return self._result
C:\ProgramData\Anaconda3\lib\concurrent\futures\thread.py in run(self)
55
56 try:
---> 57 result = self.fn(*self.args, **self.kwargs)
58 except BaseException as exc:
59 self.future.set_exception(exc)
C:\ProgramData\Anaconda3\lib\site-packages\LoopStructural\modelling\features\geological_feature.py in evaluate_value(self, evaluation_points)
124 #if evaluation_points is not a numpy array try and convert
125 #otherwise error
--> 126 self.builder.up_to_date()
127 # check if the points are within the display region
128 v = np.zeros(evaluation_points.shape[0])
C:\ProgramData\Anaconda3\lib\site-packages\LoopStructural\modelling\features\geological_feature_builder.py in up_to_date(self)
88 #has anything changed in the builder since we built the feature? if so update
89 if self._up_to_date == False:
---> 90 self.update()
91 #check if the interpolator is up to date, if not solve
92 if self._interpolator.up_to_date == False:
C:\ProgramData\Anaconda3\lib\site-packages\LoopStructural\modelling\features\geological_feature_builder.py in update(self)
77
78 def update(self):
---> 79 self.build(**self.build_arguments)
80 @property
81 def name(self):
C:\ProgramData\Anaconda3\lib\site-packages\LoopStructural\modelling\features\geological_feature_builder.py in build(self, fold, fold_weights, data_region, **kwargs)
436 # try adding very small cg
437 kwargs['cgw'] = 0.0
--> 438 self.install_gradient_constraint()
439 self.install_equality_constraints()
440 self.interpolator.setup_interpolator(**kwargs)
C:\ProgramData\Anaconda3\lib\site-packages\LoopStructural\modelling\features\geological_feature_builder.py in install_gradient_constraint(self)
249 for g in self._orthogonal_features.values():
250 feature,w,region,step,B = g
--> 251 vector = feature.evaluate_gradient(self.interpolator.support.barycentre())
252 vector /= np.linalg.norm(vector,axis=1)[:,None]
253 element_idx = np.arange(self.interpolator.support.n_elements)
C:\ProgramData\Anaconda3\lib\site-packages\LoopStructural\modelling\features\geological_feature.py in evaluate_gradient(self, evaluation_points)
158
159 """
--> 160 self.builder.up_to_date()
161 v = np.zeros(evaluation_points.shape)
162 v[:] = np.nan
C:\ProgramData\Anaconda3\lib\site-packages\LoopStructural\modelling\features\geological_feature_builder.py in up_to_date(self)
88 #has anything changed in the builder since we built the feature? if so update
89 if self._up_to_date == False:
---> 90 self.update()
91 #check if the interpolator is up to date, if not solve
92 if self._interpolator.up_to_date == False:
C:\ProgramData\Anaconda3\lib\site-packages\LoopStructural\modelling\features\geological_feature_builder.py in update(self)
77
78 def update(self):
---> 79 self.build(**self.build_arguments)
80 @property
81 def name(self):
C:\ProgramData\Anaconda3\lib\site-packages\LoopStructural\modelling\features\geological_feature_builder.py in build(self, fold, fold_weights, data_region, **kwargs)
438 self.install_gradient_constraint()
439 self.install_equality_constraints()
--> 440 self.interpolator.setup_interpolator(**kwargs)
441 self.interpolator.solve_system(**kwargs)
442 self._up_to_date = True
C:\ProgramData\Anaconda3\lib\site-packages\LoopStructural\interpolators\geological_interpolator.py in setup_interpolator(self, **kwargs)
171 Runs all of the required setting up stuff
172 """
--> 173 self._setup_interpolator(**kwargs)
174
175 def solve_system(self, **kwargs):
C:\ProgramData\Anaconda3\lib\site-packages\LoopStructural\interpolators\piecewiselinear_interpolator.py in _setup_interpolator(self, **kwargs)
85 self.n_t, self.n_i, self.propertyname))
86 self.add_gradient_ctr_pts(self.interpolation_weights['gpw'])
---> 87 self.add_norm_ctr_pts(self.interpolation_weights['npw'])
88 self.add_ctr_pts(self.interpolation_weights['cpw'])
89 self.add_tangent_ctr_pts(self.interpolation_weights['tpw'])
C:\ProgramData\Anaconda3\lib\site-packages\LoopStructural\interpolators\piecewiselinear_interpolator.py in add_norm_ctr_pts(self, w)
312 vol[outside, None],
313 idc[outside],
--> 314 name='norm')
315
316 def add_ctr_pts(self, w=1.0): # for now weight all value points the same
C:\ProgramData\Anaconda3\lib\site-packages\LoopStructural\interpolators\discrete_interpolator.py in add_constraints_to_least_squares(self, A, B, idc, name)
181
182 self.constraints[name]['A'] = np.vstack([self.constraints[name]['A'],A])
--> 183 self.constraints[name]['B'] = np.hstack([self.constraints[name]['B'], B])
184 self.constraints[name]['idc'] = np.vstack([self.constraints[name]['idc'],
185 idc])
<__array_function__ internals> in hstack(*args, **kwargs)
C:\ProgramData\Anaconda3\lib\site-packages\numpy\core\shape_base.py in hstack(tup)
341 # As a special case, dimension 0 of 1-dimensional arrays is "horizontal"
342 if arrs and arrs[0].ndim == 1:
--> 343 return _nx.concatenate(arrs, 0)
344 else:
345 return _nx.concatenate(arrs, 1)
<__array_function__ internals> in concatenate(*args, **kwargs)
ValueError: all the input arrays must have same number of dimensions, but the array at index 0 has 1 dimension(s) and the array at index 1 has 2 dimension(s)