Coverage for /opt/hostedtoolcache/Python/3.11.10/x64/lib/python3.11/site-packages/hypervehicle/geometry/surface.py: 43%
140 statements
« prev ^ index » next coverage.py v7.6.4, created at 2024-10-29 02:51 +0000
« prev ^ index » next coverage.py v7.6.4, created at 2024-10-29 02:51 +0000
1import numpy as np
2from copy import deepcopy
3from typing import Optional
4from abc import ABC, abstractmethod
5from hypervehicle.geometry.path import ClusterFunction, LinearFunction, Line, Path
6from hypervehicle.geometry.vector import Vector3, approximately_equal_vectors
9class ParametricSurface(ABC):
10 """Astract parametric surface."""
12 @abstractmethod
13 def __repr__(self):
14 pass
16 @abstractmethod
17 def __call__(self, r, s):
18 pass
21class CoonsPatch(ParametricSurface):
22 """A surface constructed by transfinite interpolation of the edges."""
24 _slots_ = [
25 "north",
26 "east",
27 "south",
28 "west",
29 "p00",
30 "p10",
31 "p11",
32 "p01",
33 "defined_by_corners",
34 "offset",
35 ]
37 def __init__(
38 self,
39 north: Path = None,
40 east: Path = None,
41 south: Path = None,
42 west: Path = None,
43 p00: Vector3 = None,
44 p10: Vector3 = None,
45 p11: Vector3 = None,
46 p01: Vector3 = None,
47 offset=Vector3(0, 0, 0),
48 ):
49 """Initialise from edges or corner points."""
50 if all([north, east, south, west]):
51 self.north = deepcopy(north)
52 self.east = deepcopy(east)
53 self.south = deepcopy(south)
54 self.west = deepcopy(west)
55 self.p00 = self.south(0.0)
56 self.p10 = self.south(1.0)
57 self.p01 = self.north(0.0)
58 self.p11 = self.north(1.0)
59 p00_alt = self.west(0.0)
60 p10_alt = self.east(0.0)
61 p01_alt = self.west(1.0)
62 p11_alt = self.east(1.0)
63 if not approximately_equal_vectors(self.p00, p00_alt):
64 raise Exception(f"CoonsPatch open corner p00={self.p00}, {p00_alt}")
65 if not approximately_equal_vectors(self.p10, p10_alt):
66 raise Exception(f"CoonsPatch open corner p10={self.p10}, {p10_alt}")
67 if not approximately_equal_vectors(self.p01, p01_alt):
68 raise Exception(f"CoonsPatch open corner p01={self.p01}, {p01_alt}")
69 if not approximately_equal_vectors(self.p11, p11_alt):
70 raise Exception(f"CoonsPatch open corner p11={self.p11}, {p11_alt}")
71 self.defined_by_corners = False
73 elif all([p00, p10, p11, p01]):
74 self.north = Line(p01, p11)
75 self.east = Line(p10, p11)
76 self.south = Line(p00, p10)
77 self.west = Line(p00, p01)
78 self.p00 = deepcopy(p00)
79 self.p10 = deepcopy(p10)
80 self.p11 = deepcopy(p11)
81 self.p01 = deepcopy(p01)
82 self.defined_by_corners = True
84 else:
85 raise Exception(
86 "CoonsPatch should be defined by four edge paths or four corners."
87 )
88 self.offset = offset
90 def __repr__(self):
91 str = "CoonsPatch("
92 if self.defined_by_corners:
93 str += f"p00={self.p00}, p10={self.p10}, p11={self.p10}, p01={self.p10}"
94 else:
95 str += f"north={self.north}, east={self.east}, south={self.south}, west={self.west}"
96 str += f", offset={self.offset})"
97 return str
99 def __call__(self, r, s):
100 """
101 Transfinite interpolation to an interior point, p.
102 """
103 south_r = self.south(r)
104 north_r = self.north(r)
105 west_s = self.west(s)
106 east_s = self.east(s)
107 p = (
108 south_r * (1.0 - s)
109 + north_r * s
110 + west_s * (1.0 - r)
111 + east_s * r
112 - (
113 self.p00 * (1.0 - r) * (1.0 - s)
114 + self.p01 * (1.0 - r) * s
115 + self.p10 * r * (1.0 - s)
116 + self.p11 * r * s
117 )
118 + self.offset
119 )
120 return p
122 def __add__(self, offset):
123 """
124 Returns a copy of the original, displaced by a Vector3 object.
125 """
126 if not isinstance(offset, Vector3):
127 raise Exception(f"Cannot add a {type(offset)} to a CoonsPatch.")
128 new_patch = deepcopy(self)
129 new_patch.offset += offset
130 return new_patch
132 def __sub__(self, offset):
133 """
134 Returns a copy of the original, displaced by a Vector3 object.
135 """
136 if not isinstance(offset, Vector3):
137 raise Exception(f"Cannot subtract a {type(offset)} from a CoonsPatch.")
138 new_patch = deepcopy(self)
139 new_patch.offset -= offset
140 return new_patch
143class StructuredGrid:
144 """Structured grid."""
146 _slots_ = ["dimensions", "niv", "njv", "nkv", "vertices", "label", "tags"]
148 def __init__(
149 self,
150 surface: ParametricSurface,
151 niv: Optional[int] = 1,
152 njv: Optional[int] = 1,
153 cf_list: Optional[list[ClusterFunction]] = None,
154 tags: Optional[list[str]] = None,
155 label: Optional[str] = "no label",
156 ):
157 cf_list = cf_list if cf_list is not None else [None, None, None, None]
158 cf_list = [
159 cf if isinstance(cf, ClusterFunction) else LinearFunction()
160 for cf in cf_list
161 ]
162 self.generate(surface, niv, njv, cf_list)
163 self.tags = (tags if tags is not None else ["", "", "", ""],)
164 self.label = label
166 def __repr__(self):
167 str = "StructuredGrid("
168 str += f"dimensions={self.dimensions}, niv={self.niv}, njv={self.njv}, nkv={self.nkv}"
169 str += f", vertices={self.vertices}"
170 str += f", tags={self.tags}"
171 str += ")"
172 return str
174 def generate(self, surface, niv, njv, cf_list):
175 if not isinstance(surface, ParametricSurface):
176 raise Exception("Need to supply a ParametricSurface to construct the grid.")
177 if niv < 2:
178 raise Exception(f"niv is too small: {niv}")
179 if njv < 2:
180 raise Exception(f"njv is too small: {njv}")
181 self.niv = niv
182 self.njv = njv
183 self.nkv = 1
184 self.dimensions = 2
186 # Start with uniformly-distributed sample points.
187 r = np.fromfunction(lambda i, j: i, (niv, njv), dtype=float) / (niv - 1)
188 s = np.fromfunction(lambda i, j: j, (niv, njv), dtype=float) / (njv - 1)
190 # Compute independent cluster function along each edge.
191 rNorth = cf_list[0](r)
192 sEast = cf_list[1](s)
193 rSouth = cf_list[2](r)
194 sWest = cf_list[3](s)
196 # Blend the clustered sample points from each edge of the rs unit square.
197 sdash = (1.0 - r) * sWest + r * sEast
198 rdash = (1.0 - s) * rSouth + s * rNorth
200 # Compute the xyz spatial coordinates of the surface.
201 self.vertices = surface(rdash, sdash)
202 return
204 def subgrid(self, i0=0, j0=0, k0=0, niv=1, njv=1, nkv=1):
205 """
206 Returns a copy of a subgrid of vertices.
207 Start at vertex i0,j0,k0 and extend for niv,njv,nkv vertices.
208 """
209 new_grid = StructuredGrid(empty=True)
210 if self.dimensions == 3:
211 new_xs = self.vertices.x[i0 : i0 + niv, j0 : j0 + njv, k0 : k0 + nkv].copy()
212 new_ys = self.vertices.y[i0 : i0 + niv, j0 : j0 + njv, k0 : k0 + nkv].copy()
213 new_zs = self.vertices.z[i0 : i0 + niv, j0 : j0 + njv, k0 : k0 + nkv].copy()
214 new_grid.vertices = Vector3(new_xs, new_ys, new_zs)
215 elif self.dimensions == 2:
216 new_xs = self.vertices.x[i0 : i0 + niv, j0 : j0 + njv].copy()
217 new_ys = self.vertices.y[i0 : i0 + niv, j0 : j0 + njv].copy()
218 new_zs = self.vertices.z[i0 : i0 + niv, j0 : j0 + njv].copy()
219 new_grid.vertices = Vector3(new_xs, new_ys, new_zs)
220 elif self.dimensions == 1:
221 new_xs = self.vertices.x[i0 : i0 + niv].copy()
222 new_ys = self.vertices.y[i0 : i0 + niv].copy()
223 new_zs = self.vertices.z[i0 : i0 + niv].copy()
224 new_grid.vertices = Vector3(new_xs, new_ys, new_zs)
225 new_grid.niv = niv
226 new_grid.njv = njv
227 new_grid.nkv = nkv
228 if niv == 1 and njv == 1 and nkv == 1:
229 new_grid.dimensions = 0
230 elif njv == 1 and nkv == 1:
231 new_grid.dimensions = 1
232 elif nkv == 1:
233 new_grid.dimensions = 2
234 else:
235 new_grid.dimensions = 3
236 return new_grid