diff --git a/ogstools/propertylib/mesh_dependent.py b/ogstools/propertylib/mesh_dependent.py index 373fa25a0245721ba4cd09f1996628e67d1dd2a8..ca5e0d2fd8acfe4a34d5af5d59e284ebcdccf8ce 100644 --- a/ogstools/propertylib/mesh_dependent.py +++ b/ogstools/propertylib/mesh_dependent.py @@ -48,11 +48,16 @@ def depth(mesh: pv.UnstructuredGrid, use_coords: bool = False) -> np.ndarray: ] edge_centers = edges.cell_centers().points adj_cells = [mesh.find_closest_cell(point) for point in edge_centers] - adj_cell_centers = mesh.extract_cells(adj_cells).cell_centers().points - are_above = ( - edge_centers[..., vertical_dim] - > adj_cell_centers[..., vertical_dim] + adj_centers = np.vstack( + [ + mesh.extract_cells(adj_cell).cell_centers().points + for adj_cell in adj_cells + ] ) + are_above = [ + edge_center[vertical_dim] > adj_center[vertical_dim] + for edge_center, adj_center in zip(edge_centers, adj_centers) + ] are_non_vertical = np.asarray(edge_horizontal_extent) > 1e-12 top_cells = are_above & are_non_vertical top = edges.extract_cells(top_cells) @@ -110,8 +115,10 @@ def fluid_pressure_criterion( """ Qty = u_reg.Quantity - sigma = -Qty(mesh[mesh_property.data_name], mesh_property.data_unit) - return p_fluid(mesh) - eigenvalues(sigma)[..., 0] + sigma = mesh[mesh_property.data_name] + return p_fluid(mesh) - Qty( + eigenvalues(-sigma)[..., 0], mesh_property.data_unit + ) def dilatancy_critescu(