Coverage for src/meshpy/core/nurbs

1# The MIT License (MIT)

5# Permission is hereby granted, free of charge, to any person obtaining a copy

6# of this software and associated documentation files (the "Software"), to deal

7# in the Software without restriction, including without limitation the rights

8# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell

9# copies of the Software, and to permit persons to whom the Software is

10# furnished to do so, subject to the following conditions:

11#

12# The above copyright notice and this permission notice shall be included in

13# all copies or substantial portions of the Software.

14#

15# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR

16# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,

17# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE

18# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER

19# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,

20# OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN

21# THE SOFTWARE.

22"""This module implements NURBS patches for the mesh."""

24import numpy as _np

26from meshpy.core.conf import mpy as _mpy

27from meshpy.core.element import Element as _Element

28from meshpy.core.material import (

29 MaterialSolidBase as _MaterialSolidBase,

30)

33class NURBSPatch(_Element):

34 """A base class for a NURBS patch."""

36 # A list of valid material types for this element

37 valid_material = [_MaterialSolidBase]

39 def __init__(

40 self,

41 knot_vectors,

42 polynomial_orders,

43 element_string,

44 material=None,

45 nodes=None,

46 element_description=None,

47 ):

48 super().__init__(nodes=nodes, material=material)

50 # Knot vectors

51 self.knot_vectors = knot_vectors

53 # Polynomial degrees

54 self.polynomial_orders = polynomial_orders

56 # Set numbers for elements

57 self.n_nurbs_patch = None

59 # Set the element definitions

60 self.element_string = element_string

61 self.element_description = element_description

63 def get_nurbs_dimension(self):

64 """Return the number of dimensions of the NURBS structure."""

65 n_knots = len(self.knot_vectors)

66 n_polynomial = len(self.polynomial_orders)

67 if not n_knots == n_polynomial:

68 raise ValueError(

69 "The variables n_knots and polynomial_orders should have "

70 f"the same length. Got {n_knots} and {n_polynomial}"

71 )

72 return n_knots

74 def dump_element_specific_section(self, input_file):

75 """Set the knot vectors of the NURBS patch in the input file."""

77 patch_data = {

78 "knot_vectors": [],

79 }

81 for dir_manifold in range(len(self.knot_vectors)):

82 num_knotvectors = len(self.knot_vectors[dir_manifold])

84 # Check the type of knot vector, in case that the multiplicity of the first and last

85 # knot vectors is not p + 1, then it is a closed (periodic) knot vector, otherwise it

86 # is an open (interpolated) knot vector.

87 knotvector_type = "Interpolated"

89 for i in range(self.polynomial_orders[dir_manifold] - 1):

90 if (

91 abs(

92 self.knot_vectors[dir_manifold][i]

93 - self.knot_vectors[dir_manifold][i + 1]

94 )

95 > _mpy.eps_knot_vector

96 ) or (

97 abs(

98 self.knot_vectors[dir_manifold][num_knotvectors - 2 - i]

99 - self.knot_vectors[dir_manifold][num_knotvectors - 1 - i]

100 )

101 > _mpy.eps_knot_vector

102 ):

103 knotvector_type = "Periodic"

104 break

105

106 patch_data["knot_vectors"].append(

107 {

108 "DEGREE": self.polynomial_orders[dir_manifold],

109 "TYPE": knotvector_type,

110 "knots": [

111 knot_vector_val

112 for knot_vector_val in self.knot_vectors[dir_manifold]

113 ],

114 }

115 )

116

117 if "STRUCTURE KNOTVECTORS" not in input_file:

118 input_file.add({"STRUCTURE KNOTVECTORS": []})

119 input_file["STRUCTURE KNOTVECTORS"] = []

120

121 patches = input_file["STRUCTURE KNOTVECTORS"]

122 patch_data["ID"] = len(patches) + 1

123 patches.append(patch_data)

124

125 def get_number_elements(self):

126 """Get the number of elements in this patch by checking the amount of

127 nonzero knot spans in the knot vector."""

128

129 num_elements_dir = _np.zeros(len(self.knot_vectors), dtype=int)

130

131 for i_dir in range(len(self.knot_vectors)):

132 for i_knot in range(len(self.knot_vectors[i_dir]) - 1):

133 if (

134 abs(

135 self.knot_vectors[i_dir][i_knot]

136 - self.knot_vectors[i_dir][i_knot + 1]

137 )

138 > _mpy.eps_knot_vector

139 ):

140 num_elements_dir[i_dir] += 1

141

142 total_num_elements = _np.prod(num_elements_dir)

143

144 return total_num_elements

145

146 def _check_material(self):

147 """Check if the linked material is valid for this type of NURBS solid

148 element."""

149 for material_type in type(self).valid_material:

150 if isinstance(self.material, material_type):

151 return

152 raise TypeError(

153 f"NURBS solid of type {type(self)} can not have a material of"

154 f" type {type(self.material)}!"

155 )

156

157

158class NURBSSurface(NURBSPatch):

159 """A patch of a NURBS surface."""

160

161 def __init__(self, *args, element_string=None, **kwargs):

162 if element_string is None:

163 element_string = "WALLNURBS"

164 super().__init__(*args, element_string, **kwargs)

165

166 def dump_to_list(self):

167 """Return a list with all the element definitions contained in this

168 patch."""

169

170 # Check the material

171 self._check_material()

172

173 # Calculate how many control points are on the u direction

174 ctrlpts_size_u = len(self.knot_vectors[0]) - self.polynomial_orders[0] - 1

175

176 def get_ids_ctrlpts_surface(knot_span_u, knot_span_v):

177 """For an interpolated patch, calculate control points involved in

178 evaluation of the surface point at the knot span (knot_span_u,

179 knot_span_v)"""

180

181 id_u = knot_span_u - self.polynomial_orders[0]

182 id_v = knot_span_v - self.polynomial_orders[1]

183

184 element_ctrlpts_ids = []

185 for j in range(self.polynomial_orders[1] + 1):

186 for i in range(self.polynomial_orders[0] + 1):

187 # Calculating the global index of the control point, based on the book

188 # "Isogeometric Analysis: toward Integration of CAD and FEA" of J. Austin

189 # Cottrell, p. 314.

190 index_global = ctrlpts_size_u * (id_v + j) + id_u + i

191 element_ctrlpts_ids.append(index_global)

192

193 return element_ctrlpts_ids

194

195 patch_elements = []

196

197 # Adding an increment j to self.global to obtain the ID of an element in the patch

198 j = 0

199

200 # Loop over the knot spans to obtain the elements inside the patch

201 for knot_span_v in range(

202 self.polynomial_orders[1],

203 len(self.knot_vectors[1]) - self.polynomial_orders[1] - 1,

204 ):

205 for knot_span_u in range(

206 self.polynomial_orders[0],

207 len(self.knot_vectors[0]) - self.polynomial_orders[0] - 1,

208 ):

209 element_cps_ids = get_ids_ctrlpts_surface(knot_span_u, knot_span_v)

210

211 connectivity = [self.nodes[i].i_global for i in element_cps_ids]

212

213 num_cp_in_element = (self.polynomial_orders[0] + 1) * (

214 self.polynomial_orders[1] + 1

215 )

216

217 patch_elements.append(

218 {

219 "id": self.i_global + j,

220 "cell": {

221 "type": f"NURBS{num_cp_in_element}",

222 "connectivity": connectivity,

223 },

224 "data": {

225 "type": "WALLNURBS",

226 "MAT": self.material.i_global,

227 **(

228 self.element_description

229 if self.element_description

230 else {}

231 ),

232 },

233 }

234 )

235 j += 1

236

237 return patch_elements

238

239

240class NURBSVolume(NURBSPatch):

241 """A patch of a NURBS volume."""

242

243 def __init__(self, *args, element_string=None, **kwargs):

244 if element_string is not None:

245 raise ValueError("element_string is not yet implemented for NURBS volumes")

246 super().__init__(*args, element_string, **kwargs)

247

248 def dump_to_list(self):

249 """Return a list with all the element definitions contained in this

250 patch."""

251

252 # Check the material

253 self._check_material()

254

255 # Calculate how many control points are on the u and v directions

256 ctrlpts_size_u = len(self.knot_vectors[0]) - self.polynomial_orders[0] - 1

257 ctrlpts_size_v = len(self.knot_vectors[1]) - self.polynomial_orders[1] - 1

258

259 def get_ids_ctrlpts_volume(knot_span_u, knot_span_v, knot_span_w):

260 """For an interpolated patch, calculate control points involved in

261 evaluation of the surface point at the knot span (knot_span_u,

262 knot_span_v, knot_span_w)"""

263

264 id_u = knot_span_u - self.polynomial_orders[0]

265 id_v = knot_span_v - self.polynomial_orders[1]

266 id_w = knot_span_w - self.polynomial_orders[2]

267

268 element_ctrlpts_ids = []

269

270 for k in range(self.polynomial_orders[2] + 1):

271 for j in range(self.polynomial_orders[1] + 1):

272 for i in range(self.polynomial_orders[0] + 1):

273 # Calculating the global index of the control point, based on the paper

274 # "Isogeometric analysis: an overview and computer implementation aspects"

275 # of Vinh-Phu Nguyen, Mathematics and Computers in Simulation, Jun-2015.

276 index_global = (

277 ctrlpts_size_u * ctrlpts_size_v * (id_w + k)

278 + ctrlpts_size_u * (id_v + j)

279 + id_u

280 + i

281 )

282 element_ctrlpts_ids.append(index_global)

283

284 return element_ctrlpts_ids

285

286 patch_elements = []

287

288 # Adding an increment to self.global to obtain the ID of an element in the patch

289 increment_ele = 0

290

291 # Loop over the knot spans to obtain the elements inside the patch

292 for knot_span_w in range(

293 self.polynomial_orders[2],

294 len(self.knot_vectors[2]) - self.polynomial_orders[2] - 1,

295 ):

296 for knot_span_v in range(

297 self.polynomial_orders[1],

298 len(self.knot_vectors[1]) - self.polynomial_orders[1] - 1,

299 ):

300 for knot_span_u in range(

301 self.polynomial_orders[0],

302 len(self.knot_vectors[0]) - self.polynomial_orders[0] - 1,

303 ):

304 element_cps_ids = get_ids_ctrlpts_volume(

305 knot_span_u, knot_span_v, knot_span_w

306 )

307

308 connectivity = [self.nodes[i].i_global for i in element_cps_ids]

309

310 num_cp_in_element = (

311 (self.polynomial_orders[0] + 1)

312 * (self.polynomial_orders[1] + 1)

313 * (self.polynomial_orders[2] + 1)

314 )

315

316 patch_elements.append(

317 {

318 "id": self.i_global + increment_ele,

319 "cell": {

320 "type": f"NURBS{num_cp_in_element}",

321 "connectivity": connectivity,

322 },

323 "data": {

324 "type": "SOLID",

325 "MAT": self.material.i_global,

326 **(

327 self.element_description

328 if self.element_description

329 else {}

330 ),

331 },

332 }

333 )

334 increment_ele += 1

335

336 return patch_elements

Coverage for src/meshpy/core/nurbs_patch.py: 92%

106 statements