564 auto nb_rows = row_data.getIndices().size();
565 auto nb_cols = col_data.getIndices().size();
568 locMat.resize(nb_rows, nb_cols,
false);
571 if (nb_cols && nb_rows) {
573 auto nb_gauss_pts = getGaussPts().size2();
574 auto t_w = getFTensor0IntegrationWeight();
575 auto t_coords = getFTensor1CoordsAtGaussPts();
576 auto t_disp = getFTensor1FromMat<3>(
commonDataPtr->contactDisp);
577 auto t_traction = getFTensor1FromMat<3>(
commonDataPtr->contactTraction);
580 auto [block_id, m_normals_at_pts, v_sdf, m_grad_sdf, m_hess_sdf] =
585 auto t_grad_sdf_v = getFTensor1FromMat<3>(m_grad_sdf);
586 auto t_hess_sdf_v = getFTensor2SymmetricFromMat<3>(m_hess_sdf);
587 auto t_normalized_normal = getFTensor1FromMat<3>(m_normals_at_pts);
597 ++t_normalized_normal;
600 auto face_data_vec_ptr =
602 auto face_gauss_pts_it = face_data_vec_ptr->begin();
604 auto t_row_base = row_data.getFTensor0N();
605 auto nb_face_functions = row_data.getN().size2() / 3;
608 for (
size_t gg = 0; gg != nb_gauss_pts; ++gg) {
610 auto face_data_ptr =
contactTreePtr->getFaceDataPtr(face_gauss_pts_it, gg,
613 auto check_face_contact = [&]() {
623 auto tn = t_traction(
i) * t_grad_sdf_v(
i);
626 constexpr
double c = 0;
629 if (!
c && check_face_contact()) {
631 t_spatial_coords(
i) = t_coords(
i) + t_disp(
i);
632 constexpr
double eps = std::numeric_limits<float>::epsilon();
633 for (
auto ii = 0; ii < 3; ++ii) {
635 t_spatial_coords(0), t_spatial_coords(1), t_spatial_coords(2)};
636 t_spatial_coords_cx(ii) +=
eps * 1
i;
640 for (
int jj = 0; jj != 3; ++jj) {
641 auto v = t_rhs_tmp(jj).imag();
642 t_res_dU(jj, ii) =
v /
eps;
650 if (ContactOps::sdfPythonWeakPtr.lock()) {
654 (-
c) * (t_hess_sdf_v(
i,
j) * t_grad_sdf_v(
k) * t_traction(
k) +
655 t_grad_sdf_v(
i) * t_hess_sdf_v(
k,
j) * t_traction(
k))
657 + (
c * inv_cn) * (t_sdf_v * t_hess_sdf_v(
i,
j) +
659 t_grad_sdf_v(
j) * t_grad_sdf_v(
i));
668 auto alpha = t_w * getMeasure();
671 for (; rr != nb_rows / 3; ++rr) {
673 auto t_mat = getFTensor2FromArray<3, 3, 3>(locMat, 3 * rr);
674 auto t_col_base = col_data.getFTensor0N(gg, 0);
676 for (
size_t cc = 0; cc != nb_cols / 3; ++cc) {
677 auto beta = alpha * t_row_base * t_col_base;
678 t_mat(
i,
j) -= beta * t_res_dU(
i,
j);
685 for (; rr < nb_face_functions; ++rr)