IgANet
IgANets - Isogeometric Analysis Networks
Loading...
Searching...
No Matches
unittest_bsplinelib.hpp
Go to the documentation of this file.
1
15#include <iganet.h>
16#include <iostream>
17
18#if __clang_major__ > 14
19#pragma clang diagnostic push
20#pragma clang diagnostic ignored "-Wunqualified-std-cast-call"
21#endif
22
23#include <BSplineLib/Splines/b_spline.hpp>
24
25#if __clang_major__ > 14
26#pragma clang diagnostic pop
27#endif
28
29#include <gtest/gtest.h>
30
31#pragma once
32
33template <typename Spline> auto to_bsplinelib_bspline(const Spline &bspline) {
34 static_assert(Spline::geoDim() < 5, "Unsupported geometric dimension");
35
36 // B-spline construction
37 using BSpline = bsplinelib::splines::BSpline<Spline::parDim()>;
38 using ParameterSpace = typename BSpline::ParameterSpace_;
39 using VectorSpace = typename BSpline::VectorSpace_;
40 using Coordinates = typename VectorSpace::Coordinates_;
41 using Degrees = typename ParameterSpace::Degrees_;
42 using Degree = typename Degrees::value_type;
43 using KnotVectors = typename ParameterSpace::KnotVectors_;
44 using KnotVector = typename KnotVectors::value_type::element_type;
45 using Knots = typename KnotVector::Knots_;
46 using Knot = typename Knots::value_type;
47
48 // B-spline evaluation
49 using ParametricCoordinate = typename BSpline::ParametricCoordinate_;
50 using ScalarParametricCoordinate = typename ParametricCoordinate::value_type;
51 using Derivative = typename ParameterSpace::Derivative_;
52 using ScalarDerivative = typename Derivative::value_type;
53
54 // Create degress structure
55 Degrees degrees;
56 for (iganet::short_t k = 0; k < bspline.parDim(); ++k)
57 degrees[k] = Degree{bspline.degree(k)};
58
59 // Create knot vectors
60 KnotVectors knot_vectors;
61 for (iganet::short_t k = 0; k < bspline.parDim(); ++k) {
62 Knots knots_;
63 for (int64_t i = 0; i < bspline.nknots(k); ++i)
64 knots_.emplace_back(Knot{
65 bspline.knots(k)[i].template item<typename Spline::value_type>()});
66 bsplinelib::SharedPointer<KnotVector> knot_vector{
67 std::make_shared<KnotVector>(knots_)};
68 knot_vectors[k] = std::move(knot_vector);
69 }
70
71 // Create parametric space
72 bsplinelib::SharedPointer<ParameterSpace> parameter_space{
73 std::make_shared<ParameterSpace>(knot_vectors, degrees)};
74
75 // Create coordinate vector(s)
76 Coordinates coordinates(bspline.ncumcoeffs(), Spline::geoDim());
77 for (int64_t i = 0; i < bspline.ncumcoeffs(); ++i) {
78 for (int64_t j = 0; j < Spline::geoDim(); ++j) {
79 coordinates(i, j) =
80 bspline.coeffs(j)[i].template item<typename Spline::value_type>();
81 }
82 }
83
84 // Create vector space
85 bsplinelib::SharedPointer<VectorSpace> vector_space{
86 std::make_shared<VectorSpace>(std::move(coordinates))};
87
88 // Create B-Spline
89 BSpline bsplinelib_bspline{parameter_space, vector_space};
90
91 return bsplinelib_bspline;
92}
93
94template <iganet::deriv deriv, bool memory_optimized, bool precompute,
95 typename Spline, typename BSplineLibSpline, typename TensorArray_t>
96void test_bspline_eval(const Spline &bspline,
97 BSplineLibSpline bsplinelib_bspline,
98 const TensorArray_t &xi,
99 typename Spline::value_type tol = 1e-12) {
100 static_assert(Spline::parDim() < 5, "Unsupported parametric dimension");
101
102 // B-spline evaluation
103 using Derivative = typename BSplineLibSpline::ParameterSpace_::Derivative_;
104 using ScalarDerivative = typename Derivative::value_type;
105 using Coordinate = typename BSplineLibSpline::Coordinate_;
106 using ScalarCoordinate = typename Coordinate::value_type;
107
108 using BSplineValue_type =
109 iganet::utils::BlockTensor<torch::Tensor, 1, Spline::geoDim()>;
110
111 BSplineValue_type bspline_val;
112 if constexpr (precompute) {
113 auto knot_indices = bspline.find_knot_indices(xi);
114 auto coeff_indices = bspline.find_coeff_indices(knot_indices);
115 auto basfunc = bspline.template eval_basfunc<deriv>(xi, knot_indices);
116 bspline_val = bspline.eval_from_precomputed(basfunc, coeff_indices,
117 xi[0].numel(), xi[0].sizes());
118 } else
119 bspline_val = bspline.template eval<deriv, memory_optimized>(xi);
120
121 Coordinate bsplinelib_query(Spline::parDim());
122 auto query = [&](const int64_t &i) -> ScalarCoordinate * {
123 for (iganet::short_t j = 0; j < Spline::parDim(); ++j) {
124 bsplinelib_query[j] =
125 (xi[j])[i].template item<typename Spline::value_type>();
126 }
127 return bsplinelib_query.data();
128 };
129
130 auto ten_pow = [](const iganet::short_t &p) {
131 iganet::short_t value{1};
132 for (iganet::short_t i = 0; i < p; ++i) {
133 value *= 10;
134 }
135 return value;
136 };
137
138 const Derivative bsplinelib_deriv = [&ten_pow]() {
139 Derivative d_query;
140 for (iganet::short_t i = 0; i < Spline::parDim(); ++i) {
141 d_query[i] = ((iganet::short_t)deriv / ((i == 0) ? 1 : ten_pow(i))) % 10;
142 }
143 return d_query;
144 }();
145
146 Coordinate bsplinelib_val(Spline::geoDim());
147 for (int64_t i = 0; i < xi[0].size(0); ++i) {
148 bsplinelib_bspline.EvaluateDerivative(query(i), bsplinelib_deriv.data(),
149 bsplinelib_val.data());
150 for (iganet::short_t j = 0; j < Spline::geoDim(); ++j) {
151 EXPECT_NEAR(
152 bspline_val(j)[i].template item<typename Spline::value_type>(),
153 bsplinelib_val[j], tol);
154 }
155 }
156}
157
158template <bool memory_optimized, bool precompute, typename Spline,
159 typename TensorArray_t>
160void test_bspline_grad(const Spline &bspline, const TensorArray_t &xi,
161 typename Spline::value_type tol = 1e-12) {
162 iganet::utils::BlockTensor<torch::Tensor, 1, Spline::parDim()>
163 bspline_grad_val;
164
165 if constexpr (precompute) {
166 auto knot_indices = bspline.find_knot_indices(xi);
167 auto coeff_indices = bspline.find_coeff_indices(knot_indices);
168 bspline_grad_val = bspline.grad(xi, knot_indices, coeff_indices);
169 } else
170 bspline_grad_val = bspline.grad(xi);
171
172 if constexpr (Spline::parDim() == 1) {
173 EXPECT_TRUE(torch::allclose(
174 bspline_grad_val(0, 0),
175 bspline.template eval<iganet::deriv::dx, memory_optimized>(xi)(0)));
176 } else if constexpr (Spline::parDim() == 2) {
177 EXPECT_TRUE(torch::allclose(
178 bspline_grad_val(0, 0),
179 bspline.template eval<iganet::deriv::dx, memory_optimized>(xi)(0)));
180 EXPECT_TRUE(torch::allclose(
181 bspline_grad_val(0, 1),
182 bspline.template eval<iganet::deriv::dy, memory_optimized>(xi)(0)));
183 } else if constexpr (Spline::parDim() == 3) {
184 EXPECT_TRUE(torch::allclose(
185 bspline_grad_val(0, 0),
186 bspline.template eval<iganet::deriv::dx, memory_optimized>(xi)(0)));
187 EXPECT_TRUE(torch::allclose(
188 bspline_grad_val(0, 1),
189 bspline.template eval<iganet::deriv::dy, memory_optimized>(xi)(0)));
190 EXPECT_TRUE(torch::allclose(
191 bspline_grad_val(0, 2),
192 bspline.template eval<iganet::deriv::dz, memory_optimized>(xi)(0)));
193 } else if constexpr (Spline::parDim() == 4) {
194 EXPECT_TRUE(torch::allclose(
195 bspline_grad_val(0, 0),
196 bspline.template eval<iganet::deriv::dx, memory_optimized>(xi)(0)));
197 EXPECT_TRUE(torch::allclose(
198 bspline_grad_val(0, 1),
199 bspline.template eval<iganet::deriv::dy, memory_optimized>(xi)(0)));
200 EXPECT_TRUE(torch::allclose(
201 bspline_grad_val(0, 2),
202 bspline.template eval<iganet::deriv::dz, memory_optimized>(xi)(0)));
203 EXPECT_TRUE(torch::allclose(
204 bspline_grad_val(0, 3),
205 bspline.template eval<iganet::deriv::dt, memory_optimized>(xi)(0)));
206 }
207}
208
209template <bool memory_optimized, bool precompute, typename Geometry_t,
210 typename Spline, typename TensorArray_t>
211void test_bspline_igrad(const Geometry_t &geometry, const Spline &bspline,
212 const TensorArray_t &xi,
213 typename Spline::value_type tol = 1e-12) {
214 iganet::utils::BlockTensor<torch::Tensor, 1, Spline::parDim()>
215 bspline_igrad_val;
216
217 if constexpr (precompute) {
218 auto knot_indices = bspline.find_knot_indices(xi);
219 auto coeff_indices = bspline.find_coeff_indices(knot_indices);
220 auto knot_indices_G = geometry.find_knot_indices(xi);
221 auto coeff_indices_G = geometry.find_coeff_indices(knot_indices_G);
222 bspline_igrad_val = bspline.igrad(geometry, xi, knot_indices, coeff_indices,
223 knot_indices_G, coeff_indices_G);
224 } else
225 bspline_igrad_val = bspline.igrad(geometry, xi);
226
227 if constexpr (Spline::parDim() == 1) {
228 EXPECT_TRUE(torch::allclose(
229 bspline_igrad_val(0, 0),
230 bspline.template eval<iganet::deriv::dx, memory_optimized>(xi)(0)));
231 } else if constexpr (Spline::parDim() == 2) {
232 EXPECT_TRUE(torch::allclose(
233 bspline_igrad_val(0, 0),
234 bspline.template eval<iganet::deriv::dx, memory_optimized>(xi)(0)));
235 EXPECT_TRUE(torch::allclose(
236 bspline_igrad_val(0, 1),
237 bspline.template eval<iganet::deriv::dy, memory_optimized>(xi)(0)));
238 } else if constexpr (Spline::parDim() == 3) {
239 EXPECT_TRUE(torch::allclose(
240 bspline_igrad_val(0, 0),
241 bspline.template eval<iganet::deriv::dx, memory_optimized>(xi)(0)));
242 EXPECT_TRUE(torch::allclose(
243 bspline_igrad_val(0, 1),
244 bspline.template eval<iganet::deriv::dy, memory_optimized>(xi)(0)));
245 EXPECT_TRUE(torch::allclose(
246 bspline_igrad_val(0, 2),
247 bspline.template eval<iganet::deriv::dz, memory_optimized>(xi)(0)));
248 } else if constexpr (Spline::parDim() == 4) {
249 EXPECT_TRUE(torch::allclose(
250 bspline_igrad_val(0, 0),
251 bspline.template eval<iganet::deriv::dx, memory_optimized>(xi)(0)));
252 EXPECT_TRUE(torch::allclose(
253 bspline_igrad_val(0, 1),
254 bspline.template eval<iganet::deriv::dy, memory_optimized>(xi)(0)));
255 EXPECT_TRUE(torch::allclose(
256 bspline_igrad_val(0, 2),
257 bspline.template eval<iganet::deriv::dz, memory_optimized>(xi)(0)));
258 EXPECT_TRUE(torch::allclose(
259 bspline_igrad_val(0, 3),
260 bspline.template eval<iganet::deriv::dt, memory_optimized>(xi)(0)));
261 }
262}
263
264template <bool memory_optimized, bool precompute, typename Spline,
265 typename TensorArray_t>
266void test_bspline_jac(const Spline &bspline, const TensorArray_t &xi,
267 typename Spline::value_type tol = 1e-12) {
268 iganet::utils::BlockTensor<torch::Tensor, Spline::geoDim(), Spline::parDim()>
269 bspline_jac_val;
270
271 if constexpr (precompute) {
272 auto knot_indices = bspline.find_knot_indices(xi);
273 auto coeff_indices = bspline.find_coeff_indices(knot_indices);
274 bspline_jac_val = bspline.jac(xi, knot_indices, coeff_indices);
275 } else
276 bspline_jac_val = bspline.jac(xi);
277
278 if constexpr (Spline::parDim() >= 1) {
279 for (iganet::short_t k = 0; k < Spline::geoDim(); ++k)
280 EXPECT_TRUE(torch::allclose(
281 bspline_jac_val(k, 0),
282 bspline.template eval<iganet::deriv::dx, memory_optimized>(xi)(k)));
283 }
284
285 if constexpr (Spline::parDim() >= 2) {
286 for (iganet::short_t k = 0; k < Spline::geoDim(); ++k)
287 EXPECT_TRUE(torch::allclose(
288 bspline_jac_val(k, 1),
289 bspline.template eval<iganet::deriv::dy, memory_optimized>(xi)(k)));
290 }
291
292 if constexpr (Spline::parDim() >= 3) {
293 for (iganet::short_t k = 0; k < Spline::geoDim(); ++k)
294 EXPECT_TRUE(torch::allclose(
295 bspline_jac_val(k, 2),
296 bspline.template eval<iganet::deriv::dz, memory_optimized>(xi)(k)));
297 }
298
299 if constexpr (Spline::parDim() >= 4) {
300 for (iganet::short_t k = 0; k < Spline::geoDim(); ++k)
301 EXPECT_TRUE(torch::allclose(
302 bspline_jac_val(k, 3),
303 bspline.template eval<iganet::deriv::dt, memory_optimized>(xi)(k)));
304 }
305}
306
307template <bool memory_optimized, bool precompute, typename Geometry_t,
308 typename Spline, typename TensorArray_t>
309void test_bspline_ijac(const Geometry_t &geometry, const Spline &bspline,
310 const TensorArray_t &xi,
311 typename Spline::value_type tol = 1e-12) {
312 iganet::utils::BlockTensor<torch::Tensor, Spline::geoDim(), Spline::parDim()>
313 bspline_ijac_val;
314
315 if constexpr (precompute) {
316 auto knot_indices = bspline.find_knot_indices(xi);
317 auto coeff_indices = bspline.find_coeff_indices(knot_indices);
318 auto knot_indices_G = geometry.find_knot_indices(xi);
319 auto coeff_indices_G = geometry.find_coeff_indices(knot_indices_G);
320 bspline_ijac_val = bspline.ijac(geometry, xi, knot_indices, coeff_indices,
321 knot_indices_G, coeff_indices_G);
322 } else
323 bspline_ijac_val = bspline.ijac(geometry, xi);
324
325 if constexpr (Spline::parDim() >= 1) {
326 for (iganet::short_t k = 0; k < Spline::geoDim(); ++k)
327 EXPECT_TRUE(torch::allclose(
328 bspline_ijac_val(k, 0),
329 bspline.template eval<iganet::deriv::dx, memory_optimized>(xi)(k)));
330 }
331
332 if constexpr (Spline::parDim() >= 2) {
333 for (iganet::short_t k = 0; k < Spline::geoDim(); ++k)
334 EXPECT_TRUE(torch::allclose(
335 bspline_ijac_val(k, 1),
336 bspline.template eval<iganet::deriv::dy, memory_optimized>(xi)(k)));
337 }
338
339 if constexpr (Spline::parDim() >= 3) {
340 for (iganet::short_t k = 0; k < Spline::geoDim(); ++k)
341 EXPECT_TRUE(torch::allclose(
342 bspline_ijac_val(k, 2),
343 bspline.template eval<iganet::deriv::dz, memory_optimized>(xi)(k)));
344 }
345
346 if constexpr (Spline::parDim() >= 4) {
347 for (iganet::short_t k = 0; k < Spline::geoDim(); ++k)
348 EXPECT_TRUE(torch::allclose(
349 bspline_ijac_val(k, 3),
350 bspline.template eval<iganet::deriv::dt, memory_optimized>(xi)(k)));
351 }
352}
353
354template <bool memory_optimized, bool precompute, typename Spline,
355 typename TensorArray_t>
356void test_bspline_hess(const Spline &bspline, const TensorArray_t &xi,
357 typename Spline::value_type tol = 1e-12) {
358 iganet::utils::BlockTensor<torch::Tensor, Spline::parDim(), Spline::parDim(),
359 Spline::geoDim()>
360 bspline_hess_val;
361
362 if constexpr (precompute) {
363 auto knot_indices = bspline.find_knot_indices(xi);
364 auto coeff_indices = bspline.find_coeff_indices(knot_indices);
365 bspline_hess_val = bspline.hess(xi, knot_indices, coeff_indices);
366 } else
367 bspline_hess_val = bspline.hess(xi);
368
369 if constexpr (Spline::parDim() >= 1) {
370 for (iganet::short_t k = 0; k < Spline::geoDim(); ++k)
371 EXPECT_TRUE(torch::allclose(
372 bspline_hess_val(0, 0, k),
373 bspline.template eval<iganet::deriv::dx ^ 2, memory_optimized>(xi)(
374 k)));
375 }
376
377 if constexpr (Spline::parDim() >= 2) {
378 for (iganet::short_t k = 0; k < Spline::geoDim(); ++k) {
379 EXPECT_TRUE(torch::allclose(
380 bspline_hess_val(0, 1, k),
381 bspline.template eval<iganet::deriv::dx + iganet::deriv::dy,
382 memory_optimized>(xi)(k)));
383 EXPECT_TRUE(torch::allclose(
384 bspline_hess_val(1, 0, k),
385 bspline.template eval<iganet::deriv::dy + iganet::deriv::dx,
386 memory_optimized>(xi)(k)));
387 EXPECT_TRUE(torch::allclose(
388 bspline_hess_val(1, 1, k),
389 bspline.template eval<iganet::deriv::dy ^ 2, memory_optimized>(xi)(
390 k)));
391 }
392 }
393
394 if constexpr (Spline::parDim() >= 3) {
395 for (iganet::short_t k = 0; k < Spline::geoDim(); ++k) {
396 EXPECT_TRUE(torch::allclose(
397 bspline_hess_val(0, 2, k),
398 bspline.template eval<iganet::deriv::dx + iganet::deriv::dz,
399 memory_optimized>(xi)(k)));
400 EXPECT_TRUE(torch::allclose(
401 bspline_hess_val(1, 2, k),
402 bspline.template eval<iganet::deriv::dy + iganet::deriv::dz,
403 memory_optimized>(xi)(k)));
404 EXPECT_TRUE(torch::allclose(
405 bspline_hess_val(2, 0, k),
406 bspline.template eval<iganet::deriv::dz + iganet::deriv::dx,
407 memory_optimized>(xi)(k)));
408 EXPECT_TRUE(torch::allclose(
409 bspline_hess_val(2, 1, k),
410 bspline.template eval<iganet::deriv::dz + iganet::deriv::dy,
411 memory_optimized>(xi)(k)));
412 EXPECT_TRUE(torch::allclose(
413 bspline_hess_val(2, 2, k),
414 bspline.template eval<iganet::deriv::dz ^ 2, memory_optimized>(xi)(
415 k)));
416 }
417 }
418
419 if constexpr (Spline::parDim() >= 4) {
420 for (iganet::short_t k = 0; k < Spline::geoDim(); ++k) {
421 EXPECT_TRUE(torch::allclose(
422 bspline_hess_val(0, 3, k),
423 bspline.template eval<iganet::deriv::dx + iganet::deriv::dt,
424 memory_optimized>(xi)(k)));
425 EXPECT_TRUE(torch::allclose(
426 bspline_hess_val(1, 3, k),
427 bspline.template eval<iganet::deriv::dy + iganet::deriv::dt,
428 memory_optimized>(xi)(k)));
429 EXPECT_TRUE(torch::allclose(
430 bspline_hess_val(2, 3, k),
431 bspline.template eval<iganet::deriv::dz + iganet::deriv::dt,
432 memory_optimized>(xi)(k)));
433 EXPECT_TRUE(torch::allclose(
434 bspline_hess_val(3, 0, k),
435 bspline.template eval<iganet::deriv::dt + iganet::deriv::dx,
436 memory_optimized>(xi)(k)));
437 EXPECT_TRUE(torch::allclose(
438 bspline_hess_val(3, 1, k),
439 bspline.template eval<iganet::deriv::dt + iganet::deriv::dy,
440 memory_optimized>(xi)(k)));
441 EXPECT_TRUE(torch::allclose(
442 bspline_hess_val(3, 2, k),
443 bspline.template eval<iganet::deriv::dt + iganet::deriv::dz,
444 memory_optimized>(xi)(k)));
445 EXPECT_TRUE(torch::allclose(
446 bspline_hess_val(3, 3, k),
447 bspline.template eval<iganet::deriv::dt ^ 2, memory_optimized>(xi)(
448 k)));
449 }
450 }
451}
452
453template <bool memory_optimized, bool precompute, typename Geometry_t,
454 typename Spline, typename TensorArray_t>
455void test_bspline_ihess(const Geometry_t &geometry, const Spline &bspline,
456 const TensorArray_t &xi,
457 typename Spline::value_type tol = 1e-12) {
458 iganet::utils::BlockTensor<torch::Tensor, Spline::parDim(), Spline::parDim(),
459 Spline::geoDim()>
460 bspline_ihess_val;
461
462 if constexpr (precompute) {
463 auto knot_indices = bspline.find_knot_indices(xi);
464 auto coeff_indices = bspline.find_coeff_indices(knot_indices);
465 auto knot_indices_G = geometry.find_knot_indices(xi);
466 auto coeff_indices_G = geometry.find_coeff_indices(knot_indices_G);
467 bspline_ihess_val = bspline.ihess(geometry, xi, knot_indices, coeff_indices,
468 knot_indices_G, coeff_indices_G);
469 } else
470 bspline_ihess_val = bspline.ihess(geometry, xi);
471
472 if constexpr (Spline::parDim() >= 1) {
473 for (iganet::short_t k = 0; k < Spline::geoDim(); ++k)
474 EXPECT_TRUE(torch::allclose(
475 bspline_ihess_val(0, 0, k),
476 bspline.template eval<iganet::deriv::dx ^ 2, memory_optimized>(xi)(
477 k)));
478 }
479
480 if constexpr (Spline::parDim() >= 2) {
481 for (iganet::short_t k = 0; k < Spline::geoDim(); ++k) {
482 EXPECT_TRUE(torch::allclose(
483 bspline_ihess_val(0, 1, k),
484 bspline.template eval<iganet::deriv::dx + iganet::deriv::dy,
485 memory_optimized>(xi)(k)));
486 EXPECT_TRUE(torch::allclose(
487 bspline_ihess_val(1, 0, k),
488 bspline.template eval<iganet::deriv::dy + iganet::deriv::dx,
489 memory_optimized>(xi)(k)));
490 EXPECT_TRUE(torch::allclose(
491 bspline_ihess_val(1, 1, k),
492 bspline.template eval<iganet::deriv::dy ^ 2, memory_optimized>(xi)(
493 k)));
494 }
495 }
496
497 if constexpr (Spline::parDim() >= 3) {
498 for (iganet::short_t k = 0; k < Spline::geoDim(); ++k) {
499 EXPECT_TRUE(torch::allclose(
500 bspline_ihess_val(0, 2, k),
501 bspline.template eval<iganet::deriv::dx + iganet::deriv::dz,
502 memory_optimized>(xi)(k)));
503 EXPECT_TRUE(torch::allclose(
504 bspline_ihess_val(1, 2, k),
505 bspline.template eval<iganet::deriv::dy + iganet::deriv::dz,
506 memory_optimized>(xi)(k)));
507 EXPECT_TRUE(torch::allclose(
508 bspline_ihess_val(2, 0, k),
509 bspline.template eval<iganet::deriv::dz + iganet::deriv::dx,
510 memory_optimized>(xi)(k)));
511 EXPECT_TRUE(torch::allclose(
512 bspline_ihess_val(2, 1, k),
513 bspline.template eval<iganet::deriv::dz + iganet::deriv::dy,
514 memory_optimized>(xi)(k)));
515 EXPECT_TRUE(torch::allclose(
516 bspline_ihess_val(2, 2, k),
517 bspline.template eval<iganet::deriv::dz ^ 2, memory_optimized>(xi)(
518 k)));
519 }
520 }
521
522 if constexpr (Spline::parDim() >= 4) {
523 for (iganet::short_t k = 0; k < Spline::geoDim(); ++k) {
524 EXPECT_TRUE(torch::allclose(
525 bspline_ihess_val(0, 3, k),
526 bspline.template eval<iganet::deriv::dx + iganet::deriv::dt,
527 memory_optimized>(xi)(k)));
528 EXPECT_TRUE(torch::allclose(
529 bspline_ihess_val(1, 3, k),
530 bspline.template eval<iganet::deriv::dy + iganet::deriv::dt,
531 memory_optimized>(xi)(k)));
532 EXPECT_TRUE(torch::allclose(
533 bspline_ihess_val(2, 3, k),
534 bspline.template eval<iganet::deriv::dz + iganet::deriv::dt,
535 memory_optimized>(xi)(k)));
536 EXPECT_TRUE(torch::allclose(
537 bspline_ihess_val(3, 0, k),
538 bspline.template eval<iganet::deriv::dt + iganet::deriv::dx,
539 memory_optimized>(xi)(k)));
540 EXPECT_TRUE(torch::allclose(
541 bspline_ihess_val(3, 1, k),
542 bspline.template eval<iganet::deriv::dt + iganet::deriv::dy,
543 memory_optimized>(xi)(k)));
544 EXPECT_TRUE(torch::allclose(
545 bspline_ihess_val(3, 2, k),
546 bspline.template eval<iganet::deriv::dt + iganet::deriv::dz,
547 memory_optimized>(xi)(k)));
548 EXPECT_TRUE(torch::allclose(
549 bspline_ihess_val(3, 3, k),
550 bspline.template eval<iganet::deriv::dt ^ 2, memory_optimized>(xi)(
551 k)));
552 }
553 }
554}
555
556template <typename Geometry_t, typename Spline, typename TensorArray_t>
557void test_bspline_eval(const Geometry_t &geometry, const Spline &bspline,
558 const TensorArray_t &xi,
559 typename Spline::value_type tol = 1e-12) {
560 // Create B-Spline
561 auto bsplinelib_bspline = to_bsplinelib_bspline(bspline);
562
563 // Evaluate function and derivatives (non-memory optimized)
564 test_bspline_eval<iganet::deriv::func, false, false>(
565 bspline, bsplinelib_bspline, xi, tol);
566
567 if constexpr (Spline::parDim() == 1) {
568 test_bspline_eval<iganet::deriv::dx, false, false>(
569 bspline, bsplinelib_bspline, xi, tol);
570 test_bspline_eval<iganet::deriv::dx ^ 2, false, false>(
571 bspline, bsplinelib_bspline, xi, tol);
572 test_bspline_eval<iganet::deriv::dx ^ 3, false, false>(
573 bspline, bsplinelib_bspline, xi, tol);
574 test_bspline_eval<iganet::deriv::dx ^ 4, false, false>(
575 bspline, bsplinelib_bspline, xi, tol);
576 }
577
578 if constexpr (Spline::parDim() == 2) {
579 test_bspline_eval<iganet::deriv::dy, false, false>(
580 bspline, bsplinelib_bspline, xi, tol);
581 test_bspline_eval<iganet::deriv::dy ^ 2, false, false>(
582 bspline, bsplinelib_bspline, xi, tol);
583 test_bspline_eval<iganet::deriv::dy ^ 3, false, false>(
584 bspline, bsplinelib_bspline, xi, tol);
585 test_bspline_eval<iganet::deriv::dy ^ 4, false, false>(
586 bspline, bsplinelib_bspline, xi, tol);
587 }
588
589 if constexpr (Spline::parDim() == 3) {
590 test_bspline_eval<iganet::deriv::dz, false, false>(
591 bspline, bsplinelib_bspline, xi, tol);
592 test_bspline_eval<iganet::deriv::dz ^ 2, false, false>(
593 bspline, bsplinelib_bspline, xi, tol);
594 test_bspline_eval<iganet::deriv::dz ^ 3, false, false>(
595 bspline, bsplinelib_bspline, xi, tol);
596 test_bspline_eval<iganet::deriv::dz ^ 4, false, false>(
597 bspline, bsplinelib_bspline, xi, tol);
598 }
599
600 if constexpr (Spline::parDim() == 4) {
601 test_bspline_eval<iganet::deriv::dt, false, false>(
602 bspline, bsplinelib_bspline, xi, tol);
603 test_bspline_eval<iganet::deriv::dt ^ 2, false, false>(
604 bspline, bsplinelib_bspline, xi, tol);
605 test_bspline_eval<iganet::deriv::dt ^ 3, false, false>(
606 bspline, bsplinelib_bspline, xi, tol);
607 test_bspline_eval<iganet::deriv::dt ^ 4, false, false>(
608 bspline, bsplinelib_bspline, xi, tol);
609 }
610
611 // Evaluate function and derivatives (memory optimized)
612 test_bspline_eval<iganet::deriv::func, true, false>(
613 bspline, bsplinelib_bspline, xi, tol);
614
615 if constexpr (Spline::parDim() == 1) {
616 test_bspline_eval<iganet::deriv::dx, true, false>(
617 bspline, bsplinelib_bspline, xi, tol);
618 test_bspline_eval<iganet::deriv::dx ^ 2, true, false>(
619 bspline, bsplinelib_bspline, xi, tol);
620 test_bspline_eval<iganet::deriv::dx ^ 3, true, false>(
621 bspline, bsplinelib_bspline, xi, tol);
622 test_bspline_eval<iganet::deriv::dx ^ 4, true, false>(
623 bspline, bsplinelib_bspline, xi, tol);
624 }
625
626 if constexpr (Spline::parDim() == 2) {
627 test_bspline_eval<iganet::deriv::dy, true, false>(
628 bspline, bsplinelib_bspline, xi, tol);
629 test_bspline_eval<iganet::deriv::dy ^ 2, true, false>(
630 bspline, bsplinelib_bspline, xi, tol);
631 test_bspline_eval<iganet::deriv::dy ^ 3, true, false>(
632 bspline, bsplinelib_bspline, xi, tol);
633 test_bspline_eval<iganet::deriv::dy ^ 4, true, false>(
634 bspline, bsplinelib_bspline, xi, tol);
635 }
636
637 if constexpr (Spline::parDim() == 3) {
638 test_bspline_eval<iganet::deriv::dz, true, false>(
639 bspline, bsplinelib_bspline, xi, tol);
640 test_bspline_eval<iganet::deriv::dz ^ 2, true, false>(
641 bspline, bsplinelib_bspline, xi, tol);
642 test_bspline_eval<iganet::deriv::dz ^ 3, true, false>(
643 bspline, bsplinelib_bspline, xi, tol);
644 test_bspline_eval<iganet::deriv::dz ^ 4, true, false>(
645 bspline, bsplinelib_bspline, xi, tol);
646 }
647
648 if constexpr (Spline::parDim() == 4) {
649 test_bspline_eval<iganet::deriv::dt, true, false>(
650 bspline, bsplinelib_bspline, xi, tol);
651 test_bspline_eval<iganet::deriv::dt ^ 2, true, false>(
652 bspline, bsplinelib_bspline, xi, tol);
653 test_bspline_eval<iganet::deriv::dt ^ 3, true, false>(
654 bspline, bsplinelib_bspline, xi, tol);
655 test_bspline_eval<iganet::deriv::dt ^ 4, true, false>(
656 bspline, bsplinelib_bspline, xi, tol);
657 }
658
659 // Evaluate function and derivatives from precomputed data (non-memory
660 // optimized)
661 test_bspline_eval<iganet::deriv::func, false, true>(
662 bspline, bsplinelib_bspline, xi, tol);
663
664 if constexpr (Spline::parDim() == 1) {
665 test_bspline_eval<iganet::deriv::dx, false, true>(
666 bspline, bsplinelib_bspline, xi, tol);
667 test_bspline_eval<iganet::deriv::dx ^ 2, false, true>(
668 bspline, bsplinelib_bspline, xi, tol);
669 test_bspline_eval<iganet::deriv::dx ^ 3, false, true>(
670 bspline, bsplinelib_bspline, xi, tol);
671 test_bspline_eval<iganet::deriv::dx ^ 4, false, true>(
672 bspline, bsplinelib_bspline, xi, tol);
673 }
674
675 if constexpr (Spline::parDim() == 2) {
676 test_bspline_eval<iganet::deriv::dy, false, true>(
677 bspline, bsplinelib_bspline, xi, tol);
678 test_bspline_eval<iganet::deriv::dy ^ 2, false, true>(
679 bspline, bsplinelib_bspline, xi, tol);
680 test_bspline_eval<iganet::deriv::dy ^ 3, false, true>(
681 bspline, bsplinelib_bspline, xi, tol);
682 test_bspline_eval<iganet::deriv::dy ^ 4, false, true>(
683 bspline, bsplinelib_bspline, xi, tol);
684 }
685
686 if constexpr (Spline::parDim() == 3) {
687 test_bspline_eval<iganet::deriv::dz, false, true>(
688 bspline, bsplinelib_bspline, xi, tol);
689 test_bspline_eval<iganet::deriv::dz ^ 2, false, true>(
690 bspline, bsplinelib_bspline, xi, tol);
691 test_bspline_eval<iganet::deriv::dz ^ 3, false, true>(
692 bspline, bsplinelib_bspline, xi, tol);
693 test_bspline_eval<iganet::deriv::dz ^ 4, false, true>(
694 bspline, bsplinelib_bspline, xi, tol);
695 }
696
697 if constexpr (Spline::parDim() == 4) {
698 test_bspline_eval<iganet::deriv::dt, false, true>(
699 bspline, bsplinelib_bspline, xi, tol);
700 test_bspline_eval<iganet::deriv::dt ^ 2, false, true>(
701 bspline, bsplinelib_bspline, xi, tol);
702 test_bspline_eval<iganet::deriv::dt ^ 3, false, true>(
703 bspline, bsplinelib_bspline, xi, tol);
704 test_bspline_eval<iganet::deriv::dt ^ 4, false, true>(
705 bspline, bsplinelib_bspline, xi, tol);
706 }
707
708 // Evaluate function and derivatives from precomputed data (memory optimized)
709 test_bspline_eval<iganet::deriv::func, true, true>(
710 bspline, bsplinelib_bspline, xi, tol);
711
712 if constexpr (Spline::parDim() == 1) {
713 test_bspline_eval<iganet::deriv::dx, true, true>(
714 bspline, bsplinelib_bspline, xi, tol);
715 test_bspline_eval<iganet::deriv::dx ^ 2, true, true>(
716 bspline, bsplinelib_bspline, xi, tol);
717 test_bspline_eval<iganet::deriv::dx ^ 3, true, true>(
718 bspline, bsplinelib_bspline, xi, tol);
719 test_bspline_eval<iganet::deriv::dx ^ 4, true, true>(
720 bspline, bsplinelib_bspline, xi, tol);
721 }
722
723 if constexpr (Spline::parDim() == 2) {
724 test_bspline_eval<iganet::deriv::dy, true, true>(
725 bspline, bsplinelib_bspline, xi, tol);
726 test_bspline_eval<iganet::deriv::dy ^ 2, true, true>(
727 bspline, bsplinelib_bspline, xi, tol);
728 test_bspline_eval<iganet::deriv::dy ^ 3, true, true>(
729 bspline, bsplinelib_bspline, xi, tol);
730 test_bspline_eval<iganet::deriv::dy ^ 4, true, true>(
731 bspline, bsplinelib_bspline, xi, tol);
732 }
733
734 if constexpr (Spline::parDim() == 3) {
735 test_bspline_eval<iganet::deriv::dz, true, true>(
736 bspline, bsplinelib_bspline, xi, tol);
737 test_bspline_eval<iganet::deriv::dz ^ 2, true, true>(
738 bspline, bsplinelib_bspline, xi, tol);
739 test_bspline_eval<iganet::deriv::dz ^ 3, true, true>(
740 bspline, bsplinelib_bspline, xi, tol);
741 test_bspline_eval<iganet::deriv::dz ^ 4, true, true>(
742 bspline, bsplinelib_bspline, xi, tol);
743 }
744
745 if constexpr (Spline::parDim() == 4) {
746 test_bspline_eval<iganet::deriv::dt, true, true>(
747 bspline, bsplinelib_bspline, xi, tol);
748 test_bspline_eval<iganet::deriv::dt ^ 2, true, true>(
749 bspline, bsplinelib_bspline, xi, tol);
750 test_bspline_eval<iganet::deriv::dt ^ 3, true, true>(
751 bspline, bsplinelib_bspline, xi, tol);
752 test_bspline_eval<iganet::deriv::dt ^ 4, true, true>(
753 bspline, bsplinelib_bspline, xi, tol);
754 }
755
756 // Evaluate gradients
757 if constexpr (Spline::geoDim() == 1) {
758 test_bspline_grad<false, false>(bspline, xi, tol);
759 test_bspline_grad<false, true>(bspline, xi, tol);
760 test_bspline_grad<true, false>(bspline, xi, tol);
761 test_bspline_grad<true, true>(bspline, xi, tol);
762
763 test_bspline_igrad<false, false>(geometry, bspline, xi, tol);
764 test_bspline_igrad<false, true>(geometry, bspline, xi, tol);
765 test_bspline_igrad<true, false>(geometry, bspline, xi, tol);
766 test_bspline_igrad<true, true>(geometry, bspline, xi, tol);
767 }
768
770 test_bspline_jac<false, false>(bspline, xi, tol);
771 test_bspline_jac<false, true>(bspline, xi, tol);
772 test_bspline_jac<true, false>(bspline, xi, tol);
773 test_bspline_jac<true, true>(bspline, xi, tol);
774
775 test_bspline_ijac<false, false>(geometry, bspline, xi, tol);
776 test_bspline_ijac<false, true>(geometry, bspline, xi, tol);
777 test_bspline_ijac<true, false>(geometry, bspline, xi, tol);
778 test_bspline_ijac<true, true>(geometry, bspline, xi, tol);
779
781 if constexpr (Spline::geoDim() == 1) {
782 test_bspline_hess<false, false>(bspline, xi, tol);
783 test_bspline_hess<false, true>(bspline, xi, tol);
784 test_bspline_hess<true, false>(bspline, xi, tol);
785 test_bspline_hess<true, true>(bspline, xi, tol);
786
787 test_bspline_ihess<false, false>(geometry, bspline, xi, tol);
788 test_bspline_ihess<false, true>(geometry, bspline, xi, tol);
789 test_bspline_ihess<true, false>(geometry, bspline, xi, tol);
790 test_bspline_ihess<true, true>(geometry, bspline, xi, tol);
791 }
792}
iganet.h
Isogeometric analysis network main header file.
iganet::utils::BlockTensor
Forward declaration of BlockTensor.
Definition blocktensor.hpp:46
iganet::deriv
deriv
Enumerator for specifying the derivative of B-spline evaluation.
Definition bspline.hpp:74
iganet::short_t
short int short_t
Definition core.hpp:74
test_bspline_jac
void test_bspline_jac(const Spline &bspline, const TensorArray_t &xi, typename Spline::value_type tol=1e-12)
Definition unittest_bsplinelib.hpp:266
test_bspline_grad
void test_bspline_grad(const Spline &bspline, const TensorArray_t &xi, typename Spline::value_type tol=1e-12)
Definition unittest_bsplinelib.hpp:160
test_bspline_igrad
void test_bspline_igrad(const Geometry_t &geometry, const Spline &bspline, const TensorArray_t &xi, typename Spline::value_type tol=1e-12)
Definition unittest_bsplinelib.hpp:211
test_bspline_ihess
void test_bspline_ihess(const Geometry_t &geometry, const Spline &bspline, const TensorArray_t &xi, typename Spline::value_type tol=1e-12)
Definition unittest_bsplinelib.hpp:455
test_bspline_hess
void test_bspline_hess(const Spline &bspline, const TensorArray_t &xi, typename Spline::value_type tol=1e-12)
Definition unittest_bsplinelib.hpp:356
test_bspline_ijac
void test_bspline_ijac(const Geometry_t &geometry, const Spline &bspline, const TensorArray_t &xi, typename Spline::value_type tol=1e-12)
Definition unittest_bsplinelib.hpp:309
to_bsplinelib_bspline
auto to_bsplinelib_bspline(const Spline &bspline)
Definition unittest_bsplinelib.hpp:33
test_bspline_eval
void test_bspline_eval(const Spline &bspline, BSplineLibSpline bsplinelib_bspline, const TensorArray_t &xi, typename Spline::value_type tol=1e-12)
Definition unittest_bsplinelib.hpp:96