From 19a58a213c1b1459a779586766bf8ac6acc4d170 Mon Sep 17 00:00:00 2001 From: Piotr Bartman Date: Thu, 25 Jul 2024 07:20:03 +0200 Subject: [PATCH] alder version of meshzoo --- conmech/dynamics/contact/contact_law.py | 4 +- .../contact/damped_normal_compliance.py | 8 +--- .../factory/_abstract_dynamics_factory.py | 8 +--- .../dynamics/factory/_dynamics_factory_2d.py | 21 ++------ .../dynamics/factory/_dynamics_factory_3d.py | 13 ++--- .../factory/dynamics_factory_method.py | 16 ++----- conmech/dynamics/statement.py | 12 +---- conmech/helpers/nph.py | 3 +- conmech/mesh/boundaries.py | 26 +++------- conmech/mesh/boundaries_factory.py | 44 +++++------------ conmech/mesh/mesh.py | 10 +--- conmech/mesh/zoo/barboteu_2008.py | 4 +- conmech/mesh/zoo/cross_for_tests.py | 19 ++------ conmech/mesh/zoo/rectangle.py | 8 ++-- conmech/mesh/zoo/sofochbar_2023.py | 3 +- conmech/plotting/drawer.py | 8 +--- conmech/plotting/membrane.py | 8 +--- conmech/scenarios/problems.py | 4 +- conmech/scene/body_forces.py | 12 ++--- conmech/simulations/problem_solver.py | 43 +++++------------ conmech/solvers/optimization/optimization.py | 8 +--- .../solvers/optimization/schur_complement.py | 24 +++------- conmech/solvers/solver.py | 8 +--- conmech/solvers/solver_methods.py | 21 ++------ conmech/state/state.py | 22 +++------ examples/BOSK_2024_example_1.py | 4 +- examples/BOSK_2024_example_2.py | 18 ++----- examples/BOSK_plot_postprocessing.py | 12 ++--- examples/BOST_2024.py | 19 ++------ examples/Jureczka_Ochal_Bartman_2023.py | 12 ++--- examples/Jureczka_and_Ochal_2019.py | 4 +- examples/Sofonea_Ochal_Bartman_2023.py | 48 +++++-------------- examples/example_dynamic_membrane.py | 4 +- examples/example_piezo_quasistatic.py | 12 ++--- examples/example_piezoelectric_dynamic.py | 8 +--- examples/example_static.py | 4 +- examples/example_temperature_dynamic.py | 16 ++----- examples/utils.py | 4 +- tests/test_conmech/regression/std_boundary.py | 4 +- .../regression/test_density_update.py | 12 ++--- .../regression/test_dynamic_membrane.py | 4 +- .../regression/test_membrane_DNC.py | 4 +- .../regression/test_membrane_boundary_DNC.py | 4 +- .../test_piezoelectric_quasistatic.py | 24 +++------- tests/test_conmech/regression/test_poisson.py | 8 +--- .../regression/test_temperature_dynamic.py | 24 +++------- tests/test_conmech/smoke/test_examples.py | 12 ++--- tests/test_conmech/unit_test/test_matrices.py | 9 +--- tests/test_conmech/unit_test/test_mesh.py | 16 ++----- .../unit_test/test_setting_mesh.py | 4 +- 50 files changed, 163 insertions(+), 484 deletions(-) diff --git a/conmech/dynamics/contact/contact_law.py b/conmech/dynamics/contact/contact_law.py index 27ca0a42..3202e96d 100644 --- a/conmech/dynamics/contact/contact_law.py +++ b/conmech/dynamics/contact/contact_law.py @@ -37,9 +37,7 @@ def normal_bound(var_nu: float, static_displacement_nu: float, dt: float) -> flo return 1.0 @staticmethod - def tangential_bound( - var_nu: float, static_displacement_nu: float, dt: float - ) -> float: + def tangential_bound(var_nu: float, static_displacement_nu: float, dt: float) -> float: """ Friction bound diff --git a/conmech/dynamics/contact/damped_normal_compliance.py b/conmech/dynamics/contact/damped_normal_compliance.py index ee58f645..995308ee 100644 --- a/conmech/dynamics/contact/damped_normal_compliance.py +++ b/conmech/dynamics/contact/damped_normal_compliance.py @@ -20,16 +20,12 @@ from conmech.dynamics.contact.interior_contact_law import InteriorContactLaw -def make_damped_norm_compl( - obstacle_level: float, kappa: float, beta: float, interior=False -): +def make_damped_norm_compl(obstacle_level: float, kappa: float, beta: float, interior=False): superclass = InteriorContactLaw if interior else PotentialOfContactLaw class DampedNormalCompliance(superclass): @staticmethod - def normal_bound( - var_nu: float, static_displacement_nu: float, dt: float - ) -> float: + def normal_bound(var_nu: float, static_displacement_nu: float, dt: float) -> float: """ Since multiply by var_nu """ diff --git a/conmech/dynamics/factory/_abstract_dynamics_factory.py b/conmech/dynamics/factory/_abstract_dynamics_factory.py index 9a10c4f8..c9e61e65 100644 --- a/conmech/dynamics/factory/_abstract_dynamics_factory.py +++ b/conmech/dynamics/factory/_abstract_dynamics_factory.py @@ -23,14 +23,10 @@ def get_relaxation_tensor(self, W: np.ndarray, coeff: np.ndarray): def calculate_acceleration(self, U: np.ndarray, density: float) -> np.ndarray: raise NotImplementedError() - def calculate_thermal_expansion( - self, V: np.ndarray, coeff: np.ndarray - ) -> np.ndarray: + def calculate_thermal_expansion(self, V: np.ndarray, coeff: np.ndarray) -> np.ndarray: raise NotImplementedError() - def calculate_thermal_conductivity( - self, W: np.ndarray, coeff: np.ndarray - ) -> np.ndarray: + def calculate_thermal_conductivity(self, W: np.ndarray, coeff: np.ndarray) -> np.ndarray: raise NotImplementedError() def get_piezoelectric_tensor(self, W: np.ndarray, coeff: np.ndarray) -> np.ndarray: diff --git a/conmech/dynamics/factory/_dynamics_factory_2d.py b/conmech/dynamics/factory/_dynamics_factory_2d.py index fee35998..b1767439 100644 --- a/conmech/dynamics/factory/_dynamics_factory_2d.py +++ b/conmech/dynamics/factory/_dynamics_factory_2d.py @@ -43,9 +43,7 @@ def get_edges_features_matrix_numba(elements, nodes): ) en0[:] = element_nodes[:, 0] en1[:] = element_nodes[:, 1] - int_xy = ( - en0 @ np.array([[2.0, 1.0, 1.0], [1.0, 2.0, 1.0], [1.0, 1.0, 2.0]]) @ en1 - ) + int_xy = en0 @ np.array([[2.0, 1.0, 1.0], [1.0, 2.0, 1.0], [1.0, 1.0, 2.0]]) @ en1 int_matrix = np.array( [ [ @@ -103,18 +101,11 @@ def get_edges_features_matrix_numba(elements, nodes): v = [INT_PH * j_d_phi for j_d_phi in j_d_phi_vec] - w = [ - [i_d_phi * j_d_phi for j_d_phi in j_d_phi_vec] - for i_d_phi in i_d_phi_vec - ] + w = [[i_d_phi * j_d_phi for j_d_phi in j_d_phi_vec] for i_d_phi in i_d_phi_vec] - local_stifness_matrices[:, :, element_index, i, j] = ( - element_volume * np.asarray(w) - ) + local_stifness_matrices[:, :, element_index, i, j] = element_volume * np.asarray(w) - edges_features_matrix[ - :, element[i], element[j] - ] += element_volume * np.array( + edges_features_matrix[:, element[i], element[j]] += element_volume * np.array( [ volume_at_nodes, u, @@ -169,9 +160,7 @@ def get_edges_features_matrix(self, elements, nodes): def dimension(self) -> int: return DIMENSION - def calculate_constitutive_matrices( - self, W: np.ndarray, mu: float, lambda_: float - ) -> SM2: + def calculate_constitutive_matrices(self, W: np.ndarray, mu: float, lambda_: float) -> SM2: A_11 = (2 * mu + lambda_) * W[0, 0] + mu * W[1, 1] A_12 = mu * W[1, 0] + lambda_ * W[0, 1] A_21 = lambda_ * W[1, 0] + mu * W[0, 1] diff --git a/conmech/dynamics/factory/_dynamics_factory_3d.py b/conmech/dynamics/factory/_dynamics_factory_3d.py index 77280b5f..65411b32 100644 --- a/conmech/dynamics/factory/_dynamics_factory_3d.py +++ b/conmech/dynamics/factory/_dynamics_factory_3d.py @@ -48,14 +48,9 @@ def get_edges_features_matrix_numba(elements, nodes): v = [INT_PH * j_d_phi for j_d_phi in j_d_phi_vec] - w = [ - [i_d_phi * j_d_phi for j_d_phi in j_d_phi_vec] - for i_d_phi in i_d_phi_vec - ] + w = [[i_d_phi * j_d_phi for j_d_phi in j_d_phi_vec] for i_d_phi in i_d_phi_vec] - edges_features_matrix[ - :, element[i], element[j] - ] += element_volume * np.array( + edges_features_matrix[:, element[i], element[j]] += element_volume * np.array( [ volume_at_nodes, u, @@ -81,9 +76,7 @@ def get_edges_features_matrix_numba(elements, nodes): @numba.njit def get_integral_parts_numba(element_nodes, element_index): x_i = element_nodes[element_index] - x_j1, x_j2, x_j3 = list( - element_nodes[np.arange(ELEMENT_NODES_COUNT) != element_index] - ) + x_j1, x_j2, x_j3 = list(element_nodes[np.arange(ELEMENT_NODES_COUNT) != element_index]) dm = denominator_numba(x_i, x_j1, x_j2, x_j3) element_volume = np.abs(dm) / VOLUME_DIVIDER diff --git a/conmech/dynamics/factory/dynamics_factory_method.py b/conmech/dynamics/factory/dynamics_factory_method.py index fac5b398..6247d642 100644 --- a/conmech/dynamics/factory/dynamics_factory_method.py +++ b/conmech/dynamics/factory/dynamics_factory_method.py @@ -60,23 +60,17 @@ def get_dynamics(elements: np.ndarray, body_prop: BodyProperties, U, V, W): poisson_operator = None if isinstance(body_prop, ElasticProperties): - elasticity = factory.calculate_constitutive_matrices( - W, body_prop.mu, body_prop.lambda_ - ) + elasticity = factory.calculate_constitutive_matrices(W, body_prop.mu, body_prop.lambda_) else: elasticity = None if isinstance(body_prop, ViscoelasticProperties): - viscosity = factory.calculate_constitutive_matrices( - W, body_prop.theta, body_prop.zeta - ) + viscosity = factory.calculate_constitutive_matrices(W, body_prop.theta, body_prop.zeta) else: viscosity = None if isinstance(body_prop, TemperatureBodyProperties): - thermal_expansion = factory.calculate_thermal_expansion( - V, body_prop.thermal_expansion - ) + thermal_expansion = factory.calculate_thermal_expansion(V, body_prop.thermal_expansion) thermal_conductivity = factory.calculate_thermal_conductivity( W, body_prop.thermal_conductivity ) @@ -85,9 +79,7 @@ def get_dynamics(elements: np.ndarray, body_prop: BodyProperties, U, V, W): thermal_conductivity = None if isinstance(body_prop, PiezoelectricBodyProperties): - piezoelectricity = factory.get_piezoelectric_tensor( - W, body_prop.piezoelectricity - ) + piezoelectricity = factory.get_piezoelectric_tensor(W, body_prop.piezoelectricity) permittivity = factory.get_permittivity_tensor(W, body_prop.permittivity) else: piezoelectricity = None diff --git a/conmech/dynamics/statement.py b/conmech/dynamics/statement.py index 3a6e890a..78a11eaf 100644 --- a/conmech/dynamics/statement.py +++ b/conmech/dynamics/statement.py @@ -99,11 +99,7 @@ def update_right_hand_side(self, var): A = -1 * self.body.dynamics.poisson_operator @ var.displacement[:ind] - A += ( - (1 / var.time_step) - * self.body.dynamics.acceleration_operator.SM1 - @ var.velocity[:ind] - ) + A += (1 / var.time_step) * self.body.dynamics.acceleration_operator.SM1 @ var.velocity[:ind] self.right_hand_side = self.body.dynamics.force.integrate(time=var.time) + A @@ -183,11 +179,7 @@ def update_right_hand_side(self, var): A = -1 * self.body.dynamics.elasticity @ var.displacement - A += ( - (1 / var.time_step) - * self.body.dynamics.acceleration_operator - @ var.velocity - ) + A += (1 / var.time_step) * self.body.dynamics.acceleration_operator @ var.velocity self.right_hand_side = self.body.dynamics.force.integrate(time=var.time) + A diff --git a/conmech/helpers/nph.py b/conmech/helpers/nph.py index 47eef0ad..18741eff 100644 --- a/conmech/helpers/nph.py +++ b/conmech/helpers/nph.py @@ -39,6 +39,5 @@ def get_node_index_numba(node, nodes): @numba.njit(inline="always") def length(edge, nodes): return np.sqrt( - (nodes[edge[0]][0] - nodes[edge[1]][0]) ** 2 - + (nodes[edge[0]][1] - nodes[edge[1]][1]) ** 2 + (nodes[edge[0]][0] - nodes[edge[1]][0]) ** 2 + (nodes[edge[0]][1] - nodes[edge[1]][1]) ** 2 ) diff --git a/conmech/mesh/boundaries.py b/conmech/mesh/boundaries.py index 022be3d2..4bd8041a 100644 --- a/conmech/mesh/boundaries.py +++ b/conmech/mesh/boundaries.py @@ -61,9 +61,7 @@ def get_boundary_surfaces_normals(nodes, boundary_surfaces, boundary_internal_in internal_nodes = nodes[boundary_internal_indices] external_orientation = (-1) * np.sign( - nph.elementwise_dot( - internal_nodes - tail_nodes, unoriented_normals, keepdims=True - ) + nph.elementwise_dot(internal_nodes - tail_nodes, unoriented_normals, keepdims=True) ) return unoriented_normals * external_orientation @@ -86,19 +84,13 @@ def __init__(self, nodes, boundary_internal_indices: np.ndarray, **kwargs): @property def boundary_surfaces(self): return np.unique( - np.vstack( - (self.contact_boundary, self.neumann_boundary, self.dirichlet_boundary) - ), + np.vstack((self.contact_boundary, self.neumann_boundary, self.dirichlet_boundary)), axis=1, ) @property def boundary_nodes_count(self): - return ( - self.contact_nodes_count - + self.neumann_nodes_count - + self.dirichlet_nodes_count - ) + return self.contact_nodes_count + self.neumann_nodes_count + self.dirichlet_nodes_count @property def boundary_indices(self): @@ -139,16 +131,12 @@ def get_all_boundary_indices(self, boundary, total_node_count, dimension): condition_start = 0 condition_stop = condition_start + i.stop - i.start yield ( - slice( - i.start + d * total_node_count, i.stop + d * total_node_count - ), + slice(i.start + d * total_node_count, i.stop + d * total_node_count), slice(condition_start, condition_stop), ) else: # assuming node_indices are sorted - discontinuities = np.concatenate( - ([0], np.nonzero(np.diff(i) - 1)[0] + 1, [len(i)]) - ) + discontinuities = np.concatenate(([0], np.nonzero(np.diff(i) - 1)[0] + 1, [len(i)])) starts = i[discontinuities[:-1]] stops = i[(discontinuities - 1)[1:]] + 1 for d in range(dimension): @@ -156,9 +144,7 @@ def get_all_boundary_indices(self, boundary, total_node_count, dimension): for start, stop in zip(starts, stops): condition_stop = condition_start + stop - start yield ( - slice( - start + d * total_node_count, stop + d * total_node_count - ), + slice(start + d * total_node_count, stop + d * total_node_count), slice(condition_start, condition_stop), ) condition_start = condition_stop diff --git a/conmech/mesh/boundaries_factory.py b/conmech/mesh/boundaries_factory.py index 4b8eb545..811c5b64 100644 --- a/conmech/mesh/boundaries_factory.py +++ b/conmech/mesh/boundaries_factory.py @@ -21,9 +21,7 @@ def identify_surfaces_numba(sorted_elements): for j in range(element_size): # exclude each node from sorted elements and get all combinations to obtain surfaces surfaces[i : i + elements_count, :j] = sorted_elements[:, :j] - surfaces[i : i + elements_count, j:dim] = sorted_elements[ - :, j + 1 : element_size - ] + surfaces[i : i + elements_count, j:dim] = sorted_elements[:, j + 1 : element_size] opposing_indices[i : i + elements_count] = sorted_elements[:, j] i += elements_count return surfaces, opposing_indices @@ -45,17 +43,13 @@ def extract_unique_indices(surfaces): def extract_unique_elements(elements, opposing_indices): - _, indices, count = np.unique( - elements, axis=0, return_index=True, return_counts=True - ) + _, indices, count = np.unique(elements, axis=0, return_index=True, return_counts=True) unique_indices = indices[count == 1] return elements[unique_indices], opposing_indices[unique_indices] def apply_predicate_to_surfaces(surfaces, nodes, predicate: Callable): - mask = [ - predicate(m) for m in np.mean(nodes[surfaces], axis=1) - ] # TODO: #65 Use numba (?) + mask = [predicate(m) for m in np.mean(nodes[surfaces], axis=1)] # TODO: #65 Use numba (?) return surfaces[mask] @@ -71,9 +65,7 @@ def reorder_boundary_nodes(nodes, elements, is_contact, is_dirichlet): Dirichlet override other """ # move boundary nodes to the top - nodes, elements, boundary_nodes_count = reorder( - nodes, elements, lambda _: True, to_top=True - ) + nodes, elements, boundary_nodes_count = reorder(nodes, elements, lambda _: True, to_top=True) # then move contact nodes to the top nodes, elements, contact_nodes_count = reorder( nodes, @@ -82,9 +74,7 @@ def reorder_boundary_nodes(nodes, elements, is_contact, is_dirichlet): to_top=True, ) # finally move dirichlet nodes to the bottom - nodes, elements, dirichlet_nodes_count = reorder( - nodes, elements, is_dirichlet, to_top=False - ) + nodes, elements, dirichlet_nodes_count = reorder(nodes, elements, is_dirichlet, to_top=False) return ( nodes, elements, @@ -173,17 +163,11 @@ def identify_boundaries_and_reorder_nodes( is_dirichlet=is_dirichlet, ) - neumann_nodes_count = ( - boundary_nodes_count - contact_nodes_count - dirichlet_nodes_count - ) + neumann_nodes_count = boundary_nodes_count - contact_nodes_count - dirichlet_nodes_count - (boundary_surfaces, boundary_internal_indices, *_) = get_boundary_surfaces( - elements - ) + (boundary_surfaces, boundary_internal_indices, *_) = get_boundary_surfaces(elements) - contact_boundary = apply_predicate_to_surfaces( - boundary_surfaces, initial_nodes, is_contact - ) + contact_boundary = apply_predicate_to_surfaces(boundary_surfaces, initial_nodes, is_contact) dirichlet_boundary = apply_predicate_to_surfaces( boundary_surfaces, initial_nodes, is_dirichlet ) @@ -200,14 +184,10 @@ def identify_boundaries_and_reorder_nodes( ) neumann_boundary = Boundary( surfaces=neumann_boundary, - node_indices=slice( - contact_nodes_count, contact_nodes_count + neumann_nodes_count - ), + node_indices=slice(contact_nodes_count, contact_nodes_count + neumann_nodes_count), node_count=neumann_nodes_count, ) - dirichlet_indices = slice( - len(initial_nodes) - dirichlet_nodes_count, len(initial_nodes) - ) + dirichlet_indices = slice(len(initial_nodes) - dirichlet_nodes_count, len(initial_nodes)) dirichlet_boundary = Boundary( surfaces=dirichlet_boundary, node_indices=dirichlet_indices, @@ -220,9 +200,7 @@ def identify_boundaries_and_reorder_nodes( other_boundaries = {} for name, indicator in boundaries_description.indicators.items(): if name not in ("contact", "dirichlet"): - surfaces = apply_predicate_to_surfaces( - boundary_surfaces, initial_nodes, indicator - ) + surfaces = apply_predicate_to_surfaces(boundary_surfaces, initial_nodes, indicator) node_indices = np.unique(surfaces) if len(node_indices): other_boundaries[name] = Boundary( diff --git a/conmech/mesh/mesh.py b/conmech/mesh/mesh.py index a359d3e4..289b6c29 100644 --- a/conmech/mesh/mesh.py +++ b/conmech/mesh/mesh.py @@ -105,9 +105,7 @@ def __init__( ) = BoundariesFactory.identify_boundaries_and_reorder_nodes( unordered_nodes, unordered_elements, boundaries_description ) - edges_matrix = get_edges_matrix( - nodes_count=len(self.nodes), elements=self.elements - ) + edges_matrix = get_edges_matrix(nodes_count=len(self.nodes), elements=self.elements) self.edges = get_edges_list_numba(edges_matrix) @property @@ -173,11 +171,7 @@ def independent_nodes_count(self): @property def free_nodes_count(self): # TODO: #65 CHECK - return ( - self.independent_nodes_count - - self.contact_nodes_count - - self.dirichlet_nodes_count - ) + return self.independent_nodes_count - self.contact_nodes_count - self.dirichlet_nodes_count @property def boundary_indices(self): diff --git a/conmech/mesh/zoo/barboteu_2008.py b/conmech/mesh/zoo/barboteu_2008.py index ac5832f6..340663c2 100644 --- a/conmech/mesh/zoo/barboteu_2008.py +++ b/conmech/mesh/zoo/barboteu_2008.py @@ -44,9 +44,7 @@ def __init__(self, mesh_descr: Barboteu2008MeshDescription): @staticmethod def _set_mesh_size(geom, mesh_descr: Barboteu2008MeshDescription): - geom.set_mesh_size_callback( - lambda dim, tag, x, y, z, *_: mesh_descr.max_element_perimeter - ) + geom.set_mesh_size_callback(lambda dim, tag, x, y, z, *_: mesh_descr.max_element_perimeter) @staticmethod def _get_nodes_and_elements(geom, dim): diff --git a/conmech/mesh/zoo/cross_for_tests.py b/conmech/mesh/zoo/cross_for_tests.py index 91973749..d44f2600 100644 --- a/conmech/mesh/zoo/cross_for_tests.py +++ b/conmech/mesh/zoo/cross_for_tests.py @@ -14,8 +14,7 @@ def __init__(self, mesh_descr: CrossMeshDescription): def _get_cross_rectangle(mesh_descr: CrossMeshDescription): scale_x, scale_y = mesh_descr.scale size_x, size_y = [ - int(np.ceil(scale / mesh_descr.max_element_perimeter)) - for scale in mesh_descr.scale + int(np.ceil(scale / mesh_descr.max_element_perimeter)) for scale in mesh_descr.scale ] min_ = np.array((0.0, 0.0)) @@ -44,21 +43,15 @@ def _set_cross_elements_numba( index = 0 for i in range(size_x): for j in range(size_y): - left_bottom = ( - np.array((i * edge_len_x, j * edge_len_y)) + left_bottom_node - ) + left_bottom = np.array((i * edge_len_x, j * edge_len_y)) + left_bottom_node lb = nph.get_node_index_numba(left_bottom, nodes) - rb = nph.get_node_index_numba( - left_bottom + np.array((edge_len_x, 0.0)), nodes - ) + rb = nph.get_node_index_numba(left_bottom + np.array((edge_len_x, 0.0)), nodes) c = nph.get_node_index_numba( left_bottom + np.array((0.5 * edge_len_x, 0.5 * edge_len_y)), nodes, ) - lt = nph.get_node_index_numba( - left_bottom + np.array((0.0, edge_len_y)), nodes - ) + lt = nph.get_node_index_numba(left_bottom + np.array((0.0, edge_len_y)), nodes) rt = nph.get_node_index_numba( left_bottom + np.array((edge_len_x, edge_len_y)), nodes ) @@ -80,9 +73,7 @@ def _set_cross_nodes_ordered_numba( index = 0 for j in range(size_y - 1, -1, -1): for i in range(size_x - 1, -1, -1): - nodes[index] = np.array( - ((i + 0.5) * edge_len_x, (j + 0.5) * edge_len_y) - ) + nodes[index] = np.array(((i + 0.5) * edge_len_x, (j + 0.5) * edge_len_y)) index += 1 for j in range(size_y - 1, 0, -1): diff --git a/conmech/mesh/zoo/rectangle.py b/conmech/mesh/zoo/rectangle.py index 41be4c02..c8e51f9c 100644 --- a/conmech/mesh/zoo/rectangle.py +++ b/conmech/mesh/zoo/rectangle.py @@ -27,14 +27,14 @@ class Rectangle(RawMesh): def __init__(self, mesh_descr: RectangleMeshDescription): scale_x, scale_y = mesh_descr.scale mesh_density = [ - int(np.ceil(scale / mesh_descr.max_element_perimeter)) - for scale in mesh_descr.scale + int(np.ceil(scale / mesh_descr.max_element_perimeter)) for scale in mesh_descr.scale ] # pylint: disable=no-member nodes, elements = meshzoo.rectangle_tri( - np.linspace(0.0, scale_x, int(mesh_density[0]) + 1), - np.linspace(0.0, scale_y, int(mesh_density[1]) + 1), + (0.0, scale_x), + (0.0, scale_y), + n=(int(mesh_density[0]) + 1, int(mesh_density[1]) + 1), variant="zigzag", ) super().__init__(nodes, elements) diff --git a/conmech/mesh/zoo/sofochbar_2023.py b/conmech/mesh/zoo/sofochbar_2023.py index 7815f462..b0f87368 100644 --- a/conmech/mesh/zoo/sofochbar_2023.py +++ b/conmech/mesh/zoo/sofochbar_2023.py @@ -61,8 +61,7 @@ def __init__(self, mesh_descr: SOB2023MeshDescription): foundation = dmsh.Path([[0.0, 0.0], [diameter + 2 * thickness, 0.0]]) dense_x, dense_y = [ - int(np.ceil(scale / mesh_descr.max_element_perimeter)) - for scale in mesh_descr.scale + int(np.ceil(scale / mesh_descr.max_element_perimeter)) for scale in mesh_descr.scale ] def target_edge_length(x): diff --git a/conmech/plotting/drawer.py b/conmech/plotting/drawer.py index 599a7e17..3d7fe7e5 100644 --- a/conmech/plotting/drawer.py +++ b/conmech/plotting/drawer.py @@ -268,9 +268,7 @@ def save_plot(self, format_, name=None): ) plt.close() - def draw_boundary( - self, edges, nodes, axes, label="", node_color="k", edge_color="k" - ): + def draw_boundary(self, edges, nodes, axes, label="", node_color="k", edge_color="k"): graph = nx.Graph() for edge in edges: graph.add_edge(edge[0], edge[1]) @@ -312,8 +310,6 @@ def draw_field(self, field, v_min, v_max, axes, fig): # from mpl_toolkits.axes_grid1 import make_axes_locatable # divider = make_axes_locatable(axes) # cax = divider.append_axes("bottom", size="5%", pad=0.15) - sm = plt.cm.ScalarMappable( - cmap=self.cmap, norm=plt.Normalize(vmin=v_min, vmax=v_max) - ) + sm = plt.cm.ScalarMappable(cmap=self.cmap, norm=plt.Normalize(vmin=v_min, vmax=v_max)) sm.set_array([]) fig.colorbar(sm, orientation="horizontal", label=self.field_label, ax=axes) diff --git a/conmech/plotting/membrane.py b/conmech/plotting/membrane.py index ae2bb05e..54982c4a 100644 --- a/conmech/plotting/membrane.py +++ b/conmech/plotting/membrane.py @@ -56,9 +56,7 @@ def plot_limit_points( label = None buff_size = 0 if buff_size > 0: - plt.scatter( - buff[0, :buff_size], buff[1, :buff_size], s=1, color=color, label=label - ) + plt.scatter(buff[0, :buff_size], buff[1, :buff_size], s=1, color=color, label=label) plt.title(title) plt.ylim(*ylim) @@ -93,9 +91,7 @@ def do_plot( mesh_node_x = state.body.mesh.nodes[:, 0] mesh_node_y = state.body.mesh.nodes[:, 1] - soltri = tri.Triangulation( - mesh_node_x, mesh_node_y, triangles=state.body.mesh.elements - ) + soltri = tri.Triangulation(mesh_node_x, mesh_node_y, triangles=state.body.mesh.elements) interpol = tri.LinearTriInterpolator # if in3d else tri.CubicTriInterpolator v = interpol(soltri, getattr(state, field)[:, 0]) u = interpol(soltri, state.displacement[:, 0]) diff --git a/conmech/scenarios/problems.py b/conmech/scenarios/problems.py index 2dd8e78f..af5d84b9 100644 --- a/conmech/scenarios/problems.py +++ b/conmech/scenarios/problems.py @@ -152,9 +152,7 @@ def initial_velocity(x: np.ndarray) -> np.ndarray: return np.zeros_like(x) -class QuasistaticDisplacementProblem( - QuasistaticProblem, TimeDependentDisplacementProblem, ABC -): +class QuasistaticDisplacementProblem(QuasistaticProblem, TimeDependentDisplacementProblem, ABC): @classmethod def statement(cls, body) -> Statement: return QuasistaticVelocityStatement(body) diff --git a/conmech/scene/body_forces.py b/conmech/scene/body_forces.py index 4b8ac461..40dd7b98 100644 --- a/conmech/scene/body_forces.py +++ b/conmech/scene/body_forces.py @@ -10,16 +10,12 @@ @numba.njit -def get_surface_per_boundary_node_numba( - boundary_surfaces, considered_nodes_count, moved_nodes -): +def get_surface_per_boundary_node_numba(boundary_surfaces, considered_nodes_count, moved_nodes): surface_per_boundary_node = np.zeros((considered_nodes_count, 1), dtype=np.float64) for boundary_surface in boundary_surfaces: face_nodes = moved_nodes[boundary_surface] - surface_per_boundary_node[boundary_surface] += element_volume_part_numba( - face_nodes - ) + surface_per_boundary_node[boundary_surface] += element_volume_part_numba(face_nodes) return surface_per_boundary_node @@ -54,9 +50,7 @@ class FieldSource: def node_source(self, nodes, time: float): # pylint: disable=not-callable if time != self.timestamp: - self.cache = np.array( - [self.source(nodes[i], time) for i in range(len(nodes))] - ) + self.cache = np.array([self.source(nodes[i], time) for i in range(len(nodes))]) self.timestamp = time return self.cache diff --git a/conmech/simulations/problem_solver.py b/conmech/simulations/problem_solver.py index 78ac06d5..13d9037d 100644 --- a/conmech/simulations/problem_solver.py +++ b/conmech/simulations/problem_solver.py @@ -182,15 +182,11 @@ def run( print(f"{self.done / self.to_do * 100:.2f}%", end="\r") def find_solution(self, state, **kwargs) -> np.ndarray: - initial_guess = ( - state[self.coordinates].T.ravel().reshape(state.body.mesh.dimension, -1) - ) + initial_guess = state[self.coordinates].T.ravel().reshape(state.body.mesh.dimension, -1) solution = self.step_solver.solve(initial_guess, **kwargs) return solution - def find_solution_uzawa( - self, solution, solution_t - ) -> Tuple[np.ndarray, np.ndarray]: + def find_solution_uzawa(self, solution, solution_t) -> Tuple[np.ndarray, np.ndarray]: # TODO #95 norm = np.inf old_solution = solution.copy().reshape(-1, 1).squeeze() @@ -254,8 +250,7 @@ def find_solution_uzawa( ) self.second_step_solver.v_vector = velocity.reshape(-1) self.second_step_solver.u_vector = ( - old_u_vector - + self.second_step_solver.time_step * self.second_step_solver.v_vector + old_u_vector + self.second_step_solver.time_step * self.second_step_solver.v_vector ) self.step_solver.p_vector = solution_t self.second_step_solver.p_vector = solution_t @@ -318,9 +313,7 @@ def __init__(self, problem: StaticDisplacementProblem, solving_method: str): # super class method takes **kwargs, so signatures are consistent # pylint: disable=arguments-differ - def solve( - self, *, initial_displacement: Callable, verbose: bool = False, **kwargs - ) -> State: + def solve(self, *, initial_displacement: Callable, verbose: bool = False, **kwargs) -> State: """ :param initial_displacement: for the solver :param verbose: show prints @@ -390,9 +383,7 @@ def solve( output_step = (0, *output_step) if output_step else (0, n_steps) # 0 for diff state = State(self.body) - state.absement[:] = initial_absement( - self.body.mesh.nodes[: self.body.mesh.nodes_count] - ) + state.absement[:] = initial_absement(self.body.mesh.nodes[: self.body.mesh.nodes_count]) state.displacement[:] = initial_displacement( self.body.mesh.nodes[: self.body.mesh.nodes_count] ) @@ -462,9 +453,7 @@ def solve( state.displacement[:] = initial_displacement( self.body.mesh.nodes[: self.body.mesh.nodes_count] ) - state.velocity[:] = initial_velocity( - self.body.mesh.nodes[: self.body.mesh.nodes_count] - ) + state.velocity[:] = initial_velocity(self.body.mesh.nodes[: self.body.mesh.nodes_count]) self.step_solver.u_vector[:] = state.displacement.T.ravel().copy() self.step_solver.v_vector[:] = state.velocity.T.ravel().copy() @@ -534,9 +523,7 @@ def solve( state.displacement[:] = initial_displacement( self.body.mesh.nodes[: self.body.mesh.nodes_count] ) - state.velocity[:] = initial_velocity( - self.body.mesh.nodes[: self.body.mesh.nodes_count] - ) + state.velocity[:] = initial_velocity(self.body.mesh.nodes[: self.body.mesh.nodes_count]) state.temperature[:] = initial_temperature( self.body.mesh.nodes[: self.body.mesh.nodes_count] ) @@ -562,9 +549,7 @@ def solve( done += 1 print(f"{done/n_steps*100:.2f}%", end="\r") self.step_solver.current_time += self.step_solver.time_step - self.second_step_solver.current_time += ( - self.second_step_solver.time_step - ) + self.second_step_solver.current_time += self.second_step_solver.time_step solution, solution_t = self.find_solution_uzawa(solution, solution_t) @@ -631,9 +616,7 @@ def solve( state.displacement[:] = initial_displacement( self.body.mesh.nodes[: self.body.mesh.nodes_count] ) - state.velocity[:] = initial_velocity( - self.body.mesh.nodes[: self.body.mesh.nodes_count] - ) + state.velocity[:] = initial_velocity(self.body.mesh.nodes[: self.body.mesh.nodes_count]) state.electric_potential[:] = initial_electric_potential( self.body.mesh.nodes[: self.body.mesh.nodes_count] ) @@ -660,9 +643,7 @@ def solve( done += 1 print(f"{done/n_steps*100:.2f}%", end="\r") self.step_solver.current_time += self.step_solver.time_step - self.second_step_solver.current_time += ( - self.second_step_solver.time_step - ) + self.second_step_solver.current_time += self.second_step_solver.time_step solution, solution_t = self.find_solution_uzawa(solution, solution_t) @@ -731,9 +712,7 @@ def solve( state.displacement[:] = initial_displacement( self.body.mesh.nodes[: self.body.mesh.nodes_count] ) - state.velocity[:] = initial_velocity( - self.body.mesh.nodes[: self.body.mesh.nodes_count] - ) + state.velocity[:] = initial_velocity(self.body.mesh.nodes[: self.body.mesh.nodes_count]) self.step_solver.current_time = state.time diff --git a/conmech/solvers/optimization/optimization.py b/conmech/solvers/optimization/optimization.py index 5b34a5ac..19f84d24 100644 --- a/conmech/solvers/optimization/optimization.py +++ b/conmech/solvers/optimization/optimization.py @@ -54,9 +54,7 @@ def __init__( if isinstance(statement, StaticPoissonStatement): self.loss = make_cost_functional( normal_condition=( - contact_law.potential_normal_direction - if contact_law is not None - else None + contact_law.potential_normal_direction if contact_law is not None else None ), variable_dimension=statement.dimension_out, problem_dimension=statement.dimension_in, @@ -118,9 +116,7 @@ def _solve_impl( # pylint: disable=import-outside-toplevel,import-error) from kosopt import qsmlm - solution = qsmlm.minimize( - self.loss, solution, args=args, maxiter=maxiter - ) + solution = qsmlm.minimize(self.loss, solution, args=args, maxiter=maxiter) sols.append(solution.copy()) elif method.lower() == "constrained": contact_nodes_count = self.body.mesh.boundaries.contact_nodes_count diff --git a/conmech/solvers/optimization/schur_complement.py b/conmech/solvers/optimization/schur_complement.py index 7bad5c56..58658815 100644 --- a/conmech/solvers/optimization/schur_complement.py +++ b/conmech/solvers/optimization/schur_complement.py @@ -128,9 +128,7 @@ def complement_free_nodes(self, truncated_solution: np.ndarray) -> np.ndarray: result = self.free_x_free_inverted @ _result return result - def merge( - self, solution_contact: np.ndarray, solution_free: np.ndarray - ) -> np.ndarray: + def merge(self, solution_contact: np.ndarray, solution_free: np.ndarray) -> np.ndarray: if self.statement.dimension_in == 2 or self.driving_vector: u_contact = solution_contact.reshape(2, -1) u_free = solution_free.reshape(2, -1) @@ -144,16 +142,12 @@ def merge( return result -@SolversRegistry.register( - "static", "schur", "schur complement", "schur complement method" -) +@SolversRegistry.register("static", "schur", "schur complement", "schur complement method") class StaticSchurOptimization(SchurComplementOptimization): pass -@SolversRegistry.register( - "quasistatic", "schur", "schur complement", "schur complement method" -) +@SolversRegistry.register("quasistatic", "schur", "schur complement", "schur complement method") class QuasistaticSchurOptimization(SchurComplementOptimization): def iterate(self): self.statement.update( @@ -183,9 +177,7 @@ def iterate(self): self.node_forces_, self.forces_free = self.recalculate_forces() -@SolversRegistry.register( - "dynamic", "schur", "schur complement", "schur complement method" -) +@SolversRegistry.register("dynamic", "schur", "schur complement", "schur complement method") class DynamicSchurOptimization(SchurComplementOptimization): def iterate(self): self.statement.update( @@ -212,9 +204,7 @@ def calculate_schur_complement_vector( vector_split = nph.unstack(vector, dimension) vector_contact = nph.stack_column(vector_split[contact_indices, :]) vector_free = nph.stack_column(vector_split[free_indices, :]) - vector_boundary = vector_contact - ( - contact_x_free @ (free_x_free_inverted @ vector_free) - ) + vector_boundary = vector_contact - (contact_x_free @ (free_x_free_inverted @ vector_free)) return vector_boundary, vector_free @@ -235,8 +225,6 @@ def get_sliced(matrix_split, indices_height, indices_width): contact_x_contact = get_sliced(matrix_split, contact_indices, contact_indices) free_x_free_inverted = np.linalg.inv(free_x_free) - matrix_boundary = contact_x_contact - contact_x_free @ ( - free_x_free_inverted @ free_x_contact - ) + matrix_boundary = contact_x_contact - contact_x_free @ (free_x_free_inverted @ free_x_contact) return matrix_boundary, free_x_contact, contact_x_free, free_x_free_inverted diff --git a/conmech/solvers/solver.py b/conmech/solvers/solver.py index 4d62648b..a0ba4aee 100644 --- a/conmech/solvers/solver.py +++ b/conmech/solvers/solver.py @@ -27,12 +27,8 @@ def __init__( self.time_step: float = time_step self.current_time: float = 0 self.b_vector = np.zeros(self.body.mesh.nodes_count * self.body.mesh.dimension) - self.u_vector: np.ndarray = np.zeros( - self.body.mesh.nodes_count * self.body.mesh.dimension - ) - self.v_vector: np.ndarray = np.zeros( - self.body.mesh.nodes_count * self.body.mesh.dimension - ) + self.u_vector: np.ndarray = np.zeros(self.body.mesh.nodes_count * self.body.mesh.dimension) + self.v_vector: np.ndarray = np.zeros(self.body.mesh.nodes_count * self.body.mesh.dimension) self.t_vector: np.ndarray = np.zeros(self.body.mesh.nodes_count) self.p_vector: np.ndarray = np.zeros(self.body.mesh.nodes_count) # TODO #23 diff --git a/conmech/solvers/solver_methods.py b/conmech/solvers/solver_methods.py index 707bef65..ee763124 100644 --- a/conmech/solvers/solver_methods.py +++ b/conmech/solvers/solver_methods.py @@ -55,8 +55,7 @@ def equation( res = ( 0.5 * np.dot(np.dot(lhs, var[:ind]), var[:ind]) - np.dot(rhs, var[:ind]) - + 0.5 - * np.dot(np.dot(volume_multiplier, response), np.ones_like(var[:ind])) + + 0.5 * np.dot(np.dot(volume_multiplier, response), np.ones_like(var[:ind])) + np.dot(var[ind:], var[ind:].T) ) @@ -77,9 +76,7 @@ def contact_part(u_vector, nodes, contact_boundary, contact_normals): normal_vector = contact_normals[ei] # ASSUMING `u_vector` and `nodes` have the same order! - um = interpolate_node_between( - edge, u_vector, u_vector, dimension=problem_dimension - ) + um = interpolate_node_between(edge, u_vector, u_vector, dimension=problem_dimension) um_normal = (um * normal_vector).sum() edge_len = nph.length(edge, nodes) @@ -158,9 +155,7 @@ def contact_cost_functional( edge = contact_boundary[ei] normal_vector = contact_normals[ei] # ASSUMING `u_vector` and `nodes` have the same order! - vm = interpolate_node_between( - edge, var, var_old, dimension=variable_dimension - ) + vm = interpolate_node_between(edge, var, var_old, dimension=variable_dimension) if variable_dimension == 1: vm_normal = vm[0] vm_tangential = np.empty(0) @@ -174,9 +169,7 @@ def contact_cost_functional( static_displacement, dimension=problem_dimension, ) - static_displacement_normal = ( - static_displacement_mean * normal_vector - ).sum() + static_displacement_normal = (static_displacement_mean * normal_vector).sum() static_displacement_tangential = ( static_displacement_mean - static_displacement_normal * normal_vector ) @@ -212,11 +205,7 @@ def cost_functional( var, var_old, u_vector, nodes, contact_boundary, contact_normals, dt ) ind = lhs.shape[0] - result = ( - 0.5 * np.dot(np.dot(lhs, var[:ind]), var[:ind]) - - np.dot(rhs, var[:ind]) - + ju - ) + result = 0.5 * np.dot(np.dot(lhs, var[:ind]), var[:ind]) - np.dot(rhs, var[:ind]) + ju result = np.asarray(result).ravel() return result diff --git a/conmech/state/state.py b/conmech/state/state.py index a2e6410f..a4646c8c 100644 --- a/conmech/state/state.py +++ b/conmech/state/state.py @@ -14,16 +14,12 @@ def __init__(self, body): self.products = {} - self.absement: np.ndarray = np.zeros( - (self.body.mesh.nodes_count, self.body.mesh.dimension) - ) + self.absement: np.ndarray = np.zeros((self.body.mesh.nodes_count, self.body.mesh.dimension)) self.displacement: np.ndarray = np.zeros( (self.body.mesh.nodes_count, self.body.mesh.dimension) ) self.displaced_nodes: np.ndarray = np.copy(self.body.mesh.nodes) - self.velocity: np.ndarray = np.zeros( - (self.body.mesh.nodes_count, self.body.mesh.dimension) - ) + self.velocity: np.ndarray = np.zeros((self.body.mesh.nodes_count, self.body.mesh.dimension)) self.setup = None self.__stress: np.ndarray = None self.constitutive_law = None @@ -37,12 +33,9 @@ def set_displacement( *, update_absement: bool = False, ): - self.displacement = displacement_vector.reshape( - (self.body.mesh.dimension, -1) - ).T + self.displacement = displacement_vector.reshape((self.body.mesh.dimension, -1)).T self.displaced_nodes[: self.body.mesh.nodes_count, :] = ( - self.body.mesh.nodes[: self.body.mesh.nodes_count, :] - + self.displacement[:, :] + self.body.mesh.nodes[: self.body.mesh.nodes_count, :] + self.displacement[:, :] ) if update_absement: dt = time - self.time @@ -50,16 +43,13 @@ def set_displacement( self.time = time self.update_products() - def set_velocity( - self, velocity_vector: np.ndarray, time: float, *, update_displacement: bool - ): + def set_velocity(self, velocity_vector: np.ndarray, time: float, *, update_displacement: bool): self.velocity = velocity_vector.reshape((self.body.mesh.dimension, -1)).T if update_displacement: dt = time - self.time self.displacement += dt * self.velocity self.displaced_nodes[: self.body.mesh.nodes_count, :] = ( - self.body.mesh.nodes[: self.body.mesh.nodes_count, :] - + self.displacement[:, :] + self.body.mesh.nodes[: self.body.mesh.nodes_count, :] + self.displacement[:, :] ) self.time = time self.update_products() diff --git a/examples/BOSK_2024_example_1.py b/examples/BOSK_2024_example_1.py index 83b2c116..c0032f8a 100644 --- a/examples/BOSK_2024_example_1.py +++ b/examples/BOSK_2024_example_1.py @@ -72,9 +72,7 @@ def main(config: Config, setup, name, steps): n_steps=steps, output_step=output_step, products=[ - VerticalIntersectionContactLimitPoints( - obstacle_level=OBSTACLE_LEVEL, x=1.0 - ), + VerticalIntersectionContactLimitPoints(obstacle_level=OBSTACLE_LEVEL, x=1.0), VerticalIntersection(x=1.0), ], initial_displacement=setup.initial_displacement, diff --git a/examples/BOSK_2024_example_2.py b/examples/BOSK_2024_example_2.py index bd3efdf2..400bca66 100644 --- a/examples/BOSK_2024_example_2.py +++ b/examples/BOSK_2024_example_2.py @@ -58,9 +58,7 @@ def outer_forces( ) -> np.ndarray: return np.array([0.0]) - boundaries: ... = BoundariesDescription( - dirichlet=lambda x: x[0] in (0, 1) or x[1] in (0, 1) - ) + boundaries: ... = BoundariesDescription(dirichlet=lambda x: x[0] in (0, 1) or x[1] in (0, 1)) def main(config: Config, setup, name, steps): @@ -96,9 +94,7 @@ def main(config: Config, setup, name, steps): n_steps=steps, output_step=output_step, products=[ - VerticalIntersectionContactLimitPoints( - obstacle_level=OBSTACLE_LEVEL, x=0.50 - ), + VerticalIntersectionContactLimitPoints(obstacle_level=OBSTACLE_LEVEL, x=0.50), VerticalIntersection(x=0.50), ], initial_displacement=setup.initial_displacement, @@ -122,14 +118,12 @@ def main(config: Config, setup, name, steps): state = State.load(path) plot_limit_points( state.products["limit points at 0.50"], - title=rf"$\kappa={setup.contact_law.KAPPA}$ " - rf"$\beta={setup.contact_law.BETA}$", + title=rf"$\kappa={setup.contact_law.KAPPA}$ " rf"$\beta={setup.contact_law.BETA}$", finish=False, ) plot_limit_points( states[-1].products["limit points at 0.50"], - title=rf"$\kappa={setup.contact_law.KAPPA}$ " - rf"$\beta={setup.contact_law.BETA}$", + title=rf"$\kappa={setup.contact_law.KAPPA}$ " rf"$\beta={setup.contact_law.BETA}$", finish=False, ) print(time.time() - start) @@ -152,9 +146,7 @@ def main(config: Config, setup, name, steps): return states_ids = list(range(len(states))) - to_plot = states_ids[ - :64:4 - ] # [2, 4, 8, 11, 13, 15, 17, 20] #states_ids[1:4] + states_ids[-1:] + to_plot = states_ids[:64:4] # [2, 4, 8, 11, 13, 15, 17, 20] #states_ids[1:4] + states_ids[-1:] vmin = np.inf vmax = -np.inf field = "velocity" diff --git a/examples/BOSK_plot_postprocessing.py b/examples/BOSK_plot_postprocessing.py index 60ff49ec..4d6c30a7 100644 --- a/examples/BOSK_plot_postprocessing.py +++ b/examples/BOSK_plot_postprocessing.py @@ -38,9 +38,7 @@ def case1(): (100.0, 0.5), )[:]: print((kappa, beta)) - states[kappa, beta] = State.load( - f"output/BOSK.pub/c1_kappa={kappa:.2f};beta={beta:.2f}" - ) + states[kappa, beta] = State.load(f"output/BOSK.pub/c1_kappa={kappa:.2f};beta={beta:.2f}") c1_reference(states, output_path="output/BOSK.pub") c1_steady_state(states, output_path="output/BOSK.pub") @@ -48,12 +46,8 @@ def case1(): def show(output_path, name): - plt.gca().yaxis.set_major_formatter( - StrMethodFormatter("{x:,.2f}") - ) # 2 decimal places - plt.gca().xaxis.set_major_formatter( - StrMethodFormatter("{x:,.2f}") - ) # 2 decimal places + plt.gca().yaxis.set_major_formatter(StrMethodFormatter("{x:,.2f}")) # 2 decimal places + plt.gca().xaxis.set_major_formatter(StrMethodFormatter("{x:,.2f}")) # 2 decimal places if output_path is None: plt.show() else: diff --git a/examples/BOST_2024.py b/examples/BOST_2024.py index 65301edc..f1598769 100644 --- a/examples/BOST_2024.py +++ b/examples/BOST_2024.py @@ -45,9 +45,7 @@ def potential_tangential_direction( return np.sum(var_tau * var_tau) ** 0.5 @staticmethod - def tangential_bound( - var_nu: float, static_displacement_nu: float, dt: float - ) -> float: + def tangential_bound(var_nu: float, static_displacement_nu: float, dt: float) -> float: return 1.0 @@ -93,9 +91,7 @@ def potential_normal_direction( elif var_nu < b: result = (a + np.exp(-b)) / (2 * b) * var_nu**2 else: - result = ( - a * var_nu - np.exp(-var_nu) + ((b + 2) * np.exp(-b) - a * b) / 2 - ) + result = a * var_nu - np.exp(-var_nu) + ((b + 2) * np.exp(-b) - a * b) / 2 return ig * result kwargs = {"method": "POWELL"} @@ -164,15 +160,12 @@ def main( if ig == ig_prev: initial_displacement = setup.initial_displacement else: - with open( - f"{config.outputs_path}/{PREFIX}_ig_{ig_prev}", "rb" - ) as output: + with open(f"{config.outputs_path}/{PREFIX}_ig_{ig_prev}", "rb") as output: state = pickle.load(output) initial_displacement = lambda _: state.displacement.copy() state.displaced_nodes[: state.body.mesh.nodes_count, :] = ( - state.body.mesh.nodes[: state.body.mesh.nodes_count, :] - + state.displacement[:, :] + state.body.mesh.nodes[: state.body.mesh.nodes_count, :] + state.displacement[:, :] ) with open(f"{config.outputs_path}/{PREFIX}_ig_{ig}", "wb+") as output: state.body.dynamics.force.outer.source = None @@ -216,9 +209,7 @@ def main( # Y = -np.asarray(Y) print(f"{X=}") print(f"{Y=}") - plot_errors( - X, Y, highlighted_id=None, save=f"{config.outputs_path}/convergence.pdf" - ) + plot_errors(X, Y, highlighted_id=None, save=f"{config.outputs_path}/convergence.pdf") def plot_errors(X, Y, highlighted_id, save: Optional[str] = None): diff --git a/examples/Jureczka_Ochal_Bartman_2023.py b/examples/Jureczka_Ochal_Bartman_2023.py index c795fbaa..a0af98d9 100644 --- a/examples/Jureczka_Ochal_Bartman_2023.py +++ b/examples/Jureczka_Ochal_Bartman_2023.py @@ -24,9 +24,7 @@ def make_contact_law(limit_value, limit, friction_bound): class JureczkaOchalBartman2023(PotentialOfContactLaw): @staticmethod - def tangential_bound( - var_nu: float, static_displacement_nu: float, dt: float - ) -> float: + def tangential_bound(var_nu: float, static_displacement_nu: float, dt: float) -> float: n = friction_bound b = 0.1 if var_nu <= 0: @@ -138,9 +136,7 @@ def main(config: Config): else: simulate = True - mesh_descr = JOB2023MeshDescription( - initial_position=None, max_element_perimeter=1 / h - ) + mesh_descr = JOB2023MeshDescription(initial_position=None, max_element_perimeter=1 / h) if simulate: for name in names: @@ -198,9 +194,7 @@ def main(config: Config): drawer.node_size = 1 drawer.original_mesh_color = "k" drawer.deformed_mesh_color = None - drawer.draw( - show=config.show, field_min=0, field_max=40, save=config.save - ) + drawer.draw(show=config.show, field_min=0, field_max=40, save=config.save) state.setup = setup state.constitutive_law = viscoelastic_constitutive_law drawer = Drawer(state=state, config=config) diff --git a/examples/Jureczka_and_Ochal_2019.py b/examples/Jureczka_and_Ochal_2019.py index 2f53c10f..21fa6a34 100644 --- a/examples/Jureczka_and_Ochal_2019.py +++ b/examples/Jureczka_and_Ochal_2019.py @@ -46,9 +46,7 @@ def potential_tangential_direction( return np.log(np.sum(var_tau * var_tau) ** 0.5 + 1) @staticmethod - def tangential_bound( - var_nu: float, static_displacement_nu: float, dt: float - ) -> float: + def tangential_bound(var_nu: float, static_displacement_nu: float, dt: float) -> float: if static_displacement_nu < 0: return 0 if static_displacement_nu < 0.1: diff --git a/examples/Sofonea_Ochal_Bartman_2023.py b/examples/Sofonea_Ochal_Bartman_2023.py index 87f3fae7..306c139d 100644 --- a/examples/Sofonea_Ochal_Bartman_2023.py +++ b/examples/Sofonea_Ochal_Bartman_2023.py @@ -150,15 +150,11 @@ def zero_relaxation(t=None): with open(f"{config.outputs_path}/{name}_h_{h}_global", "rb") as output: _ = pickle.load(output) - with open( - f"{config.outputs_path}/{name}_h_{h}_penetration", "rb" - ) as output: + with open(f"{config.outputs_path}/{name}_h_{h}_penetration", "rb") as output: _ = pickle.load(output) for time_step in steps: - with open( - f"{config.outputs_path}/{name}_t_{time_step}_h_{h}", "rb" - ) as output: + with open(f"{config.outputs_path}/{name}_t_{time_step}_h_{h}", "rb") as output: _ = pickle.load(output) except IOError: simulate = True @@ -223,9 +219,7 @@ def zero_relaxation(t=None): f_limits = pickle.load(output) for time_step in steps: - with open( - f"{config.outputs_path}/{name}_t_{time_step}_h_{h}", "rb" - ) as output: + with open(f"{config.outputs_path}/{name}_t_{time_step}_h_{h}", "rb") as output: state = pickle.load(output) # Workaround state.body.dynamics.force.outer.source = examples[name]["outer_forces"] @@ -266,9 +260,7 @@ def zero_relaxation(t=None): ) drawer.ylabel = "y" axes[1].axis("on") - axes[1].tick_params( - left=True, bottom=True, labelleft=True, labelbottom=True - ) + axes[1].tick_params(left=True, bottom=True, labelleft=True, labelbottom=True) axes[1].set_aspect("equal", adjustable="box") zoom_outside(axes[0], [3, -1.5, 6, 2], axes[1], color="gray") @@ -285,27 +277,15 @@ def zero_relaxation(t=None): save=False, ) if time_step == 0: - axes[0].annotate( - "$\Gamma_1$", xy=(0, 0), xytext=(0.33, -0.50), fontsize=18 - ) + axes[0].annotate("$\Gamma_1$", xy=(0, 0), xytext=(0.33, -0.50), fontsize=18) position = (3.66, 4.5) - axes[0].annotate( - "$\Gamma_2$", xy=(0, 0), xytext=position, fontsize=18 - ) + axes[0].annotate("$\Gamma_2$", xy=(0, 0), xytext=position, fontsize=18) axes[0].add_patch(Rectangle(position, 0.3, 0.3, color="white")) - axes[0].annotate( - "$\Gamma_2$", xy=(0, 0), xytext=(2.33, 3.0), fontsize=18 - ) - axes[0].annotate( - "$\mathbf{f}_2$", xy=(0, 0), xytext=(2.5, 5.00), fontsize=15 - ) - axes[0].annotate( - "$\Gamma_3$", xy=(0, 0), xytext=(4.33, -0.50), fontsize=18 - ) + axes[0].annotate("$\Gamma_2$", xy=(0, 0), xytext=(2.33, 3.0), fontsize=18) + axes[0].annotate("$\mathbf{f}_2$", xy=(0, 0), xytext=(2.5, 5.00), fontsize=15) + axes[0].annotate("$\Gamma_3$", xy=(0, 0), xytext=(4.33, -0.50), fontsize=18) axes[0].axis("on") - axes[0].tick_params( - left=True, bottom=True, labelleft=True, labelbottom=True - ) + axes[0].tick_params(left=True, bottom=True, labelleft=True, labelbottom=True) axes[0].set_aspect("equal", adjustable="box") fig.tight_layout(rect=[0, 0, 1, 1.2]) @@ -352,9 +332,7 @@ def plots(setup, h, examples, config): old_p = p plot_outer_force(axes[0, col], frc, t, vertical_line=pnt_sig_change) - plot_displacement_normal_direction( - axes[1, col], pnt, t, vertical_line=pnt_sig_change - ) + plot_displacement_normal_direction(axes[1, col], pnt, t, vertical_line=pnt_sig_change) format_ = "pdf" path = Drawer.get_output_path(config, format_, name="force_penetration") @@ -426,9 +404,7 @@ def zoom_outside( + list(arrow_kwargs.items()) ) # draw a rectangle on original chart - src_ax.add_patch( - Rectangle([roi[0], roi[1]], roi[2] - roi[0], roi[3] - roi[1], **roi_kwargs) - ) + src_ax.add_patch(Rectangle([roi[0], roi[1]], roi[2] - roi[0], roi[3] - roi[1], **roi_kwargs)) if not draw_lines: return diff --git a/examples/example_dynamic_membrane.py b/examples/example_dynamic_membrane.py index fad0ca05..ac97d1ba 100644 --- a/examples/example_dynamic_membrane.py +++ b/examples/example_dynamic_membrane.py @@ -43,9 +43,7 @@ def outer_forces( ) -> np.ndarray: return np.array([0.0]) - boundaries: ... = BoundariesDescription( - dirichlet=lambda x: x[0] in (0, 1) or x[1] in (0, 1) - ) + boundaries: ... = BoundariesDescription(dirichlet=lambda x: x[0] in (0, 1) or x[1] in (0, 1)) def main(config: Config): diff --git a/examples/example_piezo_quasistatic.py b/examples/example_piezo_quasistatic.py index 3d1de131..d4153548 100644 --- a/examples/example_piezo_quasistatic.py +++ b/examples/example_piezo_quasistatic.py @@ -19,9 +19,7 @@ class PPSlopeContactLaw(PotentialOfContactLaw): @staticmethod - def tangential_bound( - var_nu: float, static_displacement_nu: float, dt: float - ) -> float: + def tangential_bound(var_nu: float, static_displacement_nu: float, dt: float) -> float: return -1.0 @staticmethod @@ -61,9 +59,7 @@ class PQuasistaticSetup(PiezoelectricQuasistaticProblem): ) ) permittivity: ... = field( - default_factory=lambda: np.array( - [[8.3, 0.0, 0.0], [0.0, 8.8, 0.0], [0.0, 0.0, -8]] - ) + default_factory=lambda: np.array([[8.3, 0.0, 0.0], [0.0, 8.8, 0.0], [0.0, 0.0, -8]]) ) @staticmethod @@ -98,9 +94,7 @@ def main(config: Config): To see result of simulation you need to call from python `main(Config().init())`. """ - mesh_descr = Barboteu2008MeshDescription( - initial_position=None, max_element_perimeter=0.5 - ) + mesh_descr = Barboteu2008MeshDescription(initial_position=None, max_element_perimeter=0.5) setup = PQuasistaticSetup(mesh_descr) runner = PiezoelectricTimeDependentSolver(setup, solving_method="global") steps = 100 if not config.test else 2 diff --git a/examples/example_piezoelectric_dynamic.py b/examples/example_piezoelectric_dynamic.py index e3bdeb42..086c81e9 100644 --- a/examples/example_piezoelectric_dynamic.py +++ b/examples/example_piezoelectric_dynamic.py @@ -19,9 +19,7 @@ # TODO # 48 class PPSlopeContactLaw(PotentialOfContactLaw): @staticmethod - def tangential_bound( - var_nu: float, static_displacement_nu: float, dt: float - ) -> float: + def tangential_bound(var_nu: float, static_displacement_nu: float, dt: float) -> float: return -1.0 @staticmethod @@ -61,9 +59,7 @@ class PDynamicSetup(PiezoelectricDynamicProblem): ) ) permittivity: ... = field( - default_factory=lambda: np.array( - [[8.3, 0.0, 0.0], [0.0, 8.8, 0.0], [0.0, 0.0, -8]] - ) + default_factory=lambda: np.array([[8.3, 0.0, 0.0], [0.0, 8.8, 0.0], [0.0, 0.0, -8]]) ) @staticmethod diff --git a/examples/example_static.py b/examples/example_static.py index b0dbb0ef..7bb84571 100644 --- a/examples/example_static.py +++ b/examples/example_static.py @@ -71,9 +71,7 @@ def main(config: Config, dimension=2): from mpl_toolkits.mplot3d.art3d import Poly3DCollection axs.add_collection3d( - Poly3DCollection( - faces, facecolors="cyan", linewidths=1, edgecolors="r", alpha=0.25 - ) + Poly3DCollection(faces, facecolors="cyan", linewidths=1, edgecolors="r", alpha=0.25) ) plt.show() diff --git a/examples/example_temperature_dynamic.py b/examples/example_temperature_dynamic.py index dfcc4302..24c3574d 100644 --- a/examples/example_temperature_dynamic.py +++ b/examples/example_temperature_dynamic.py @@ -49,14 +49,10 @@ class TDynamicSetup(TemperatureDynamicProblem): contact_law: ... = make_slope_contact_law(slope=1e1) contact_law_2: ... = TPSlopeContactLaw thermal_expansion: ... = field( - default_factory=lambda: np.array( - [[0.5, 0.0, 0.0], [0.0, 0.5, 0.0], [0.0, 0.0, 0.5]] - ) + default_factory=lambda: np.array([[0.5, 0.0, 0.0], [0.0, 0.5, 0.0], [0.0, 0.0, 0.5]]) ) thermal_conductivity: ... = field( - default_factory=lambda: np.array( - [[0.1, 0.0, 0.0], [0.0, 0.1, 0.0], [0.0, 0.0, 0.1]] - ) + default_factory=lambda: np.array([[0.1, 0.0, 0.0], [0.0, 0.1, 0.0], [0.0, 0.0, 0.1]]) ) @staticmethod @@ -75,9 +71,7 @@ def outer_forces(x, t=None): return np.array([-48.0 * (0.25 - (x[1] - 0.5) ** 2), 0]) return np.array([0, 0]) - boundaries: ... = BoundariesDescription( - contact=lambda x: x[1] == 0, dirichlet=lambda x: False - ) + boundaries: ... = BoundariesDescription(contact=lambda x: x[1] == 0, dirichlet=lambda x: False) def main(config: Config): @@ -111,9 +105,7 @@ def main(config: Config): drawer = Drawer(state=state, config=config) drawer.cmap = "plasma" drawer.field_name = "temperature" - drawer.draw( - field_max=T_max, field_min=T_min, show=config.show, save=config.save - ) + drawer.draw(field_max=T_max, field_min=T_min, show=config.show, save=config.save) if __name__ == "__main__": diff --git a/examples/utils.py b/examples/utils.py index 296afe5a..04a09f6b 100644 --- a/examples/utils.py +++ b/examples/utils.py @@ -78,9 +78,7 @@ def gradient(elements, initial_nodes, f): return result -def viscoelastic_constitutive_law( - displacement, velocity, setup, elements, nodes, **_kwargs -): +def viscoelastic_constitutive_law(displacement, velocity, setup, elements, nodes, **_kwargs): grad_x = gradient(elements, nodes, displacement[:, 0]) grad_y = gradient(elements, nodes, displacement[:, 1]) grad_u = np.concatenate((grad_x, grad_y), axis=1).reshape(-1, 2, 2) diff --git a/tests/test_conmech/regression/std_boundary.py b/tests/test_conmech/regression/std_boundary.py index 35f6219b..1cf4d44f 100644 --- a/tests/test_conmech/regression/std_boundary.py +++ b/tests/test_conmech/regression/std_boundary.py @@ -60,9 +60,7 @@ def extract_boundary_paths_from_elements(elements): visited_path = extract_boundary_path(boundary_surfaces, start_node=start_node) visited_path = np.append(visited_path, visited_path[0]) boundary_paths.append(visited_path) - boundary_indices_to_visit = list( - set(boundary_indices_to_visit) - set(visited_path) - ) + boundary_indices_to_visit = list(set(boundary_indices_to_visit) - set(visited_path)) return boundary_paths diff --git a/tests/test_conmech/regression/test_density_update.py b/tests/test_conmech/regression/test_density_update.py index 2bdc7579..80afb3c1 100644 --- a/tests/test_conmech/regression/test_density_update.py +++ b/tests/test_conmech/regression/test_density_update.py @@ -168,9 +168,7 @@ def inner_forces(x, t=None): @staticmethod def outer_forces(x, t=None): - return ( - np.array([0, 0.1]) if x[1] < 0.2 and x[0] >= 2.0 else np.array([0, 0]) - ) + return np.array([0, 0.1]) if x[1] < 0.2 and x[0] >= 2.0 else np.array([0, 0]) boundaries: ... = BoundariesDescription( contact=lambda x: x[1] == 0 and x[0] <= 1.0, dirichlet=lambda x: x[0] == 0 @@ -238,12 +236,8 @@ def outer_forces(x, t=None): return test_suites -@pytest.mark.parametrize( - "setup, density_func, expected_displacement_vector", generate_test_suits() -) -def test_nonhomogenous_solver( - solving_method, setup, density_func, expected_displacement_vector -): +@pytest.mark.parametrize("setup, density_func, expected_displacement_vector", generate_test_suits()) +def test_nonhomogenous_solver(solving_method, setup, density_func, expected_displacement_vector): runner = NonHomogenousSolver(setup, solving_method) elem_centers = get_elem_centers(runner) elements_density = np.asarray([density_func(x) for x in elem_centers]) diff --git a/tests/test_conmech/regression/test_dynamic_membrane.py b/tests/test_conmech/regression/test_dynamic_membrane.py index 9cdfdd16..97f8f7b2 100644 --- a/tests/test_conmech/regression/test_dynamic_membrane.py +++ b/tests/test_conmech/regression/test_dynamic_membrane.py @@ -256,6 +256,4 @@ def test_dynamic_membrane(solving_method, setup, expected_displacement_vector): np.set_printoptions(precision=8, suppress=True) print(repr(displacement)) - np.testing.assert_array_almost_equal( - displacement, expected_displacement_vector, decimal=3 - ) + np.testing.assert_array_almost_equal(displacement, expected_displacement_vector, decimal=3) diff --git a/tests/test_conmech/regression/test_membrane_DNC.py b/tests/test_conmech/regression/test_membrane_DNC.py index 3d66c8cb..b26d04c2 100644 --- a/tests/test_conmech/regression/test_membrane_DNC.py +++ b/tests/test_conmech/regression/test_membrane_DNC.py @@ -196,6 +196,4 @@ def test_membrane_DNC(solving_method, setup, expected_displacement_vector): np.set_printoptions(precision=7, suppress=True) print(repr(displacement)) - np.testing.assert_array_almost_equal( - displacement, expected_displacement_vector, decimal=2 - ) + np.testing.assert_array_almost_equal(displacement, expected_displacement_vector, decimal=2) diff --git a/tests/test_conmech/regression/test_membrane_boundary_DNC.py b/tests/test_conmech/regression/test_membrane_boundary_DNC.py index 2e99222b..797717af 100644 --- a/tests/test_conmech/regression/test_membrane_boundary_DNC.py +++ b/tests/test_conmech/regression/test_membrane_boundary_DNC.py @@ -246,6 +246,4 @@ def test_boundary_DNC(solving_method, setup, expected_displacement_vector): print(repr(displacement)) precision = 2 if solving_method != "global optimization" else 3 - np.testing.assert_array_almost_equal( - displacement, expected_displacement_vector, decimal=2 - ) + np.testing.assert_array_almost_equal(displacement, expected_displacement_vector, decimal=2) diff --git a/tests/test_conmech/regression/test_piezoelectric_quasistatic.py b/tests/test_conmech/regression/test_piezoelectric_quasistatic.py index a7dae6a9..17417b96 100644 --- a/tests/test_conmech/regression/test_piezoelectric_quasistatic.py +++ b/tests/test_conmech/regression/test_piezoelectric_quasistatic.py @@ -25,9 +25,7 @@ def solving_method(request): class PPSlopeContactLaw(PotentialOfContactLaw): @staticmethod - def tangential_bound( - var_nu: float, static_displacement_nu: float, dt: float - ) -> float: + def tangential_bound(var_nu: float, static_displacement_nu: float, dt: float) -> float: return -1.0 @staticmethod @@ -72,9 +70,7 @@ class QuasistaticSetup_1(PiezoelectricQuasistaticProblem): ) ) permittivity: ... = field( - default_factory=lambda: np.array( - [[0.1, 0.0, 0.0], [0.0, 0.1, 0.0], [0.0, 0.0, 0.1]] - ) + default_factory=lambda: np.array([[0.1, 0.0, 0.0], [0.0, 0.1, 0.0], [0.0, 0.0, 0.1]]) ) @staticmethod @@ -203,12 +199,8 @@ def inner_forces(x, time=None): return np.array([0, -0.2]) setup_0_m02_p_0.inner_forces = inner_forces - expected_displacement_vector_0_m02_p_0 = [ - -v for v in expected_displacement_vector_0_02_p_0 - ] - expected_temperature_vector_0_m02_p_0 = [ - -v for v in expected_temperature_vector_0_02_p_0 - ] + expected_displacement_vector_0_m02_p_0 = [-v for v in expected_displacement_vector_0_02_p_0] + expected_temperature_vector_0_m02_p_0 = [-v for v in expected_temperature_vector_0_02_p_0] test_suites.append( ( setup_0_m02_p_0, @@ -237,9 +229,7 @@ class QuasistaticSetup_2(PiezoelectricQuasistaticProblem): ) ) permittivity: ... = field( - default_factory=lambda: np.array( - [[0.1, 0.0, 0.0], [0.0, 0.1, 0.0], [0.0, 0.0, 0.1]] - ) + default_factory=lambda: np.array([[0.1, 0.0, 0.0], [0.0, 0.1, 0.0], [0.0, 0.0, 0.1]]) ) @staticmethod @@ -324,9 +314,7 @@ def test_piezoelectric_time_dependent_solver( std_ids = standard_boundary_nodes(runner.body.mesh.nodes, runner.body.mesh.elements) displacement = results[-1].body.mesh.nodes[:] - results[-1].displaced_nodes[:] electric_potential = np.zeros(len(results[-1].body.mesh.nodes)) - electric_potential[: len(results[-1].electric_potential)] = results[ - -1 - ].electric_potential + electric_potential[: len(results[-1].electric_potential)] = results[-1].electric_potential # print result np.set_printoptions(precision=8, suppress=True) diff --git a/tests/test_conmech/regression/test_poisson.py b/tests/test_conmech/regression/test_poisson.py index b08f4317..911a5e9e 100644 --- a/tests/test_conmech/regression/test_poisson.py +++ b/tests/test_conmech/regression/test_poisson.py @@ -31,9 +31,7 @@ def internal_temperature(x: np.ndarray, t: float = None) -> np.ndarray: def outer_temperature(x: np.ndarray, t: float = None) -> np.ndarray: return np.array([3.0]) - boundaries: ... = BoundariesDescription( - dirichlet=lambda x: x[0] == 0 or x[0] == 1 - ) + boundaries: ... = BoundariesDescription(dirichlet=lambda x: x[0] == 0 or x[0] == 1) mesh_descr = CrossMeshDescription( initial_position=None, max_element_perimeter=0.25, scale=[1, 1] @@ -164,6 +162,4 @@ def test_poisson_solver(solving_method, setup, expected_temperature_vector): np.set_printoptions(precision=8, suppress=True) - np.testing.assert_array_almost_equal( - temperature, expected_temperature_vector, decimal=3 - ) + np.testing.assert_array_almost_equal(temperature, expected_temperature_vector, decimal=3) diff --git a/tests/test_conmech/regression/test_temperature_dynamic.py b/tests/test_conmech/regression/test_temperature_dynamic.py index 8bca6557..a7a4f014 100644 --- a/tests/test_conmech/regression/test_temperature_dynamic.py +++ b/tests/test_conmech/regression/test_temperature_dynamic.py @@ -59,14 +59,10 @@ class DynamicSetup(TemperatureDynamicProblem): contact_law: ... = make_slope_contact_law(1e1) contact_law_2: ... = TPSlopeContactLaw() thermal_expansion: ... = field( - default_factory=lambda: np.array( - [[0.5, 0.0, 0.0], [0.0, 0.5, 0.0], [0.0, 0.0, 0.5]] - ) + default_factory=lambda: np.array([[0.5, 0.0, 0.0], [0.0, 0.5, 0.0], [0.0, 0.0, 0.5]]) ) thermal_conductivity: ... = field( - default_factory=lambda: np.array( - [[0.1, 0.0, 0.0], [0.0, 0.1, 0.0], [0.0, 0.0, 0.1]] - ) + default_factory=lambda: np.array([[0.1, 0.0, 0.0], [0.0, 0.1, 0.0], [0.0, 0.0, 0.1]]) ) @staticmethod @@ -195,12 +191,8 @@ def inner_forces(x, time=None): return np.array([0, -0.2]) setup_0_m02_p_0.inner_forces = inner_forces - expected_displacement_vector_0_m02_p_0 = [ - -v for v in expected_displacement_vector_0_02_p_0 - ] - expected_temperature_vector_0_m02_p_0 = [ - -v for v in expected_temperature_vector_0_02_p_0 - ] + expected_displacement_vector_0_m02_p_0 = [-v for v in expected_displacement_vector_0_02_p_0] + expected_temperature_vector_0_m02_p_0 = [-v for v in expected_temperature_vector_0_02_p_0] test_suites.append( ( setup_0_m02_p_0, @@ -221,14 +213,10 @@ class DynamicSetup(TemperatureDynamicProblem): contact_law: ... = make_slope_contact_law(2.71) contact_law_2: ... = TPSlopeContactLaw thermal_expansion: ... = field( - default_factory=lambda: np.array( - [[0.5, 0.0, 0.0], [0.0, 0.5, 0.0], [0.0, 0.0, 0.5]] - ) + default_factory=lambda: np.array([[0.5, 0.0, 0.0], [0.0, 0.5, 0.0], [0.0, 0.0, 0.5]]) ) thermal_conductivity: ... = field( - default_factory=lambda: np.array( - [[0.1, 0.0, 0.0], [0.0, 0.1, 0.0], [0.0, 0.0, 0.1]] - ) + default_factory=lambda: np.array([[0.1, 0.0, 0.0], [0.0, 0.1, 0.0], [0.0, 0.0, 0.1]]) ) @staticmethod diff --git a/tests/test_conmech/smoke/test_examples.py b/tests/test_conmech/smoke/test_examples.py index 0c2b5714..06ec2470 100644 --- a/tests/test_conmech/smoke/test_examples.py +++ b/tests/test_conmech/smoke/test_examples.py @@ -36,24 +36,18 @@ "dynamic": lambda: dynamic(Config(**default_args).init()), "dynamic_piezo": lambda: dynamic_piezo(Config(**default_args).init()), "dynamic_temp": lambda: dynamic_temp(Config(**default_args).init()), - "Jureczka_and_Ochal_2019": lambda: Jureczka_and_Ochal_2019( - Config(**default_args).init() - ), + "Jureczka_and_Ochal_2019": lambda: Jureczka_and_Ochal_2019(Config(**default_args).init()), "Jureczka_Ochal_Bartman_2023": lambda: Jureczka_Ochal_Bartman_2023( Config(outputs_path="./output/JOB2023", **default_args).init() ), "Sofonea_Ochal_Bartman_2023": lambda: Sofonea_Ochal_Bartman_2023( Config(outputs_path="./output/SOB2023", **default_args).init() ), - "BOST_2024": lambda: BOST_2024( - Config(outputs_path="./output/BOST2024", **default_args).init() - ), + "BOST_2024": lambda: BOST_2024(Config(outputs_path="./output/BOST2024", **default_args).init()), "Tarzia_problem": lambda: Tarzia_problem( Config(outputs_path="./output/BOST2024", **default_args).init() ), - "nonhomogenous_density": lambda: nonhomogenous_density( - Config(**default_args).init() - ), + "nonhomogenous_density": lambda: nonhomogenous_density(Config(**default_args).init()), "imported_mesh": lambda: imported_mesh( Config(**default_args).init(), mesh_path="examples/meshes/example_mesh.msh", diff --git a/tests/test_conmech/unit_test/test_matrices.py b/tests/test_conmech/unit_test/test_matrices.py index bad8c8dd..e455a440 100644 --- a/tests/test_conmech/unit_test/test_matrices.py +++ b/tests/test_conmech/unit_test/test_matrices.py @@ -91,15 +91,10 @@ def test_local_stiff_mats_assembly(): scale=[scale_x, scale_y], ) ) - edges_features_matrix, _, local_stiff_mats = sut_2d( - elements=elements, nodes=initial_nodes - ) + edges_features_matrix, _, local_stiff_mats = sut_2d(elements=elements, nodes=initial_nodes) expected_w_matrix = np.asarray( [ - [ - edges_features_matrix[2 + dimension * (k + 1) + j] - for j in range(dimension) - ] + [edges_features_matrix[2 + dimension * (k + 1) + j] for j in range(dimension)] for k in range(dimension) ] ) diff --git a/tests/test_conmech/unit_test/test_mesh.py b/tests/test_conmech/unit_test/test_mesh.py index 9193b89d..820d1dfa 100644 --- a/tests/test_conmech/unit_test/test_mesh.py +++ b/tests/test_conmech/unit_test/test_mesh.py @@ -173,9 +173,7 @@ def is_contact(x): ) -unordered_nodes = np.asarray( - [[1.0, 1.0], [0.0, 0.0], [0.0, 2.0], [2.0, 2.0], [2.0, 0.0]] -) +unordered_nodes = np.asarray([[1.0, 1.0], [0.0, 0.0], [0.0, 2.0], [2.0, 2.0], [2.0, 0.0]]) unordered_elements = np.asarray([[1, 2, 0], [2, 3, 0], [3, 4, 0], [4, 1, 0]]) @@ -211,9 +209,7 @@ def test_condition_boundaries(_test_name_, params): expected_neumann_boundary, expected_dirichlet_boundary, ) = params - boundaries_description = BoundariesDescription( - contact=is_contact, dirichlet=is_dirichlet - ) + boundaries_description = BoundariesDescription(contact=is_contact, dirichlet=is_dirichlet) # Act ( @@ -228,9 +224,7 @@ def test_condition_boundaries(_test_name_, params): # Assert def unify_edges(boundary): - return frozenset( - [frozenset([str(np.sort(node)) for node in edge]) for edge in boundary] - ) + return frozenset([frozenset([str(np.sort(node)) for node in edge]) for edge in boundary]) def compare_surfaces(actual_surfaces, expected_surfaces): return unify_edges(initial_nodes[actual_surfaces]) == unify_edges( @@ -239,6 +233,4 @@ def compare_surfaces(actual_surfaces, expected_surfaces): assert compare_surfaces(boundaries_data.contact_boundary, expected_contact_boundary) assert compare_surfaces(boundaries_data.neumann_boundary, expected_neumann_boundary) - assert compare_surfaces( - boundaries_data.dirichlet_boundary, expected_dirichlet_boundary - ) + assert compare_surfaces(boundaries_data.dirichlet_boundary, expected_dirichlet_boundary) diff --git a/tests/test_conmech/unit_test/test_setting_mesh.py b/tests/test_conmech/unit_test/test_setting_mesh.py index 8f243d45..b7c85b59 100644 --- a/tests/test_conmech/unit_test/test_setting_mesh.py +++ b/tests/test_conmech/unit_test/test_setting_mesh.py @@ -71,9 +71,7 @@ def test_remove_unconnected_nodes(): elements = np.array([[4, 2], [2, 5], [4, 5]]) # Act - cleaned_nodes, cleaned_elements = mesh.remove_unconnected_nodes_numba( - nodes, elements - ) + cleaned_nodes, cleaned_elements = mesh.remove_unconnected_nodes_numba(nodes, elements) # Assert expected_nodes = np.array([[2.1, 1.2], [4.1, 1.2], [5.1, 1.2]])