Merge remote-tracking branch 'fea2/main'

Author: tomvanmele

src/compas_fea2/UI/__init__.py

@@ -1,7 +1,3 @@
-from __future__ import absolute_import
-from __future__ import division
-from __future__ import print_function
-
 from .viewer import FEA2Viewer
 
 __all__ = [

src/compas_fea2/UI/viewer/scene.py

@@ -1,7 +1,3 @@
-from __future__ import absolute_import
-from __future__ import division
-from __future__ import print_function
-
 from compas.colors import Color
 from compas.colors import ColorMap
 from compas.geometry import Vector

src/compas_fea2/cli.py

@@ -62,7 +62,7 @@ def change_setting(backend, setting, value):
             The new value for the setting.
 
     Example usage:\n
-        fea2 change-setting opensees exe "Applications/OpenSees3.5.0/bin/OpenSees"
+        fea2 change-setting opensees exe "Applications/OpenSees3.7.0/bin/OpenSees"
     """
     m = importlib.import_module("compas_fea2_" + backend.lower())
     env = os.path.join(m.HOME, "src", "compas_fea2_" + backend.lower(), ".env")

src/compas_fea2/model/elements.py

@@ -206,8 +206,19 @@ def rigid(self) -> bool:
     def mass(self) -> float:
         return self.volume * self.section.material.density
 
-    def weight(self, g: float) -> float:
-        return self.mass * g
+    @property
+    def g(self) -> float:
+        if self.model:
+            return self.model.g
+        else:
+            return self.g
+
+    @property
+    def weight(self) -> float:
+        if self.model.g:
+            return self.mass * self.g
+        else:
+            raise "Gravity constant not defined"
 
     @property
     def nodal_mass(self) -> List[float]:
@@ -356,9 +367,9 @@ def volume(self) -> float:
     def plot_section(self):
         self.section.plot()
 
-    def plot_stress_distribution(self, step: "Step", end: str = "end_1", nx: int = 100, ny: int = 100, *args, **kwargs):
-        """Plot the stress distribution along the element.
-
+    def plot_stress_distribution(self, step: "Step", end: str = "end_1", nx: int = 100, ny: int = 100, *args, **kwargs):  # noqa: F821
+        """ Plot the stress distribution along the element.
+        
         Parameters
         ----------
         step : :class:`compas_fea2.model.Step`
@@ -412,9 +423,9 @@ def forces(self, step: "Step") -> "Result":
         r = self.section_forces_result(step)
         return r.forces
 
-    def moments(self, step: "Step") -> "Result":
-        """Get the moments result for the element.
-
+    def moments(self, step: "_Step") -> "Result":  # noqa: F821
+        """ Get the moments result for the element.
+        
         Parameters
         ----------
         step : :class:`compas_fea2.model.Step`
@@ -491,8 +502,7 @@ def __init__(self, nodes: List["Node"], tag: str, element: Optional["_Element"]
         self._nodes = nodes
         self._tag = tag
         self._plane = Plane.from_three_points(*[node.xyz for node in nodes[:3]])  # TODO check when more than 3 nodes
-        self._registration = element
-        self._centroid = centroid_points([node.xyz for node in nodes])
+        self._registration = element  # FIXME: not updated when copying parts
 
     @property
     def __data__(self):
@@ -549,7 +559,7 @@ def area(self) -> float:
 
     @property
     def centroid(self) -> "Point":
-        return self._centroid
+        return centroid_points([node.xyz for node in self.nodes])
 
     @property
     def nodes_key(self) -> List:
@@ -841,7 +851,7 @@ def __init__(
         if len(nodes) not in {4, 10}:
             raise ValueError("TetrahedronElement must have either 4 (C3D4) or 10 (C3D10) nodes.")
 
-        self.element_type = "C3D10" if len(nodes) == 10 else "C3D4"
+        # self.element_type = "C3D10" if len(nodes) == 10 else "C3D4"
 
         super().__init__(
             nodes=nodes,

src/compas_fea2/model/groups.py

@@ -136,10 +136,6 @@ def group_by(self, key: Callable[[T], Any]) -> Dict[Any, "_Group"]:
             A dictionary where keys are the grouping values and values are `_Group` instances.
         """
         sorted_members = self._members
-        # try:
-        #     sorted_members = sorted(self._members, key=key)
-        # except TypeError:
-        #     sorted_members = sorted(self._members, key=lambda x: x.key)
         grouped_members = {k: set(v) for k, v in groupby(sorted_members, key=key)}
         return {k: self.__class__(v, name=f"{self.name}") for k, v in grouped_members.items()}
 

src/compas_fea2/model/materials/concrete.py

@@ -1,8 +1,5 @@
 from math import log
 
-from compas_fea2.units import UnitRegistry
-from compas_fea2.units import units as u
-
 from .material import _Material
 
 
@@ -139,11 +136,6 @@ def __from_data__(cls, data):
     # FIXME: this is only working for the basic material properties.
     @classmethod
     def C20_25(cls, units, **kwargs):
-        if not units:
-            units = u(system="SI_mm")
-        elif not isinstance(units, UnitRegistry):
-            units = u(units)
-
         return cls(fck=25 * units.MPa, E=30 * units.GPa, v=0.17, density=2400 * units("kg/m**3"), name="C20/25", **kwargs)
 
 

src/compas_fea2/model/materials/material.py

@@ -303,7 +303,6 @@ def __init__(self, *, density: float, expansion: Optional[float] = None, name: O
 # ==============================================================================
 # non-linear general
 # ==============================================================================
-# @extend_docstring(_Material)
 class ElasticPlastic(ElasticIsotropic):
     """Elastic and plastic, isotropic and homogeneous material.
 

src/compas_fea2/model/model.py

@@ -9,6 +9,7 @@
 from typing import Optional
 from typing import Set
 from typing import Union
+from typing import Tuple
 
 from compas.datastructures import Graph
 from compas.geometry import Box
@@ -18,6 +19,7 @@
 from compas.geometry import Transformation
 from compas.geometry import bounding_box
 from compas.geometry import centroid_points
+from compas.geometry import centroid_points_weighted
 from pint import UnitRegistry
 
 import compas_fea2
@@ -107,6 +109,8 @@ def __init__(self, description: Optional[str] = None, author: Optional[str] = No
         self._groups: Set[_Group] = set()
         self._problems: Set[Problem] = set()
 
+        self._constants: dict = {"g": None}
+
     @property
     def __data__(self):
         return {
@@ -121,6 +125,7 @@ def __data__(self):
             "sections": [section.__data__ for section in self.sections],
             "problems": [problem.__data__ for problem in self.problems],
             "path": str(self.path) if self.path else None,
+            "constants": self._constants,
         }
 
     @classmethod
@@ -161,6 +166,7 @@ def __from_data__(cls, data):
         problem_classes = {cls.__name__: cls for cls in Problem.__subclasses__()}
         model._problems = {problem_classes[problem_data["class"]].__from_data__(problem_data) for problem_data in data.get("problems", [])}
         model._path = Path(data.get("path")) if data.get("path") else None
+        model._constants = data.get("constants")
         return model
 
     @classmethod
@@ -217,6 +223,18 @@ def constraints(self) -> Set[_Constraint]:
     def connectors(self) -> Set[Connector]:
         return self._connectors
 
+    @property
+    def constants(self) -> dict:
+        return self._constants
+
+    @property
+    def g(self) -> float:
+        return self.constants["g"]
+
+    @g.setter
+    def g(self, value):
+        self._constants["g"] = value
+
     @property
     def materials_dict(self) -> dict[Union[_Part, "Model"], list[_Material]]:
         materials = {part: part.materials for part in filter(lambda p: not isinstance(p, RigidPart), self.parts)}
@@ -301,11 +319,24 @@ def bounding_box(self) -> Optional[Box]:
         return Box.from_bounding_box(bb)
 
     @property
-    def center(self) -> Point:
+    def bb_center(self) -> Point:
         if self.bounding_box:
-            return centroid_points(self.bounding_box.points)
+            return Point(*centroid_points(self.bounding_box.points))
         else:
-            return centroid_points(self.points)
+            raise AttributeError("The model has no bounding box")
+
+    @property
+    def center(self) -> Point:
+        return Point(*centroid_points(self.points))
+
+    @property
+    def centroid(self) -> Point:
+        weights = []
+        points = []
+        for part in self.parts:
+            points.append(part.centroid)
+            weights.append(part.weight)
+        return Point(*centroid_points_weighted(points=points, weights=weights))
 
     @property
     def bottom_plane(self) -> Plane:
@@ -855,7 +886,12 @@ def find_node_by_name(self, name: str) -> Node:
 
     @get_docstring(_Part)
     @part_method
-    def find_closest_nodes_to_node(self, node: Node, distance: float, number_of_nodes: int = 1, plane: Optional[Plane] = None) -> NodesGroup:
+    def find_closest_nodes_to_node(self, node: Node, number_of_nodes: int = 1, plane: Optional[Plane] = None) -> NodesGroup:
+        pass
+
+    @get_docstring(_Part)
+    @part_method
+    def find_closest_nodes_to_point(self, point: Point, number_of_nodes: int = 1, plane: Optional[Plane] = None) -> NodesGroup:
         pass
 
     @get_docstring(_Part)
@@ -1433,6 +1469,42 @@ def add_interfaces(self, interfaces):
         """
         return [self.add_interface(interface) for interface in interfaces]
 
+    def extract_interfaces(self, max_origin_distance=0.01, tol=1e-6) -> list[Tuple[Plane, Plane]]:
+        """Extract interfaces from the model.
+
+        Returns
+        -------
+        list[:class:`compas_fea2.model.Interface`]
+            List of interfaces extracted from the model.
+
+        """
+        from compas.geometry import distance_point_plane_signed
+        import itertools
+        from typing import List
+
+        # --- Nested helper function to check coplanarity ---
+        def _are_two_planes_coplanar(pln1: Plane, pln2: Plane, tol: float) -> bool:
+            """Checks if two COMPAS planes are geometrically the same."""
+            normals_collinear = abs(abs(pln1.normal.dot(pln2.normal)) - 1.0) < tol
+            if not normals_collinear:
+                return False
+            point_on_plane = abs(distance_point_plane_signed(pln1.point, pln2)) < tol
+            return point_on_plane
+
+        list_of_plane_lists = [p.extract_clustered_planes() for p in self.parts]
+        # 1. Flatten the list of lists into a single list of all planes.
+        all_planes_raw: List[Plane] = list(itertools.chain.from_iterable(list_of_plane_lists))
+
+        # 2. Get unique plane instances from the raw list.
+        coplanar_pairs: List[Tuple[Plane, Plane]] = []
+        for p1, p2 in itertools.combinations(all_planes_raw, 2):
+            if _are_two_planes_coplanar(p1, p2, tol):
+                # If planes are coplanar, then check the distance between their origin points
+                if max_origin_distance is None or p1.point.distance_to_point(p2.point) < max_origin_distance:
+                    coplanar_pairs.append((p1, p2))
+
+        return coplanar_pairs
+
     # ==============================================================================
     # Summary
     # ==============================================================================