FrontISTR  5.2.0
Large-scale structural analysis program with finit element method
mcontactdef Module Reference

This module manage the data structure for contact calculation. More...

Data Types

type  tcontact
 Structure to includes all info needed by contact calculation. More...
 
type  fstr_info_contactchange
 

Functions/Subroutines

subroutine fstr_write_contact (file, contact)
 Write out contact definition. More...
 
subroutine fstr_contact_finalize (contact)
 Finalizer. More...
 
logical function fstr_contact_check (contact, hecMESH)
 Check the consistency with given mesh of contact defintiion. More...
 
logical function fstr_contact_init (contact, hecMESH, myrank)
 Initializer of tContactState. More...
 
subroutine clear_contact_state (contact)
 Reset contact state all to free. More...
 
subroutine scan_contact_state (flag_ctAlgo, contact, currpos, currdisp, ndforce, infoCTChange, nodeID, elemID, is_init, active, mu, B)
 This subroutine update contact states, which include. More...
 
subroutine calcu_contact_force0 (contact, coord, disp, ddisp, fcoeff, mu, mut, B)
 This subroutine update contact condition as follows: More...
 
subroutine update_contact_multiplier (contact, coord, disp, ddisp, fcoeff, mu, mut, gnt, ctchanged)
 This subroutine update lagrangian multiplier and the distance between contacting nodes. More...
 
subroutine ass_contact_force (contact, coord, disp, B)
 This subroutine assemble contact force into contacing nodes. More...
 
subroutine set_contact_state_vector (contact, dt, relvel_vec, state_vec)
 This subroutine setup contact output nodal vectors. More...
 
subroutine update_contact_tangentforce (contact)
 
subroutine track_contact_position_exp (nslave, contact, currpos, currdisp, infoCTChange, nodeID, elemID)
 This subroutine tracks down next contact position after a finite slide. More...
 
subroutine scan_contact_state_exp (contact, currpos, currdisp, infoCTChange, nodeID, elemID, is_init, active)
 This subroutine update contact states, which include. More...
 

Variables

real(kind=kreal), parameter distclr_init = 1.d-6
 dist clearance for initial scan More...
 

Detailed Description

This module manage the data structure for contact calculation.

Contact calculation takes into act after calling the following three subrotuines provided in this module

  1. Reading contact definition with subroutine: fstr_ctrl_get_CONTACT
  2. Check its consistency with mesh definition: fstr_contact_check
  3. Initilizing the contact calculation : fstr_contact_init

Function/Subroutine Documentation

◆ ass_contact_force()

subroutine mcontactdef::ass_contact_force ( type( tcontact ), intent(in)  contact,
real(kind=kreal), dimension(:), intent(in)  coord,
real(kind=kreal), dimension(:), intent(in)  disp,
real(kind=kreal), dimension(:), intent(inout)  B 
)

This subroutine assemble contact force into contacing nodes.

Parameters
[in]contactcontact info
[in]coordmesh coordinate
[in]dispdisp till current now
[in,out]bnodal force residual

Definition at line 961 of file fstr_contact_def.F90.

Here is the call graph for this function:

◆ calcu_contact_force0()

subroutine mcontactdef::calcu_contact_force0 ( type( tcontact ), intent(inout)  contact,
real(kind=kreal), dimension(:), intent(in)  coord,
real(kind=kreal), dimension(:), intent(in)  disp,
real(kind=kreal), dimension(:), intent(in)  ddisp,
real(kind=kreal), intent(in)  fcoeff,
real(kind=kreal), intent(in)  mu,
real(kind=kreal), intent(in)  mut,
real(kind=kreal), dimension(:), intent(inout)  B 
)

This subroutine update contact condition as follows:

  1. Contact force from multipler and disp increment
  2. Update nodal force residual
    Parameters
    [in,out]contactcontact info
    [in]coordmesh coordinate
    [in]dispdisp till current step
    [in]ddispdisp till current substep
    [in]fcoefffrictional coeff
    [in]mutpenalty
    [in,out]bnodal force residual

Definition at line 802 of file fstr_contact_def.F90.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ clear_contact_state()

subroutine mcontactdef::clear_contact_state ( type(tcontact), intent(inout)  contact)

Reset contact state all to free.

Parameters
[in,out]contactcontact definition

Definition at line 286 of file fstr_contact_def.F90.

Here is the caller graph for this function:

◆ fstr_contact_check()

logical function mcontactdef::fstr_contact_check ( type(tcontact), intent(inout)  contact,
type(hecmwst_local_mesh), pointer  hecMESH 
)

Check the consistency with given mesh of contact defintiion.

Parameters
[in,out]contactcontact definition
hecmeshmesh definition

Definition at line 140 of file fstr_contact_def.F90.

◆ fstr_contact_finalize()

subroutine mcontactdef::fstr_contact_finalize ( type(tcontact), intent(inout)  contact)

Finalizer.

Parameters
[in,out]contactcontact definition

Definition at line 124 of file fstr_contact_def.F90.

Here is the call graph for this function:

◆ fstr_contact_init()

logical function mcontactdef::fstr_contact_init ( type(tcontact), intent(inout)  contact,
type(hecmwst_local_mesh), pointer  hecMESH,
integer(kint), intent(in), optional  myrank 
)

Initializer of tContactState.

Parameters
[in,out]contactcontact definition
hecmeshmesh definition

Definition at line 169 of file fstr_contact_def.F90.

Here is the call graph for this function:

◆ fstr_write_contact()

subroutine mcontactdef::fstr_write_contact ( integer(kind=kint), intent(in)  file,
type(tcontact), intent(in)  contact 
)

Write out contact definition.

Parameters
[in]filefile number
[in]contactcontact definition

Definition at line 102 of file fstr_contact_def.F90.

Here is the call graph for this function:

◆ scan_contact_state()

subroutine mcontactdef::scan_contact_state ( character(len=9), intent(in)  flag_ctAlgo,
type( tcontact ), intent(inout)  contact,
real(kind=kreal), dimension(:), intent(in)  currpos,
real(kind=kreal), dimension(:), intent(in)  currdisp,
real(kind=kreal), dimension(:), intent(in)  ndforce,
type( fstr_info_contactchange ), intent(inout)  infoCTChange,
integer(kind=kint), dimension(:), intent(in)  nodeID,
integer(kind=kint), dimension(:), intent(in)  elemID,
logical, intent(in)  is_init,
logical, intent(out)  active,
real(kind=kreal), intent(in)  mu,
real(kind=kreal), dimension(:), optional, target  B 
)

This subroutine update contact states, which include.

  1. Free to contact or contact to free state changes
  2. Clear lagrangian multipliers when free to contact
    Parameters
    [in]flag_ctalgocontact analysis algorithm flag
    [in,out]contactcontact info
    [in,out]infoctchangecontact change info
    [in]currposcurrent coordinate of each nodes
    [in]currdispcurrent displacement of each nodes
    [in]ndforcenodal force
    [in]nodeidglobal nodal ID, just for print out
    [in]elemidglobal elemental ID, just for print out
    [in]is_initwheather initial scan or not
    [out]activeif any in contact
    [in]mupenalty
    bnodal force residual

Definition at line 309 of file fstr_contact_def.F90.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ scan_contact_state_exp()

subroutine mcontactdef::scan_contact_state_exp ( type( tcontact ), intent(inout)  contact,
real(kind=kreal), dimension(:), intent(in)  currpos,
real(kind=kreal), dimension(:), intent(in)  currdisp,
type( fstr_info_contactchange ), intent(inout)  infoCTChange,
integer(kind=kint), dimension(:), intent(in)  nodeID,
integer(kind=kint), dimension(:), intent(in)  elemID,
logical, intent(in)  is_init,
logical, intent(out)  active 
)

This subroutine update contact states, which include.

  1. Free to contact or contact to free state changes
  2. Clear lagrangian multipliers when free to contact
    Parameters
    [in,out]contactcontact info
    [in,out]infoctchangecontact change info
    [in]currposcurrent coordinate of each nodes
    [in]currdispcurrent displacement of each nodes
    [in]nodeidglobal nodal ID, just for print out
    [in]elemidglobal elemental ID, just for print out
    [in]is_initwheather initial scan or not
    [out]activeif any in contact

Definition at line 1164 of file fstr_contact_def.F90.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ set_contact_state_vector()

subroutine mcontactdef::set_contact_state_vector ( type( tcontact ), intent(in)  contact,
real(kind=kreal), intent(in)  dt,
real(kind=kreal), dimension(:), intent(inout)  relvel_vec,
real(kind=kreal), dimension(:), intent(inout)  state_vec 
)

This subroutine setup contact output nodal vectors.

Parameters
[in]contactcontact info
[in,out]relvel_vecmesh coordinate
[in,out]state_vecdisp till current now

Definition at line 1015 of file fstr_contact_def.F90.

Here is the caller graph for this function:

◆ track_contact_position_exp()

subroutine mcontactdef::track_contact_position_exp ( integer, intent(in)  nslave,
type( tcontact ), intent(inout)  contact,
real(kind=kreal), dimension(:), intent(in)  currpos,
real(kind=kreal), dimension(:), intent(in)  currdisp,
type( fstr_info_contactchange ), intent(inout)  infoCTChange,
integer(kind=kint), dimension(:), intent(in)  nodeID,
integer(kind=kint), dimension(:), intent(in)  elemID 
)

This subroutine tracks down next contact position after a finite slide.

Parameters
[in]nslaveslave node
[in,out]contactcontact info
[in,out]infoctchangecontact change info
[in]currposcurrent coordinate of each nodes
[in]currdispcurrent displacement of each nodes
[in]nodeidglobal nodal ID, just for print out
[in]elemidglobal elemental ID, just for print out

checking the contact element of last step

Definition at line 1053 of file fstr_contact_def.F90.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ update_contact_multiplier()

subroutine mcontactdef::update_contact_multiplier ( type( tcontact ), intent(inout)  contact,
real(kind=kreal), dimension(:), intent(in)  coord,
real(kind=kreal), dimension(:), intent(in)  disp,
real(kind=kreal), dimension(:), intent(in)  ddisp,
real(kind=kreal), intent(in)  fcoeff,
real(kind=kreal), intent(in)  mu,
real(kind=kreal), intent(in)  mut,
real(kind=kreal), dimension(2), intent(out)  gnt,
logical, intent(inout)  ctchanged 
)

This subroutine update lagrangian multiplier and the distance between contacting nodes.

Parameters
[in,out]contactcontact info
[in]coordmesh coordinate
[in]dispdisp till current step
[in]ddispdisp till current substep
[in]fcoefffrictional coeff
[in]mutpenalty
[out]gntconvergency information
[in,out]ctchangedif contact state changes

Definition at line 877 of file fstr_contact_def.F90.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ update_contact_tangentforce()

subroutine mcontactdef::update_contact_tangentforce ( type( tcontact ), intent(inout)  contact)
Parameters
[in,out]contactcontact info

Definition at line 1035 of file fstr_contact_def.F90.

Here is the caller graph for this function:

Variable Documentation

◆ distclr_init

real(kind=kreal), parameter mcontactdef::distclr_init = 1.d-6

dist clearance for initial scan

Definition at line 70 of file fstr_contact_def.F90.