fstr_mat_con_contact.f90

Go to the documentation of this file.

!-------------------------------------------------------------------------------
! Copyright (c) 2019 FrontISTR Commons
! This software is released under the MIT License, see LICENSE.txt
!-------------------------------------------------------------------------------
!         stiffness matrix structure for the contact analysis
!         employing standard Lagrange multiplier algorithm
 
module fstr_matrix_con_contact
 
  use m_fstr
  use elementinfo
 
  implicit none
  private
  public :: fstr_get_num_lagrange_pernode
  public :: hecmwst_matrix_lagrange
  public :: fstr_save_originalmatrixstructure
  public :: fstr_mat_con_contact
  public :: fstr_is_matrixstruct_symmetric
 
  integer(kind=kint), save         :: NPL_org, NPU_org
  type(nodeRelated), pointer, save :: list_nodeRelated_org(:) => null() 
 
  type(nodeRelated), pointer       :: list_nodeRelated(:) => null() 
 
  logical                          :: permission = .false.
 
contains
 
  integer(kind=kint) function fstr_get_num_lagrange_pernode(algtype)
    integer(kind=kint) :: algtype
    if( algtype == contactsslid .or. algtype == contactfslid ) then
      fstr_get_num_lagrange_pernode = 1
    else if( algtype == contacttied ) then
      fstr_get_num_lagrange_pernode = 3
    endif
  end function
 
  subroutine fstr_save_originalmatrixstructure(hecMAT)
 
    type(hecmwst_matrix) :: hecmat
 
    if( associated(list_noderelated_org) ) return
    call hecmw_construct_noderelated_from_hecmat(hecmat, npl_org, npu_org, list_noderelated_org)
 
  end subroutine fstr_save_originalmatrixstructure
 
  subroutine fstr_mat_con_contact(cstep,contact_algo,hecMAT,fstrSOLID,hecLagMAT,infoCTChange,conMAT,is_contact_active)
 
    integer(kind=kint)                   :: cstep
    integer(kind=kint)                   :: contact_algo
    type(hecmwst_matrix)                 :: hecmat
    type(fstr_solid)                     :: fstrsolid
    type(hecmwst_matrix_lagrange) :: heclagmat
    type(fstr_info_contactchange)        :: infoctchange
 
    integer(kind=kint)                   :: num_lagrange
    integer(kind=kint)                   :: countnon0lu_node, countnon0lu_lagrange
    integer(kind=kint)                   :: numnon0_node, numnon0_lagrange
    type (hecmwst_matrix)                :: conmat
    logical, intent(in)                  :: is_contact_active
 
    integer(kind=kint)                   :: i, j, grpid
    integer(kind=kint)                   :: nlag
 
    num_lagrange = 0
    if( contact_algo == kcaslagrange ) then
      do i = 1, fstrsolid%n_contacts
        grpid = fstrsolid%contacts(i)%group
        if( .not. fstr_iscontactactive( fstrsolid, grpid, cstep ) ) cycle
        nlag = fstr_get_num_lagrange_pernode(fstrsolid%contacts(i)%algtype)
        do j = 1, size(fstrsolid%contacts(i)%slave)
          if( fstrsolid%contacts(i)%states(j)%state == contactfree ) cycle
          num_lagrange = num_lagrange + nlag
        enddo
      enddo
 
      do i = 1, fstrsolid%n_embeds
        grpid = fstrsolid%embeds(i)%group
        if( .not. fstr_iscontactactive( fstrsolid, grpid, cstep ) ) cycle
        nlag = 3
        do j = 1, size(fstrsolid%embeds(i)%slave)
          if( fstrsolid%embeds(i)%states(j)%state == contactfree ) cycle
          num_lagrange = num_lagrange + nlag
        enddo
      enddo
  endif
 
    ! Get original list of related nodes
    call hecmw_init_noderelated_from_org(hecmat%NP,num_lagrange,is_contact_active,list_noderelated_org,list_noderelated)
 
    ! Construct new list of related nodes and Lagrange multipliers
    countnon0lu_node = npl_org + npu_org
    countnon0lu_lagrange = 0
    if( is_contact_active ) call getnewlistofrelatednodesandlagrangemultipliers(cstep,contact_algo, &
       &  hecmat%NP,fstrsolid,countnon0lu_node,countnon0lu_lagrange,list_noderelated)
 
    ! Construct new matrix structure(hecMAT&hecLagMAT)
    numnon0_node = countnon0lu_node/2
    numnon0_lagrange = countnon0lu_lagrange/2
    call hecmw_construct_hecmat_from_noderelated(hecmat%N, hecmat%NP, hecmat%NDOF, &
      & numnon0_node, num_lagrange, list_noderelated, hecmat)
    call hecmw_construct_hecmat_from_noderelated(hecmat%N, hecmat%NP, hecmat%NDOF, &
      & numnon0_node, num_lagrange, list_noderelated, conmat)
    if( contact_algo == kcaslagrange ) call hecmw_construct_heclagmat_from_noderelated(hecmat%NP, &
      & hecmat%NDOF, num_lagrange, numnon0_lagrange, is_contact_active, list_noderelated, heclagmat)
    call hecmw_finalize_noderelated(list_noderelated)
 
    ! Copy Lagrange multipliers
    if( is_contact_active .and. contact_algo == kcaslagrange ) &
      call fstr_copy_lagrange_contact(fstrsolid,heclagmat)
 
  end subroutine fstr_mat_con_contact
 
  subroutine getnewlistofrelatednodesandlagrangemultipliers( &
      & cstep, contact_algo, np, fstrSOLID, countNon0LU_node, countNon0LU_lagrange, list_nodeRelated )
    integer(kind=kint),intent(in)             :: cstep
    integer(kind=kint),intent(in)             :: contact_algo
    integer(kind=kint),intent(in)             :: np
    type(fstr_solid),intent(in)               :: fstrsolid
    integer(kind=kint), intent(inout)         :: countnon0lu_node, countnon0lu_lagrange
    type(noderelated), pointer, intent(inout) :: list_noderelated(:)
 
    integer(kind=kint)            :: grpid
    integer(kind=kint)            :: count_lagrange
    integer(kind=kint)            :: ctsurf, etype, nnode, ndlocal(l_max_surface_node + 1)
    integer(kind=kint)            :: i, j, k, nlag, algtype
    real(kind=kreal)              :: fcoeff 
    logical                       :: necessary_to_insert_node
 
    count_lagrange = 0
    do i = 1, fstrsolid%n_contacts
 
      grpid = fstrsolid%contacts(i)%group
      if( .not. fstr_iscontactactive( fstrsolid, grpid, cstep ) ) cycle
 
      fcoeff = fstrsolid%contacts(i)%fcoeff
      necessary_to_insert_node = ( fcoeff /= 0.0d0 .or. contact_algo == kcaalagrange )
 
      algtype = fstrsolid%contacts(i)%algtype
      nlag = fstr_get_num_lagrange_pernode(algtype)
      if( contact_algo == kcaalagrange ) nlag = 1
      if( algtype == contacttied ) permission = .true.
 
      do j = 1, size(fstrsolid%contacts(i)%slave)
 
        if( fstrsolid%contacts(i)%states(j)%state == contactfree ) cycle
        ctsurf = fstrsolid%contacts(i)%states(j)%surface
        etype = fstrsolid%contacts(i)%master(ctsurf)%etype
        if( etype/=fe_tri3n .and. etype/=fe_quad4n ) &
          stop " ##Error: This element type is not supported in contact analysis !!! "
        nnode = size(fstrsolid%contacts(i)%master(ctsurf)%nodes)
        ndlocal(1) = fstrsolid%contacts(i)%slave(j)
        ndlocal(2:nnode+1) = fstrsolid%contacts(i)%master(ctsurf)%nodes(1:nnode)
 
        do k=1,nlag
          if( contact_algo == kcaslagrange ) count_lagrange = count_lagrange + 1
          call hecmw_ass_noderelated_from_contact_pair(np, nnode, ndlocal, count_lagrange, permission, &
            & necessary_to_insert_node, list_noderelated_org, list_noderelated, countnon0lu_node, countnon0lu_lagrange )
        enddo
      enddo
 
    enddo
 
    do i = 1, fstrsolid%n_embeds
 
      grpid = fstrsolid%embeds(i)%group
      if( .not. fstr_isembedactive( fstrsolid, grpid, cstep ) ) cycle
 
      necessary_to_insert_node = ( contact_algo == kcaalagrange )
 
      nlag = 3
      if( contact_algo == kcaalagrange ) nlag = 1
      permission = .true.
 
      do j = 1, size(fstrsolid%embeds(i)%slave)
 
        if( fstrsolid%embeds(i)%states(j)%state == contactfree ) cycle
        ctsurf = fstrsolid%embeds(i)%states(j)%surface
        etype = fstrsolid%embeds(i)%master(ctsurf)%etype
        nnode = size(fstrsolid%embeds(i)%master(ctsurf)%nodes)
        ndlocal(1) = fstrsolid%embeds(i)%slave(j)
        ndlocal(2:nnode+1) = fstrsolid%embeds(i)%master(ctsurf)%nodes(1:nnode)
 
        do k=1,nlag
          if( contact_algo == kcaslagrange ) count_lagrange = count_lagrange + 1
          call hecmw_ass_noderelated_from_contact_pair(np, nnode, ndlocal, count_lagrange, permission, &
            & necessary_to_insert_node, list_noderelated_org, list_noderelated, countnon0lu_node, countnon0lu_lagrange )
        enddo
      enddo
 
    enddo
 
  end subroutine getnewlistofrelatednodesandlagrangemultipliers
 
  subroutine fstr_copy_lagrange_contact(fstrSOLID,hecLagMAT)
 
    type(fstr_solid)                        :: fstrsolid
    type(hecmwst_matrix_lagrange)          :: heclagmat
    integer (kind=kint)                    :: id_lagrange, algtype, i, j, k, nlag
 
    id_lagrange = 0
 
    do i = 1, fstrsolid%n_contacts
 
      algtype = fstrsolid%contacts(i)%algtype
      nlag = fstr_get_num_lagrange_pernode(algtype)
 
      do j = 1, size(fstrsolid%contacts(i)%slave)
        if( fstrsolid%contacts(i)%states(j)%state == contactfree ) cycle
        do k=1,nlag
          id_lagrange = id_lagrange + 1
          heclagmat%Lagrange(id_lagrange)=fstrsolid%contacts(i)%states(j)%multiplier(k)
        enddo
      enddo
    enddo
 
    do i = 1, fstrsolid%n_embeds
      nlag = 3
      do j = 1, size(fstrsolid%embeds(i)%slave)
        if( fstrsolid%embeds(i)%states(j)%state == contactfree ) cycle
        do k=1,nlag
          id_lagrange = id_lagrange + 1
          heclagmat%Lagrange(id_lagrange)=fstrsolid%embeds(i)%states(j)%multiplier(k)
        enddo
      enddo
    enddo
 
  end subroutine fstr_copy_lagrange_contact
 
  logical function fstr_is_matrixstruct_symmetric(fstrSOLID,hecMESH)
 
    type(fstr_solid )        :: fstrsolid
    type(hecmwst_local_mesh) :: hecmesh
    integer (kind=kint)      :: is_in_contact
 
    is_in_contact = 0
    if( fstrsolid%n_contacts>0 ) then
      if( any(fstrsolid%contacts(:)%fcoeff /= 0.0d0) )  is_in_contact = 1
    endif
    call hecmw_allreduce_i1(hecmesh, is_in_contact, hecmw_max)
    fstr_is_matrixstruct_symmetric = (is_in_contact == 0)
 
  end function fstr_is_matrixstruct_symmetric
 
 
end module fstr_matrix_con_contact