455 {
457
458 const size_t nb_gauss_pts = getGaussPts().size2();
459 auto &locMat = AssemblyBoundaryEleOp::locMat;
460
461 auto t_normal = getFTensor1Normal();
462 t_normal(
i) /= sqrt(t_normal(
j) * t_normal(
j));
463
464 auto t_disp = getFTensor1FromMat<SPACE_DIM>(
commonDataPtr->contactDisp);
465 auto t_traction =
466 getFTensor1FromMat<SPACE_DIM>(
commonDataPtr->contactTraction);
467 auto t_coords = getFTensor1CoordsAtGaussPts();
468
469 auto t_w = getFTensor0IntegrationWeight();
470 auto t_row_base = row_data.getFTensor1N<3>();
471 size_t nb_face_functions = row_data.getN().size2() / 3;
472
473 for (size_t gg = 0; gg != nb_gauss_pts; ++gg) {
474
475 const double alpha = t_w * getMeasure();
476
478 t_spatial_coords(
i) = t_coords(
i) + t_disp(
i);
479
483 getTStime(), t_spatial_coords, t_traction);
484 auto tn = -t_traction(
i) * t_grad_sdf(
i);
486
488 t_cP(
i,
j) = (
c * t_grad_sdf(
i)) * t_grad_sdf(
j);
491
493 t_res_dt(
i,
j) = -
cn * t_cQ(
i,
j);
494
495 size_t rr = 0;
496 for (; rr != AssemblyBoundaryEleOp::nbRows /
SPACE_DIM; ++rr) {
497
498 auto t_mat = getFTensor2FromArray<SPACE_DIM, SPACE_DIM, SPACE_DIM>(
500 const double row_base = t_row_base(
i) * t_normal(
i);
501
502 auto t_col_base = col_data.getFTensor1N<3>(gg, 0);
503 for (
size_t cc = 0; cc != AssemblyBoundaryEleOp::nbCols /
SPACE_DIM;
504 ++cc) {
505 const double col_base = t_col_base(
i) * t_normal(
i);
506 const double beta =
alpha * row_base * col_base;
507
508 t_mat(
i,
j) -= beta * t_res_dt(
i,
j);
509
510 ++t_col_base;
511 ++t_mat;
512 }
513
514 ++t_row_base;
515 }
516 for (; rr < nb_face_functions; ++rr)
517 ++t_row_base;
518
519 ++t_disp;
520 ++t_traction;
521 ++t_coords;
522 ++t_w;
523 }
524
526}
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
const double c
speed of light (cm/ns)
Tensor2_Expr< Kronecker_Delta< T >, T, Dim0, Dim1, i, j > kronecker_delta(const Index< i, Dim0 > &, const Index< j, Dim1 > &)
Rank 2.