| 1647 | |
| 1648 | template<int DIM> |
| 1649 | void MPContactElement<DIM>::CalcLinearizedElementMatrix(FlatVector<double> elx, |
| 1650 | FlatMatrix<double> elmat, |
| 1651 | LocalHeap& lh) const |
| 1652 | { |
| 1653 | // NETGEN_TIMER_FROM_HERE("MPContactElement::CalcLinearizedElementMatrix"); |
| 1654 | |
| 1655 | elmat = 0.; |
| 1656 | |
| 1657 | auto& primary_trafo = trial_fes->GetMeshAccess()->GetTrafo(primary_ei, lh); |
| 1658 | auto& secondary_trafo = trial_fes->GetMeshAccess()->GetTrafo(secondary_ei, lh); |
| 1659 | auto& primary_deformed_trafo = primary_trafo.AddDeformation(deformation, lh); |
| 1660 | auto& secondary_deformed_trafo = secondary_trafo.AddDeformation(deformation, lh); |
| 1661 | |
| 1662 | // IntegrationRule primary_ir(1, &const_cast<IntegrationPoint&>(primary_e)); |
| 1663 | // IntegrationRule secondary_ir(1, &const_cast<IntegrationPoint&>(pair.secondary_ip)); |
| 1664 | |
| 1665 | auto& primary_trial_fel = trial_fes->GetFE(primary_ei, lh); |
| 1666 | auto& secondary_trial_fel = trial_fes->GetFE(secondary_ei, lh); |
| 1667 | auto& primary_test_fel = test_fes->GetFE(primary_ei, lh); |
| 1668 | auto& secondary_test_fel = test_fes->GetFE(secondary_ei, lh); |
| 1669 | |
| 1670 | if (primary_ei.IsBoundary()) |
| 1671 | for (bool def : { false, true }) |
| 1672 | { |
| 1673 | auto& energies = cb->GetEnergies(def); |
| 1674 | auto& integrators = cb->GetIntegrators(def); |
| 1675 | if (energies.Size() || integrators.Size()) |
| 1676 | { |
| 1677 | MappedIntegrationRule<DIM-1,DIM> primary_mir(primary_ir, def ? primary_deformed_trafo : primary_trafo, lh); |
| 1678 | MappedIntegrationRule<DIM-1,DIM> secondary_mir(secondary_ir, def ? secondary_deformed_trafo : secondary_trafo, lh); |
| 1679 | primary_mir.SetOtherMIR(&secondary_mir); |
| 1680 | |
| 1681 | for(const auto& ci : integrators) |
| 1682 | ci->CalcLinearizedAdd(primary_trial_fel, secondary_trial_fel, primary_test_fel, secondary_test_fel, primary_mir, elx, elmat, lh); |
| 1683 | for(const auto& ce : energies) |
| 1684 | ce->CalcLinearizedAdd(primary_trial_fel, secondary_trial_fel, primary_mir, elx, elmat, lh); |
| 1685 | } |
| 1686 | } |
| 1687 | else |
| 1688 | { |
| 1689 | // cout << "volume contact integrator" << endl; |
| 1690 | for (bool def : { false, true }) |
| 1691 | { |
| 1692 | auto& energies = cb->GetEnergies(def); |
| 1693 | auto& integrators = cb->GetIntegrators(def); |
| 1694 | if (energies.Size() || integrators.Size()) |
| 1695 | { |
| 1696 | MappedIntegrationRule<DIM,DIM> primary_mir(primary_ir, def ? primary_deformed_trafo : primary_trafo, lh); |
| 1697 | primary_mir.ComputeNormalsAndMeasure (primary_trafo.GetElementType(), primary_ir[0].FacetNr()); |
| 1698 | |
| 1699 | MappedIntegrationRule<DIM-1,DIM> secondary_mir(secondary_ir, def ? secondary_deformed_trafo : secondary_trafo, lh); |
| 1700 | |
| 1701 | /* |
| 1702 | cout << "primary mir: " << endl; |
| 1703 | for (const auto & mip : primary_mir) |
| 1704 | cout << mip << endl; |
| 1705 | cout << "secondary mir: " << endl; |
| 1706 | for (const auto & mip : secondary_mir) |
no test coverage detected