1"""A design: the material, the bounding box, the outline, the thickness field and
2the stress-free profile of a plate, as numbers or vectorized functions of position.
4The functions take four ``(ny, nx)`` arrays, ``(x, y, r, theta)``: the point
5relative to the center of the bounding box, its distance from the center and
6its angle from the x axis, all in meters and radians, and return an array of
7the same shape (or a scalar). ``Design.evaluate`` samples them on a grid.
8"""
10from __future__ import annotations
12from dataclasses import dataclass
13from typing import Callable, Union
15import numpy as np
17from .constants import ROOM_TEMPERATURE
18from .grid import InsideFn, make_grid
19from .material import Material
20from .plate import Plate, remove_plane
22Field = Union[float, Callable[[np.ndarray, np.ndarray, np.ndarray, np.ndarray], np.ndarray]]
23"""A number, or a function ``f(x, y, r, theta)`` of grid matrices."""
26class DesignError(ValueError):
27 """A design that cannot be evaluated."""
31class Design:
32 """A plate design. Lengths are meters.
34 Example::
36 inch, mm = 0.0254, 1e-3
37 R = 8 * inch
38 crash = Design(
39 name="16 inch flat crash",
40 material=BRONZE,
41 width=2 * R, height=2 * R,
42 inside=lambda x, y, r, theta: r <= R,
43 thickness=lambda x, y, r, theta: 1.2 * mm + (0.7 * mm - 1.2 * mm) * np.clip((r - 0.04) / (R - 0.04), 0, 1),
44 )
45 plate = crash.evaluate(spacing=0.008)
46 """
48 name: str
49 material: Material
50 width: float
51 """Bounding box [m]."""
52 height: float
53 thickness: Field
54 """Thickness [m]: a number or a function of ``(x, y, r, theta)``; must be positive on the plate."""
55 inside: InsideFn | None = None
56 """Outline: ``inside(x, y, r, theta) -> bool``; ``None`` for the whole bounding box."""
57 profile: Field = 0.0
58 """Stress-free height [m]: 0 for a flat plate, or a function for a bow or a bell. Its mean plane is removed."""
59 temperature: float = ROOM_TEMPERATURE
60 """Plate temperature [K], for the thermoelastic damping."""
62 def __post_init__(self) -> None:
63 if not (self.width >= 1e-3 and self.height >= 1e-3):
64 raise DesignError("width and height must be at least 1 mm")
65 if not callable(self.thickness) and not float(self.thickness) > 0:
66 raise DesignError("thickness must be positive")
68 def evaluate(self, spacing: float = 0.008) -> Plate:
69 """Samples the design on a grid of about the given spacing: the mask from
70 the outline at the cell centers, the thickness and the profile at the
71 nodes. Only the slopes of the profile enter the physics, so its mean is
72 immaterial, but its tilt is not (the von Karman strains are not invariant
73 under rotation), so the mass-weighted least-squares plane is removed."""
74 g = make_grid(self.width, self.height, spacing, self.inside)
75 if g.count == 0:
76 raise DesignError("the outline contains no grid cells; enlarge it or reduce the grid spacing")
77 X, Y = g.X, g.Y
78 R, T = np.hypot(X, Y), np.arctan2(Y, X)
79 h = sample_field(self.thickness, "thickness", X, Y, R, T)
80 w0 = sample_field(self.profile, "profile", X, Y, R, T)
81 m = g.mask
82 bad = m & ~np.isfinite(h)
83 if bad.any():
84 j, i = np.argwhere(bad)[0]
85 raise DesignError(f"thickness is not finite at ({X[j, i]:.4f}, {Y[j, i]:.4f})")
86 bad = m & (h <= 0)
87 if bad.any():
88 j, i = np.argwhere(bad)[0]
89 raise DesignError(f"thickness must be positive; it is {h[j, i]} at ({X[j, i]:.4f}, {Y[j, i]:.4f})")
90 bad = m & ~np.isfinite(w0)
91 if bad.any():
92 j, i = np.argwhere(bad)[0]
93 raise DesignError(f"profile is not finite at ({X[j, i]:.4f}, {Y[j, i]:.4f})")
94 w0 = np.where(m, remove_plane(np.where(m, w0, 0.0), g), 0.0)
95 return Plate(g, self.material, np.where(m, h, 0.0), w0, self.width, self.height, self.temperature, self.name)
98def sample_field(f: Field, what: str, X: np.ndarray, Y: np.ndarray, R: np.ndarray, T: np.ndarray) -> np.ndarray:
99 """Evaluates a field (a number or a function of ``(x, y, r, theta)``) on grid matrices."""
100 if callable(f):
101 try:
102 v = np.asarray(f(X, Y, R, T), dtype=float)
103 except Exception as e: # noqa: BLE001
104 raise DesignError(f"{what}: {e}") from e
105 if v.ndim == 0:
106 return np.full(X.shape, float(v))
107 if v.shape != X.shape:
108 raise DesignError(f"{what} returned shape {v.shape} for a {X.shape[0]} x {X.shape[1]} grid; use elementwise numpy operations")
109 return v
110 try:
111 return np.full(X.shape, float(f))
112 except (TypeError, ValueError) as e:
113 raise DesignError(f"{what} must be a number or a function of (x, y, r, theta)") from e