/home/runner/work/kynema/kynema/kynema/src/constraints/calculate_constraint_residual_gradient.hpp Source File

Kynema API: /home/runner/work/kynema/kynema/kynema/src/constraints/calculate_constraint_residual_gradient.hpp Source File
Kynema API
A flexible multibody structural dynamics code for wind turbines
Loading...
Searching...
No Matches
calculate_constraint_residual_gradient.hpp
Go to the documentation of this file.
1#pragma once
2
3#include <KokkosBatched_Copy_Decl.hpp>
4#include <KokkosBatched_Gemm_Decl.hpp>
5#include <KokkosBlas.hpp>
6#include <Kokkos_Core.hpp>
7
8#include "calculate_fixed_bc_3DOF_constraint.hpp"
9#include "calculate_fixed_bc_constraint.hpp"
10#include "calculate_prescribed_bc_3DOF_constraint.hpp"
11#include "calculate_prescribed_bc_constraint.hpp"
12#include "calculate_revolute_joint_constraint.hpp"
13#include "calculate_revolute_joint_force.hpp"
14#include "calculate_rigid_joint_3DOF_constraint.hpp"
15#include "calculate_rigid_joint_constraint.hpp"
16#include "calculate_rotation_control_constraint.hpp"
17#include "constraint_type.hpp"
18
19namespace kynema::constraints {
20
29template <typename DeviceType>
30struct CalculateConstraintResidualGradient {
31 using TransposeMatrix = KokkosBatched::SerialCopy<KokkosBatched::Trans::Transpose>;
32 using Gemm = KokkosBatched::SerialGemm<
33 KokkosBatched::Trans::NoTranspose, KokkosBatched::Trans::NoTranspose,
34 KokkosBatched::Algo::Gemm::Default>;
35 using CopyVector = KokkosBatched::SerialCopy<KokkosBatched::Trans::NoTranspose, 1>;
36 using CopyMatrix = KokkosBatched::SerialCopy<KokkosBatched::Trans::NoTranspose, 2>;
37 template <typename ValueType>
38 using View = Kokkos::View<ValueType, DeviceType>;
39 template <typename ValueType>
40 using ConstView = typename View<ValueType>::const_type;
41
42 ConstView<ConstraintType*> type_;
43 ConstView<size_t*> base_node_index_;
44 ConstView<size_t*> target_node_index_;
45 ConstView<double* [3]> X0_;
46 ConstView<double* [3][3]> axes_;
47 ConstView<double* [7]> constraint_inputs_;
48 ConstView<double* [6]> lambda_;
49 ConstView<double* [6][6]> tangent_;
50 ConstView<double* [7]> node_u_;
51 View<double* [6]> res_;
52 View<double* [6]> b_lambda_res_;
53 View<double* [6]> t_lambda_res_;
54 View<double* [6]> system_res_;
55 View<double* [6][6]> b_grad_;
56 View<double* [6][6]> t_grad_;
57 View<double* [6][6]> b_grad_trans_;
58 View<double* [6][6]> t_grad_trans_;
59
60 KOKKOS_FUNCTION
61 void FixedBC(size_t constraint) const {
62 using Kokkos::ALL;
63 using Kokkos::Array;
64 using Kokkos::subview;
65 using KokkosBlas::Experimental::serial_gemv;
66
67 auto X0_data = Array<double, 3>{};
68 auto t_node_u_data = Array<double, 7>{};
69 auto lambda_data = Array<double, 6>{};
70 auto res_data = Array<double, 6>{};
71 auto t_lambda_res_data = Array<double, 6>{};
72 auto t_grad_data = Array<double, 36>{};
73 auto t_grad_trans_data = Array<double, 36>{};
74 auto target_tangent_data = Array<double, 36>{};
75 auto t_grad_tan_data = Array<double, 36>{};
76
77 const auto X0 = View<double[3]>(X0_data.data());
78 const auto t_node_u = View<double[7]>(t_node_u_data.data());
79 const auto lambda = View<double[6]>(lambda_data.data());
80 const auto res = View<double[6]>(res_data.data());
81 const auto t_lambda_res = View<double[6]>(t_lambda_res_data.data());
82 const auto t_grad = View<double[6][6]>(t_grad_data.data());
83 const auto t_grad_trans = View<double[6][6]>(t_grad_trans_data.data());
84 const auto target_tangent = View<double[6][6]>(target_tangent_data.data());
85 const auto t_grad_tan = View<double[6][6]>(t_grad_tan_data.data());
86
87 const auto target_node_index = target_node_index_(constraint);
88
89 CopyVector::invoke(subview(X0_, constraint, ALL), X0);
90 CopyVector::invoke(subview(node_u_, target_node_index, ALL), t_node_u);
91 CopyVector::invoke(subview(lambda_, constraint, ALL), lambda);
92 CopyMatrix::invoke(subview(tangent_, target_node_index, ALL, ALL), target_tangent);
93
94 CalculateFixedBCConstraint<DeviceType>::invoke(X0, t_node_u, res, t_grad);
95
96 TransposeMatrix::invoke(t_grad, t_grad_trans);
97 Gemm::invoke(1., t_grad, target_tangent, 0., t_grad_tan);
98 serial_gemv('N', 1., t_grad_trans, lambda, 0., t_lambda_res);
99
100 CopyVector::invoke(res, subview(res_, constraint, ALL));
101 CopyVector::invoke(t_lambda_res, subview(t_lambda_res_, constraint, ALL));
102 CopyMatrix::invoke(t_grad_tan, subview(t_grad_, constraint, ALL, ALL));
103 CopyMatrix::invoke(t_grad_trans, subview(t_grad_trans_, constraint, ALL, ALL));
104 }
105
106 KOKKOS_FUNCTION
107 void FixedBC3DOF(size_t constraint) const {
108 using Kokkos::ALL;
109 using Kokkos::Array;
110 using Kokkos::subview;
111 using KokkosBlas::Experimental::serial_gemv;
112
113 auto X0_data = Array<double, 3>{};
114 auto t_node_u_data = Array<double, 7>{};
115 auto lambda_data = Array<double, 6>{};
116 auto res_data = Array<double, 6>{};
117 auto t_lambda_res_data = Array<double, 6>{};
118 auto t_grad_data = Array<double, 36>{};
119 auto t_grad_trans_data = Array<double, 36>{};
120 auto target_tangent_data = Array<double, 36>{};
121 auto t_grad_tan_data = Array<double, 36>{};
122
123 const auto X0 = View<double[3]>(X0_data.data());
124 const auto t_node_u = View<double[7]>(t_node_u_data.data());
125 const auto lambda = View<double[6]>(lambda_data.data());
126 const auto res = View<double[6]>(res_data.data());
127 const auto t_lambda_res = View<double[6]>(t_lambda_res_data.data());
128 const auto t_grad = View<double[6][6]>(t_grad_data.data());
129 const auto t_grad_trans = View<double[6][6]>(t_grad_trans_data.data());
130 const auto target_tangent = View<double[6][6]>(target_tangent_data.data());
131 const auto t_grad_tan = View<double[6][6]>(t_grad_tan_data.data());
132
133 const auto target_node_index = target_node_index_(constraint);
134
135 CopyVector::invoke(subview(X0_, constraint, ALL), X0);
136 CopyVector::invoke(subview(node_u_, target_node_index, ALL), t_node_u);
137 CopyVector::invoke(subview(lambda_, constraint, ALL), lambda);
138 CopyMatrix::invoke(subview(tangent_, target_node_index, ALL, ALL), target_tangent);
139
140 CalculateFixedBC3DOFConstraint<DeviceType>::invoke(X0, t_node_u, res, t_grad);
141
142 TransposeMatrix::invoke(t_grad, t_grad_trans);
143 Gemm::invoke(1., t_grad, target_tangent, 0., t_grad_tan);
144 serial_gemv('N', 1., t_grad_trans, lambda, 0., t_lambda_res);
145
146 CopyVector::invoke(res, subview(res_, constraint, ALL));
147 CopyVector::invoke(t_lambda_res, subview(t_lambda_res_, constraint, ALL));
148 CopyMatrix::invoke(t_grad_tan, subview(t_grad_, constraint, ALL, ALL));
149 CopyMatrix::invoke(t_grad_trans, subview(t_grad_trans_, constraint, ALL, ALL));
150 }
151
152 KOKKOS_FUNCTION
153 void PrescribedBC(size_t constraint) const {
154 using Kokkos::ALL;
155 using Kokkos::Array;
156 using Kokkos::subview;
157 using KokkosBlas::Experimental::serial_gemv;
158
159 auto X0_data = Array<double, 3>{};
160 auto inputs_data = Array<double, 7>{};
161 auto t_node_u_data = Array<double, 7>{};
162 auto lambda_data = Array<double, 6>{};
163 auto res_data = Array<double, 6>{};
164 auto t_lambda_res_data = Array<double, 6>{};
165 auto t_grad_data = Array<double, 36>{};
166 auto t_grad_trans_data = Array<double, 36>{};
167 auto target_tangent_data = Array<double, 36>{};
168 auto t_grad_tan_data = Array<double, 36>{};
169
170 const auto X0 = View<double[3]>(X0_data.data());
171 const auto inputs = View<double[7]>(inputs_data.data());
172 const auto t_node_u = View<double[7]>(t_node_u_data.data());
173 const auto lambda = View<double[6]>(lambda_data.data());
174 const auto res = View<double[6]>(res_data.data());
175 const auto t_lambda_res = View<double[6]>(t_lambda_res_data.data());
176 const auto t_grad = View<double[6][6]>(t_grad_data.data());
177 const auto t_grad_trans = View<double[6][6]>(t_grad_trans_data.data());
178 const auto target_tangent = View<double[6][6]>(target_tangent_data.data());
179 const auto t_grad_tan = View<double[6][6]>(t_grad_tan_data.data());
180
181 const auto target_node_index = target_node_index_(constraint);
182
183 CopyVector::invoke(subview(X0_, constraint, ALL), X0);
184 CopyVector::invoke(subview(constraint_inputs_, constraint, ALL), inputs);
185 CopyVector::invoke(subview(node_u_, target_node_index, ALL), t_node_u);
186 CopyVector::invoke(subview(lambda_, constraint, ALL), lambda);
187 CopyMatrix::invoke(subview(tangent_, target_node_index, ALL, ALL), target_tangent);
188
189 CalculatePrescribedBCConstraint<DeviceType>::invoke(X0, inputs, t_node_u, res, t_grad);
190
191 TransposeMatrix::invoke(t_grad, t_grad_trans);
192 Gemm::invoke(1., t_grad, target_tangent, 0., t_grad_tan);
193 serial_gemv('N', 1., t_grad_trans, lambda, 0., t_lambda_res);
194
195 CopyVector::invoke(res, subview(res_, constraint, ALL));
196 CopyVector::invoke(t_lambda_res, subview(t_lambda_res_, constraint, ALL));
197 CopyMatrix::invoke(t_grad_tan, subview(t_grad_, constraint, ALL, ALL));
198 CopyMatrix::invoke(t_grad_trans, subview(t_grad_trans_, constraint, ALL, ALL));
199 }
200
201 KOKKOS_FUNCTION
202 void PrescribedBC3DOF(size_t constraint) const {
203 using Kokkos::ALL;
204 using Kokkos::Array;
205 using Kokkos::subview;
206 using KokkosBlas::Experimental::serial_gemv;
207
208 auto X0_data = Array<double, 3>{};
209 auto inputs_data = Array<double, 7>{};
210 auto t_node_u_data = Array<double, 7>{};
211 auto lambda_data = Array<double, 6>{};
212 auto res_data = Array<double, 6>{};
213 auto t_lambda_res_data = Array<double, 6>{};
214 auto t_grad_data = Array<double, 36>{};
215 auto t_grad_trans_data = Array<double, 36>{};
216 auto target_tangent_data = Array<double, 36>{};
217 auto t_grad_tan_data = Array<double, 36>{};
218
219 const auto X0 = View<double[3]>(X0_data.data());
220 const auto inputs = View<double[7]>(inputs_data.data());
221 const auto t_node_u = View<double[7]>(t_node_u_data.data());
222 const auto lambda = View<double[6]>(lambda_data.data());
223 const auto res = View<double[6]>(res_data.data());
224 const auto t_lambda_res = View<double[6]>(t_lambda_res_data.data());
225 const auto t_grad = View<double[6][6]>(t_grad_data.data());
226 const auto t_grad_trans = View<double[6][6]>(t_grad_trans_data.data());
227 const auto target_tangent = View<double[6][6]>(target_tangent_data.data());
228 const auto t_grad_tan = View<double[6][6]>(t_grad_tan_data.data());
229
230 const auto target_node_index = target_node_index_(constraint);
231
232 CopyVector::invoke(subview(X0_, constraint, ALL), X0);
233 CopyVector::invoke(subview(constraint_inputs_, constraint, ALL), inputs);
234 CopyVector::invoke(subview(node_u_, target_node_index, ALL), t_node_u);
235 CopyVector::invoke(subview(lambda_, constraint, ALL), lambda);
236 CopyMatrix::invoke(subview(tangent_, target_node_index, ALL, ALL), target_tangent);
237
238 CalculatePrescribedBC3DOFConstraint<DeviceType>::invoke(X0, inputs, t_node_u, res, t_grad);
239
240 TransposeMatrix::invoke(t_grad, t_grad_trans);
241 Gemm::invoke(1., t_grad, target_tangent, 0., t_grad_tan);
242 serial_gemv('N', 1., t_grad_trans, lambda, 0., t_lambda_res);
243
244 CopyVector::invoke(res, subview(res_, constraint, ALL));
245 CopyVector::invoke(t_lambda_res, subview(t_lambda_res_, constraint, ALL));
246 CopyMatrix::invoke(t_grad_tan, subview(t_grad_, constraint, ALL, ALL));
247 CopyMatrix::invoke(t_grad_trans, subview(t_grad_trans_, constraint, ALL, ALL));
248 }
249
250 KOKKOS_FUNCTION
251 void RigidJoint(size_t constraint) const {
252 using Kokkos::ALL;
253 using Kokkos::Array;
254 using Kokkos::subview;
255 using KokkosBlas::Experimental::serial_gemv;
256
257 auto X0_data = Array<double, 3>{};
258 auto b_node_u_data = Array<double, 7>{};
259 auto t_node_u_data = Array<double, 7>{};
260 auto lambda_data = Array<double, 6>{};
261 auto res_data = Array<double, 6>{};
262 auto b_lambda_res_data = Array<double, 6>{};
263 auto t_lambda_res_data = Array<double, 6>{};
264 auto b_grad_data = Array<double, 36>{};
265 auto t_grad_data = Array<double, 36>{};
266 auto b_grad_trans_data = Array<double, 36>{};
267 auto t_grad_trans_data = Array<double, 36>{};
268 auto base_tangent_data = Array<double, 36>{};
269 auto target_tangent_data = Array<double, 36>{};
270 auto b_grad_tan_data = Array<double, 36>{};
271 auto t_grad_tan_data = Array<double, 36>{};
272
273 const auto X0 = View<double[3]>(X0_data.data());
274 const auto b_node_u = View<double[7]>(b_node_u_data.data());
275 const auto t_node_u = View<double[7]>(t_node_u_data.data());
276 const auto lambda = View<double[6]>(lambda_data.data());
277 const auto res = View<double[6]>(res_data.data());
278 const auto b_lambda_res = View<double[6]>(b_lambda_res_data.data());
279 const auto t_lambda_res = View<double[6]>(t_lambda_res_data.data());
280 const auto b_grad = View<double[6][6]>(b_grad_data.data());
281 const auto t_grad = View<double[6][6]>(t_grad_data.data());
282 const auto b_grad_trans = View<double[6][6]>(b_grad_trans_data.data());
283 const auto t_grad_trans = View<double[6][6]>(t_grad_trans_data.data());
284 const auto base_tangent = View<double[6][6]>(base_tangent_data.data());
285 const auto target_tangent = View<double[6][6]>(target_tangent_data.data());
286 const auto b_grad_tan = View<double[6][6]>(b_grad_tan_data.data());
287 const auto t_grad_tan = View<double[6][6]>(t_grad_tan_data.data());
288
289 const auto base_node_index = base_node_index_(constraint);
290 const auto target_node_index = target_node_index_(constraint);
291
292 CopyVector::invoke(subview(X0_, constraint, ALL), X0);
293 CopyVector::invoke(subview(node_u_, base_node_index, ALL), b_node_u);
294 CopyVector::invoke(subview(node_u_, target_node_index, ALL), t_node_u);
295 CopyVector::invoke(subview(lambda_, constraint, ALL), lambda);
296 CopyMatrix::invoke(subview(tangent_, base_node_index, ALL, ALL), base_tangent);
297 CopyMatrix::invoke(subview(tangent_, target_node_index, ALL, ALL), target_tangent);
298
299 CalculateRigidJointConstraint<DeviceType>::invoke(
300 X0, b_node_u, t_node_u, res, b_grad, t_grad
301 );
302
303 TransposeMatrix::invoke(b_grad, b_grad_trans);
304 TransposeMatrix::invoke(t_grad, t_grad_trans);
305 Gemm::invoke(1., b_grad, base_tangent, 0., b_grad_tan);
306 Gemm::invoke(1., t_grad, target_tangent, 0., t_grad_tan);
307 serial_gemv('N', 1., b_grad_trans, lambda, 0., b_lambda_res);
308 serial_gemv('N', 1., t_grad_trans, lambda, 0., t_lambda_res);
309
310 CopyVector::invoke(res, subview(res_, constraint, ALL));
311 CopyVector::invoke(b_lambda_res, subview(b_lambda_res_, constraint, ALL));
312 CopyVector::invoke(t_lambda_res, subview(t_lambda_res_, constraint, ALL));
313 CopyMatrix::invoke(b_grad_tan, subview(b_grad_, constraint, ALL, ALL));
314 CopyMatrix::invoke(t_grad_tan, subview(t_grad_, constraint, ALL, ALL));
315 CopyMatrix::invoke(b_grad_trans, subview(b_grad_trans_, constraint, ALL, ALL));
316 CopyMatrix::invoke(t_grad_trans, subview(t_grad_trans_, constraint, ALL, ALL));
317 }
318
319 KOKKOS_FUNCTION
320 void RigidJoint3DOF(size_t constraint) const {
321 using Kokkos::ALL;
322 using Kokkos::Array;
323 using Kokkos::subview;
324 using KokkosBlas::Experimental::serial_gemv;
325
326 auto X0_data = Array<double, 3>{};
327 auto b_node_u_data = Array<double, 7>{};
328 auto t_node_u_data = Array<double, 7>{};
329 auto lambda_data = Array<double, 6>{};
330 auto res_data = Array<double, 6>{};
331 auto b_lambda_res_data = Array<double, 6>{};
332 auto t_lambda_res_data = Array<double, 6>{};
333 auto b_grad_data = Array<double, 36>{};
334 auto t_grad_data = Array<double, 36>{};
335 auto b_grad_trans_data = Array<double, 36>{};
336 auto t_grad_trans_data = Array<double, 36>{};
337 auto base_tangent_data = Array<double, 36>{};
338 auto target_tangent_data = Array<double, 36>{};
339 auto b_grad_tan_data = Array<double, 36>{};
340 auto t_grad_tan_data = Array<double, 36>{};
341
342 const auto X0 = View<double[3]>(X0_data.data());
343 const auto b_node_u = View<double[7]>(b_node_u_data.data());
344 const auto t_node_u = View<double[7]>(t_node_u_data.data());
345 const auto lambda = View<double[6]>(lambda_data.data());
346 const auto res = View<double[6]>(res_data.data());
347 const auto b_lambda_res = View<double[6]>(b_lambda_res_data.data());
348 const auto t_lambda_res = View<double[6]>(t_lambda_res_data.data());
349 const auto b_grad = View<double[6][6]>(b_grad_data.data());
350 const auto t_grad = View<double[6][6]>(t_grad_data.data());
351 const auto b_grad_trans = View<double[6][6]>(b_grad_trans_data.data());
352 const auto t_grad_trans = View<double[6][6]>(t_grad_trans_data.data());
353 const auto base_tangent = View<double[6][6]>(base_tangent_data.data());
354 const auto target_tangent = View<double[6][6]>(target_tangent_data.data());
355 const auto b_grad_tan = View<double[6][6]>(b_grad_tan_data.data());
356 const auto t_grad_tan = View<double[6][6]>(t_grad_tan_data.data());
357
358 const auto base_node_index = base_node_index_(constraint);
359 const auto target_node_index = target_node_index_(constraint);
360
361 CopyVector::invoke(subview(X0_, constraint, ALL), X0);
362 CopyVector::invoke(subview(node_u_, base_node_index, ALL), b_node_u);
363 CopyVector::invoke(subview(node_u_, target_node_index, ALL), t_node_u);
364 CopyVector::invoke(subview(lambda_, constraint, ALL), lambda);
365 CopyMatrix::invoke(subview(tangent_, base_node_index, ALL, ALL), base_tangent);
366 CopyMatrix::invoke(subview(tangent_, target_node_index, ALL, ALL), target_tangent);
367
368 CalculateRigidJoint3DOFConstraint<DeviceType>::invoke(
369 X0, b_node_u, t_node_u, res, b_grad, t_grad
370 );
371
372 TransposeMatrix::invoke(b_grad, b_grad_trans);
373 TransposeMatrix::invoke(t_grad, t_grad_trans);
374 Gemm::invoke(1., b_grad, base_tangent, 0., b_grad_tan);
375 Gemm::invoke(1., t_grad, target_tangent, 0., t_grad_tan);
376 serial_gemv('N', 1., b_grad_trans, lambda, 0., b_lambda_res);
377 serial_gemv('N', 1., t_grad_trans, lambda, 0., t_lambda_res);
378
379 CopyVector::invoke(res, subview(res_, constraint, ALL));
380 CopyVector::invoke(b_lambda_res, subview(b_lambda_res_, constraint, ALL));
381 CopyVector::invoke(t_lambda_res, subview(t_lambda_res_, constraint, ALL));
382 CopyMatrix::invoke(b_grad_tan, subview(b_grad_, constraint, ALL, ALL));
383 CopyMatrix::invoke(t_grad_tan, subview(t_grad_, constraint, ALL, ALL));
384 CopyMatrix::invoke(b_grad_trans, subview(b_grad_trans_, constraint, ALL, ALL));
385 CopyMatrix::invoke(t_grad_trans, subview(t_grad_trans_, constraint, ALL, ALL));
386 }
387
388 KOKKOS_FUNCTION
389 void RevoluteJoint(size_t constraint) const {
390 using Kokkos::ALL;
391 using Kokkos::Array;
392 using Kokkos::subview;
393 using KokkosBlas::Experimental::serial_gemv;
394
395 auto X0_data = Array<double, 3>{};
396 auto inputs_data = Array<double, 7>{};
397 auto b_node_u_data = Array<double, 7>{};
398 auto t_node_u_data = Array<double, 7>{};
399 auto lambda_data = Array<double, 6>{};
400 auto axes_data = Array<double, 9>{};
401 auto res_data = Array<double, 6>{};
402 auto b_lambda_res_data = Array<double, 6>{};
403 auto t_lambda_res_data = Array<double, 6>{};
404 auto system_res_data = Array<double, 6>{};
405 auto b_grad_data = Array<double, 36>{};
406 auto t_grad_data = Array<double, 36>{};
407 auto b_grad_trans_data = Array<double, 36>{};
408 auto t_grad_trans_data = Array<double, 36>{};
409 auto base_tangent_data = Array<double, 36>{};
410 auto target_tangent_data = Array<double, 36>{};
411 auto b_grad_tan_data = Array<double, 36>{};
412 auto t_grad_tan_data = Array<double, 36>{};
413
414 const auto X0 = View<double[3]>(X0_data.data());
415 const auto inputs = View<double[7]>(inputs_data.data());
416 const auto b_node_u = View<double[7]>(b_node_u_data.data());
417 const auto t_node_u = View<double[7]>(t_node_u_data.data());
418 const auto lambda = View<double[6]>(lambda_data.data());
419 const auto axes = View<double[3][3]>(axes_data.data());
420 const auto res = View<double[6]>(res_data.data());
421 const auto b_lambda_res = View<double[6]>(b_lambda_res_data.data());
422 const auto t_lambda_res = View<double[6]>(t_lambda_res_data.data());
423 const auto system_res = View<double[6]>(system_res_data.data());
424 const auto b_grad = View<double[6][6]>(b_grad_data.data());
425 const auto t_grad = View<double[6][6]>(t_grad_data.data());
426 const auto b_grad_trans = View<double[6][6]>(b_grad_trans_data.data());
427 const auto t_grad_trans = View<double[6][6]>(t_grad_trans_data.data());
428 const auto base_tangent = View<double[6][6]>(base_tangent_data.data());
429 const auto target_tangent = View<double[6][6]>(target_tangent_data.data());
430 const auto b_grad_tan = View<double[6][6]>(b_grad_tan_data.data());
431 const auto t_grad_tan = View<double[6][6]>(t_grad_tan_data.data());
432
433 const auto base_node_index = base_node_index_(constraint);
434 const auto target_node_index = target_node_index_(constraint);
435
436 CopyVector::invoke(subview(X0_, constraint, ALL), X0);
437 CopyVector::invoke(subview(constraint_inputs_, constraint, ALL), inputs);
438 CopyVector::invoke(subview(node_u_, base_node_index, ALL), b_node_u);
439 CopyVector::invoke(subview(node_u_, target_node_index, ALL), t_node_u);
440 CopyVector::invoke(subview(lambda_, constraint, ALL), lambda);
441 CopyMatrix::invoke(subview(axes_, constraint, ALL, ALL), axes);
442 CopyMatrix::invoke(subview(tangent_, base_node_index, ALL, ALL), base_tangent);
443 CopyMatrix::invoke(subview(tangent_, target_node_index, ALL, ALL), target_tangent);
444
445 CalculateRevoluteJointConstraint<DeviceType>::invoke(
446 X0, axes, b_node_u, t_node_u, res, b_grad, t_grad
447 );
448 CalculateRevoluteJointForce<DeviceType>::invoke(axes, inputs, t_node_u, system_res);
449
450 TransposeMatrix::invoke(b_grad, b_grad_trans);
451 TransposeMatrix::invoke(t_grad, t_grad_trans);
452 Gemm::invoke(1., b_grad, base_tangent, 0., b_grad_tan);
453 Gemm::invoke(1., t_grad, target_tangent, 0., t_grad_tan);
454 serial_gemv('N', 1., b_grad_trans, lambda, 0., b_lambda_res);
455 serial_gemv('N', 1., t_grad_trans, lambda, 0., t_lambda_res);
456
457 CopyVector::invoke(res, subview(res_, constraint, ALL));
458 CopyVector::invoke(b_lambda_res, subview(b_lambda_res_, constraint, ALL));
459 CopyVector::invoke(t_lambda_res, subview(t_lambda_res_, constraint, ALL));
460 CopyVector::invoke(system_res, subview(system_res_, constraint, ALL));
461 CopyMatrix::invoke(b_grad_tan, subview(b_grad_, constraint, ALL, ALL));
462 CopyMatrix::invoke(t_grad_tan, subview(t_grad_, constraint, ALL, ALL));
463 CopyMatrix::invoke(b_grad_trans, subview(b_grad_trans_, constraint, ALL, ALL));
464 CopyMatrix::invoke(t_grad_trans, subview(t_grad_trans_, constraint, ALL, ALL));
465 }
466
467 KOKKOS_FUNCTION
468 void RotationControl(size_t constraint) const {
469 using Kokkos::ALL;
470 using Kokkos::Array;
471 using Kokkos::subview;
472 using KokkosBlas::Experimental::serial_gemv;
473
474 auto X0_data = Array<double, 3>{};
475 auto inputs_data = Array<double, 7>{};
476 auto b_node_u_data = Array<double, 7>{};
477 auto t_node_u_data = Array<double, 7>{};
478 auto lambda_data = Array<double, 6>{};
479 auto axes_data = Array<double, 9>{};
480 auto res_data = Array<double, 6>{};
481 auto b_lambda_res_data = Array<double, 6>{};
482 auto t_lambda_res_data = Array<double, 6>{};
483 auto b_grad_data = Array<double, 36>{};
484 auto t_grad_data = Array<double, 36>{};
485 auto b_grad_trans_data = Array<double, 36>{};
486 auto t_grad_trans_data = Array<double, 36>{};
487 auto base_tangent_data = Array<double, 36>{};
488 auto target_tangent_data = Array<double, 36>{};
489 auto b_grad_tan_data = Array<double, 36>{};
490 auto t_grad_tan_data = Array<double, 36>{};
491
492 const auto X0 = View<double[3]>(X0_data.data());
493 const auto inputs = View<double[7]>(inputs_data.data());
494 const auto b_node_u = View<double[7]>(b_node_u_data.data());
495 const auto t_node_u = View<double[7]>(t_node_u_data.data());
496 const auto lambda = View<double[6]>(lambda_data.data());
497 const auto axes = View<double[3][3]>(axes_data.data());
498 const auto res = View<double[6]>(res_data.data());
499 const auto b_lambda_res = View<double[6]>(b_lambda_res_data.data());
500 const auto t_lambda_res = View<double[6]>(t_lambda_res_data.data());
501 const auto b_grad = View<double[6][6]>(b_grad_data.data());
502 const auto t_grad = View<double[6][6]>(t_grad_data.data());
503 const auto b_grad_trans = View<double[6][6]>(b_grad_trans_data.data());
504 const auto t_grad_trans = View<double[6][6]>(t_grad_trans_data.data());
505 const auto base_tangent = View<double[6][6]>(base_tangent_data.data());
506 const auto target_tangent = View<double[6][6]>(target_tangent_data.data());
507 const auto b_grad_tan = View<double[6][6]>(b_grad_tan_data.data());
508 const auto t_grad_tan = View<double[6][6]>(t_grad_tan_data.data());
509
510 const auto base_node_index = base_node_index_(constraint);
511 const auto target_node_index = target_node_index_(constraint);
512
513 CopyVector::invoke(subview(X0_, constraint, ALL), X0);
514 CopyVector::invoke(subview(constraint_inputs_, constraint, ALL), inputs);
515 CopyVector::invoke(subview(node_u_, base_node_index, ALL), b_node_u);
516 CopyVector::invoke(subview(node_u_, target_node_index, ALL), t_node_u);
517 CopyVector::invoke(subview(lambda_, constraint, ALL), lambda);
518 CopyMatrix::invoke(subview(axes_, constraint, ALL, ALL), axes);
519 CopyMatrix::invoke(subview(tangent_, base_node_index, ALL, ALL), base_tangent);
520 CopyMatrix::invoke(subview(tangent_, target_node_index, ALL, ALL), target_tangent);
521
522 CalculateRotationControlConstraint<DeviceType>::invoke(
523 X0, axes, inputs, b_node_u, t_node_u, res, b_grad, t_grad
524 );
525
526 TransposeMatrix::invoke(b_grad, b_grad_trans);
527 TransposeMatrix::invoke(t_grad, t_grad_trans);
528 Gemm::invoke(1., b_grad, base_tangent, 0., b_grad_tan);
529 Gemm::invoke(1., t_grad, target_tangent, 0., t_grad_tan);
530 serial_gemv('N', 1., b_grad_trans, lambda, 0., b_lambda_res);
531 serial_gemv('N', 1., t_grad_trans, lambda, 0., t_lambda_res);
532
533 CopyVector::invoke(res, subview(res_, constraint, ALL));
534 CopyVector::invoke(b_lambda_res, subview(b_lambda_res_, constraint, ALL));
535 CopyVector::invoke(t_lambda_res, subview(t_lambda_res_, constraint, ALL));
536 CopyMatrix::invoke(b_grad_tan, subview(b_grad_, constraint, ALL, ALL));
537 CopyMatrix::invoke(t_grad_tan, subview(t_grad_, constraint, ALL, ALL));
538 CopyMatrix::invoke(b_grad_trans, subview(b_grad_trans_, constraint, ALL, ALL));
539 CopyMatrix::invoke(t_grad_trans, subview(t_grad_trans_, constraint, ALL, ALL));
540 }
541
542 KOKKOS_FUNCTION
543 void operator()(size_t constraint) const {
544 const auto constraint_type = type_(constraint);
545
546 if (constraint_type == ConstraintType::FixedBC) {
547 FixedBC(constraint);
548 } else if (constraint_type == ConstraintType::FixedBC3DOFs) {
549 FixedBC3DOF(constraint);
550 } else if (constraint_type == ConstraintType::PrescribedBC) {
551 PrescribedBC(constraint);
552 } else if (constraint_type == ConstraintType::PrescribedBC3DOFs) {
553 PrescribedBC3DOF(constraint);
554 } else if (constraint_type == ConstraintType::RigidJoint) {
555 RigidJoint(constraint);
556 } else if (constraint_type == ConstraintType::RigidJoint6DOFsTo3DOFs) {
557 RigidJoint3DOF(constraint);
558 } else if (constraint_type == ConstraintType::RevoluteJoint) {
559 RevoluteJoint(constraint);
560 } else if (constraint_type == ConstraintType::RotationControl) {
561 RotationControl(constraint);
562 }
563 };
564};
565
566} // namespace kynema::constraints
calculate_fixed_bc_3DOF_constraint.hpp
calculate_fixed_bc_constraint.hpp
calculate_prescribed_bc_3DOF_constraint.hpp
calculate_prescribed_bc_constraint.hpp
calculate_revolute_joint_constraint.hpp
calculate_revolute_joint_force.hpp
calculate_rigid_joint_3DOF_constraint.hpp
calculate_rigid_joint_constraint.hpp
calculate_rotation_control_constraint.hpp
constraint_type.hpp
kynema::constraints
Definition calculate_constraint_output.hpp:8
kynema::constraints::CalculateConstraintResidualGradient
Top level kernel which calculates the residual and gradient contributions of a constraint.
Definition calculate_constraint_residual_gradient.hpp:30
kynema::constraints::CalculateConstraintResidualGradient::b_grad_
View< double *[6][6]> b_grad_
Definition calculate_constraint_residual_gradient.hpp:55
kynema::constraints::CalculateConstraintResidualGradient::type_
ConstView< ConstraintType * > type_
Definition calculate_constraint_residual_gradient.hpp:42
kynema::constraints::CalculateConstraintResidualGradient::t_grad_trans_
View< double *[6][6]> t_grad_trans_
Definition calculate_constraint_residual_gradient.hpp:58
kynema::constraints::CalculateConstraintResidualGradient::RevoluteJoint
KOKKOS_FUNCTION void RevoluteJoint(size_t constraint) const
Definition calculate_constraint_residual_gradient.hpp:389
kynema::constraints::CalculateConstraintResidualGradient::system_res_
View< double *[6]> system_res_
Definition calculate_constraint_residual_gradient.hpp:54
kynema::constraints::CalculateConstraintResidualGradient::axes_
ConstView< double *[3][3]> axes_
Definition calculate_constraint_residual_gradient.hpp:46
kynema::constraints::CalculateConstraintResidualGradient::PrescribedBC
KOKKOS_FUNCTION void PrescribedBC(size_t constraint) const
Definition calculate_constraint_residual_gradient.hpp:153
kynema::constraints::CalculateConstraintResidualGradient::t_lambda_res_
View< double *[6]> t_lambda_res_
Definition calculate_constraint_residual_gradient.hpp:53
kynema::constraints::CalculateConstraintResidualGradient::node_u_
ConstView< double *[7]> node_u_
Definition calculate_constraint_residual_gradient.hpp:50
kynema::constraints::CalculateConstraintResidualGradient::b_lambda_res_
View< double *[6]> b_lambda_res_
Definition calculate_constraint_residual_gradient.hpp:52
kynema::constraints::CalculateConstraintResidualGradient::ConstView
typename View< ValueType >::const_type ConstView
Definition calculate_constraint_residual_gradient.hpp:40
kynema::constraints::CalculateConstraintResidualGradient::RotationControl
KOKKOS_FUNCTION void RotationControl(size_t constraint) const
Definition calculate_constraint_residual_gradient.hpp:468
kynema::constraints::CalculateConstraintResidualGradient::PrescribedBC3DOF
KOKKOS_FUNCTION void PrescribedBC3DOF(size_t constraint) const
Definition calculate_constraint_residual_gradient.hpp:202
kynema::constraints::CalculateConstraintResidualGradient::target_node_index_
ConstView< size_t * > target_node_index_
Definition calculate_constraint_residual_gradient.hpp:44
kynema::constraints::CalculateConstraintResidualGradient::X0_
ConstView< double *[3]> X0_
Definition calculate_constraint_residual_gradient.hpp:45
kynema::constraints::CalculateConstraintResidualGradient::FixedBC3DOF
KOKKOS_FUNCTION void FixedBC3DOF(size_t constraint) const
Definition calculate_constraint_residual_gradient.hpp:107
kynema::constraints::CalculateConstraintResidualGradient::CopyVector
KokkosBatched::SerialCopy< KokkosBatched::Trans::NoTranspose, 1 > CopyVector
Definition calculate_constraint_residual_gradient.hpp:35
kynema::constraints::CalculateConstraintResidualGradient::lambda_
ConstView< double *[6]> lambda_
Definition calculate_constraint_residual_gradient.hpp:48
kynema::constraints::CalculateConstraintResidualGradient::constraint_inputs_
ConstView< double *[7]> constraint_inputs_
Definition calculate_constraint_residual_gradient.hpp:47
kynema::constraints::CalculateConstraintResidualGradient::b_grad_trans_
View< double *[6][6]> b_grad_trans_
Definition calculate_constraint_residual_gradient.hpp:57
kynema::constraints::CalculateConstraintResidualGradient::View
Kokkos::View< ValueType, DeviceType > View
Definition calculate_constraint_residual_gradient.hpp:38
kynema::constraints::CalculateConstraintResidualGradient::base_node_index_
ConstView< size_t * > base_node_index_
Definition calculate_constraint_residual_gradient.hpp:43
kynema::constraints::CalculateConstraintResidualGradient::RigidJoint
KOKKOS_FUNCTION void RigidJoint(size_t constraint) const
Definition calculate_constraint_residual_gradient.hpp:251
kynema::constraints::CalculateConstraintResidualGradient::FixedBC
KOKKOS_FUNCTION void FixedBC(size_t constraint) const
Definition calculate_constraint_residual_gradient.hpp:61
kynema::constraints::CalculateConstraintResidualGradient::t_grad_
View< double *[6][6]> t_grad_
Definition calculate_constraint_residual_gradient.hpp:56
kynema::constraints::CalculateConstraintResidualGradient::tangent_
ConstView< double *[6][6]> tangent_
Definition calculate_constraint_residual_gradient.hpp:49
kynema::constraints::CalculateConstraintResidualGradient::Gemm
KokkosBatched::SerialGemm< KokkosBatched::Trans::NoTranspose, KokkosBatched::Trans::NoTranspose, KokkosBatched::Algo::Gemm::Default > Gemm
Definition calculate_constraint_residual_gradient.hpp:34
kynema::constraints::CalculateConstraintResidualGradient::res_
View< double *[6]> res_
Definition calculate_constraint_residual_gradient.hpp:51
kynema::constraints::CalculateConstraintResidualGradient::TransposeMatrix
KokkosBatched::SerialCopy< KokkosBatched::Trans::Transpose > TransposeMatrix
Definition calculate_constraint_residual_gradient.hpp:31
kynema::constraints::CalculateConstraintResidualGradient::operator()
KOKKOS_FUNCTION void operator()(size_t constraint) const
Definition calculate_constraint_residual_gradient.hpp:543
kynema::constraints::CalculateConstraintResidualGradient::RigidJoint3DOF
KOKKOS_FUNCTION void RigidJoint3DOF(size_t constraint) const
Definition calculate_constraint_residual_gradient.hpp:320
kynema::constraints::CalculateConstraintResidualGradient::CopyMatrix
KokkosBatched::SerialCopy< KokkosBatched::Trans::NoTranspose, 2 > CopyMatrix
Definition calculate_constraint_residual_gradient.hpp:36
kynema::constraints::CalculateFixedBC3DOFConstraint::invoke
static KOKKOS_FUNCTION void invoke(const ConstView< double[3]> &X0, const ConstView< double[7]> &node_u, const View< double[6]> &residual_terms, const View< double[6][6]> &target_gradient_terms)
Definition calculate_fixed_bc_3DOF_constraint.hpp:20
kynema::constraints::CalculateFixedBCConstraint::invoke
static KOKKOS_FUNCTION void invoke(const ConstView< double[3]> &X0, const ConstView< double[7]> &node_u, const View< double[6]> &residual_terms, const View< double[6][6]> &target_gradient_terms)
Definition calculate_fixed_bc_constraint.hpp:22
kynema::constraints::CalculatePrescribedBC3DOFConstraint::invoke
static KOKKOS_FUNCTION void invoke(const ConstView< double[3]> &X0, const ConstView< double[7]> &inputs, const ConstView< double[7]> &node_u, const View< double[6]> &residual_terms, const View< double[6][6]> &target_gradient_terms)
Definition calculate_prescribed_bc_3DOF_constraint.hpp:20
kynema::constraints::CalculatePrescribedBCConstraint::invoke
static KOKKOS_FUNCTION void invoke(const ConstView< double[3]> &X0, const ConstView< double[7]> &inputs, const ConstView< double[7]> &node_u, const View< double[6]> &residual_terms, const View< double[6][6]> &target_gradient_terms)
Definition calculate_prescribed_bc_constraint.hpp:22
kynema::constraints::CalculateRevoluteJointConstraint::invoke
static KOKKOS_FUNCTION void invoke(const ConstView< double[3]> &X0, const ConstView< double[3][3]> &axes, const ConstView< double[7]> &base_node_u, const ConstView< double[7]> &target_node_u, const View< double[6]> &residual_terms, const View< double[6][6]> &base_gradient_terms, const View< double[6][6]> &target_gradient_terms)
Definition calculate_revolute_joint_constraint.hpp:20
kynema::constraints::CalculateRevoluteJointForce::invoke
static KOKKOS_FUNCTION void invoke(const ConstView< double[3][3]> &axes, const ConstView< double[7]> &inputs, const ConstView< double[7]> &node_u, const View< double[6]> &system_residual_terms)
Definition calculate_revolute_joint_force.hpp:20
kynema::constraints::CalculateRigidJoint3DOFConstraint::invoke
static KOKKOS_FUNCTION void invoke(const ConstView< double[3]> &X0, const ConstView< double[7]> &base_node_u, const ConstView< double[7]> &target_node_u, const View< double[6]> &residual_terms, const View< double[6][6]> &base_gradient_terms, const View< double[6][6]> &target_gradient_terms)
Definition calculate_rigid_joint_3DOF_constraint.hpp:22
kynema::constraints::CalculateRigidJointConstraint::invoke
static KOKKOS_FUNCTION void invoke(const ConstView< double[3]> &X0, const ConstView< double[7]> &base_node_u, const ConstView< double[7]> &target_node_u, const View< double[6]> &residual_terms, const View< double[6][6]> &base_gradient_terms, const View< double[6][6]> &target_gradient_terms)
Definition calculate_rigid_joint_constraint.hpp:23
kynema::constraints::CalculateRotationControlConstraint::invoke
static KOKKOS_FUNCTION void invoke(const ConstView< double[3]> &X0, const ConstView< double[3][3]> &axes, const ConstView< double[7]> &inputs, const ConstView< double[7]> &base_node_u, const ConstView< double[7]> &target_node_u, const View< double[6]> &residual_terms, const View< double[6][6]> &base_gradient_terms, const View< double[6][6]> &target_gradient_terms)
Definition calculate_rotation_control_constraint.hpp:23
Generated by doxygen 1.9.8