fstr_dynamic_setup.f90

Go to the documentation of this file.

!-------------------------------------------------------------------------------
! Copyright (c) 2019 FrontISTR Commons
! This software is released under the MIT License, see LICENSE.txt
!-------------------------------------------------------------------------------
module m_fstr_dynamic_setup
  use m_fstr
  use hecmw_setup_util, only &
    : node_grp_name_to_id_ex
  use fstr_setup_util, only &
    : fstr_param_pack &
    , fstr_ctrl_err_stop
  use fstr_ctrl_dynamic
  implicit none
 
contains
 
  subroutine fstr_dynamic_init( fstrDYNAMIC )
    use m_fstr
    type(fstr_dynamic) :: fstrDYNAMIC
    fstrdynamic%idx_eqa  = 1
    fstrdynamic%idx_resp = 1
    fstrdynamic%n_step   = 1
    fstrdynamic%t_start  = 0.0
    fstrdynamic%t_curr   = 0.0d0
    fstrdynamic%t_end    = 1.0
    fstrdynamic%t_delta  = 1.0
    fstrdynamic%gamma    = 0.5
    fstrdynamic%beta     = 0.25
    fstrdynamic%alpha    = 0.0
    fstrdynamic%idx_mas  = kmasslumped
    fstrdynamic%idx_dmp  = 1
    fstrdynamic%ray_m    = 0.0
    fstrdynamic%ray_k    = 0.0
    fstrdynamic%restart_nout = 0
    fstrdynamic%nout         = 100
    fstrdynamic%ngrp_monit   = 0
    fstrdynamic%nout_monit   = 1
    fstrdynamic%iout_list(1) = 0
    fstrdynamic%iout_list(2) = 0
    fstrdynamic%iout_list(3) = 0
    fstrdynamic%iout_list(4) = 0
    fstrdynamic%iout_list(5) = 0
    fstrdynamic%iout_list(6) = 0
 
  end subroutine fstr_dynamic_init
 
  subroutine fstr_dynamic_alloc( hecMESH, fstrDYNAMIC )
    use m_fstr
    type(hecmwst_local_mesh),target :: hecMESH
    type(fstr_dynamic) :: fstrDYNAMIC
 
    integer :: ierror, ndof,nnod
 
    ndof=hecmesh%n_dof
    nnod=hecmesh%n_node
    if(fstrdynamic%idx_eqa == 11) then
      allocate( fstrdynamic%DISP(ndof*nnod,3)  ,stat=ierror )
      if( ierror /= 0 ) then
        write(idbg,*) 'stop due to allocation error <fstr_solve_LINEAR_DYNAMIC, DISP>'
        write(idbg,*) '  rank = ', hecmesh%my_rank,'  ierror = ',ierror
        call flush(idbg)
        call hecmw_abort( hecmw_comm_get_comm())
      end if
      allocate( fstrdynamic%VEL (ndof*nnod,1)  ,stat=ierror )
      if( ierror /= 0 ) then
        write(idbg,*) 'stop due to allocation error <fstr_solve_LINEAR_DYNAMIC, VEL>'
        write(idbg,*) '  rank = ', hecmesh%my_rank,'  ierror = ',ierror
        call flush(idbg)
        call hecmw_abort( hecmw_comm_get_comm())
      end if
      allocate( fstrdynamic%ACC (ndof*nnod,1)  ,stat=ierror )
      if( ierror /= 0 ) then
        write(idbg,*) 'stop due to allocation error <fstr_solve_LINEAR_DYNAMIC, ACC>'
        write(idbg,*) '  rank = ', hecmesh%my_rank,'  ierror = ',ierror
        call flush(idbg)
        call hecmw_abort( hecmw_comm_get_comm())
      end if
    else
      allocate( fstrdynamic%DISP(ndof*nnod,2)  ,stat=ierror )
      if( ierror /= 0 ) then
        write(idbg,*) 'stop due to allocation error <fstr_solve_LINEAR_DYNAMIC, DISP>'
        write(idbg,*) '  rank = ', hecmesh%my_rank,'  ierror = ',ierror
        call flush(idbg)
        call hecmw_abort( hecmw_comm_get_comm())
      end if
      allocate( fstrdynamic%VEL (ndof*nnod,2)  ,stat=ierror )
      if( ierror /= 0 ) then
        write(idbg,*) 'stop due to allocation error <fstr_solve_LINEAR_DYNAMIC, VEL>'
        write(idbg,*) '  rank = ', hecmesh%my_rank,'  ierror = ',ierror
        call flush(idbg)
        call hecmw_abort( hecmw_comm_get_comm())
      end if
      allocate( fstrdynamic%ACC (ndof*nnod,2)  ,stat=ierror )
      if( ierror /= 0 ) then
        write(idbg,*) 'stop due to allocation error <fstr_solve_LINEAR_DYNAMIC, ACC>'
        write(idbg,*) '  rank = ', hecmesh%my_rank,'  ierror = ',ierror
        call flush(idbg)
        call hecmw_abort( hecmw_comm_get_comm())
      end if
    endif
 
 
    allocate( fstrdynamic%VEC1(ndof*nnod)    ,stat=ierror )
    if( ierror /= 0 ) then
      write(idbg,*) 'stop due to allocation error <fstr_solve_LINEAR_DYNAMIC, VEC1>'
      write(idbg,*) '  rank = ', hecmesh%my_rank,'  ierror = ',ierror
      call flush(idbg)
      call hecmw_abort( hecmw_comm_get_comm())
    end if
    allocate( fstrdynamic%VEC2(ndof*nnod)    ,stat=ierror )
    if( ierror /= 0 ) then
      write(idbg,*) 'stop due to allocation error <fstr_solve_LINEAR_DYNAMIC, VEC2>'
      write(idbg,*) '  rank = ', hecmesh%my_rank,'  ierror = ',ierror
      call flush(idbg)
      call hecmw_abort( hecmw_comm_get_comm())
    end if
    allocate( fstrdynamic%VEC3(ndof*nnod)    ,stat=ierror )
    if( ierror /= 0 ) then
      write(idbg,*) 'stop due to allocation error <fstr_solve_LINEAR_DYNAMIC, VEC3>'
      write(idbg,*) '  rank = ', hecmesh%my_rank,'  ierror = ',ierror
      call flush(idbg)
      call hecmw_abort( hecmw_comm_get_comm())
    end if
    if( fstrdynamic%idx_eqa == 2 ) then
      allocate( fstrdynamic%HHT_FORCE(ndof*nnod)    ,stat=ierror )
      if( ierror /= 0 ) then
        write(idbg,*) 'stop due to allocation error <fstr_solve_LINEAR_DYNAMIC, HHT_FORCE>'
        write(idbg,*) '  rank = ', hecmesh%my_rank,'  ierror = ',ierror
        call flush(idbg)
        call hecmw_abort( hecmw_comm_get_comm())
      end if
      fstrdynamic%HHT_FORCE(:) = 0.0d0
    endif
 
  end subroutine fstr_dynamic_alloc
 
  subroutine fstr_dynamic_finalize( fstrDYNAMIC )
    type(fstr_dynamic) :: fstrDYNAMIC
 
    integer :: ierror
    if( associated(fstrdynamic%DISP) )        &
      deallocate( fstrdynamic%DISP    ,stat=ierror )
    if( ierror /= 0 ) then
      write(idbg,*) 'stop due to deallocation error <fstr_solve_LINEAR_DYNAMIC, DISP>'
      call flush(idbg)
      call hecmw_abort( hecmw_comm_get_comm())
    end if
    if( associated(fstrdynamic%VEL) )        &
      deallocate( fstrdynamic%VEL     ,stat=ierror )
    if( ierror /= 0 ) then
      write(idbg,*) 'stop due to deallocation error <fstr_solve_LINEAR_DYNAMIC, VEL>'
      call flush(idbg)
      call hecmw_abort( hecmw_comm_get_comm())
    end if
    if( associated(fstrdynamic%ACC) )        &
      deallocate( fstrdynamic%ACC     ,stat=ierror )
    if( ierror /= 0 ) then
      write(idbg,*) 'stop due to deallocation error <fstr_solve_LINEAR_DYNAMIC, ACC>'
      call flush(idbg)
      call hecmw_abort( hecmw_comm_get_comm())
    end if
    if( associated(fstrdynamic%HHT_FORCE) ) then
      deallocate( fstrdynamic%HHT_FORCE    ,stat=ierror )
      if( ierror /= 0 ) then
        write(idbg,*) 'stop due to deallocation error <fstr_solve_LINEAR_DYNAMIC, HHT_FORCE>'
        call flush(idbg)
        call hecmw_abort( hecmw_comm_get_comm())
      end if
    endif
    if( associated(fstrdynamic%VEC1) )        &
      deallocate( fstrdynamic%VEC1    ,stat=ierror )
    if( ierror /= 0 ) then
      write(idbg,*) 'stop due to deallocation error <fstr_solve_LINEAR_DYNAMIC, VEC1>'
      call flush(idbg)
      call hecmw_abort( hecmw_comm_get_comm())
    end if
    if( associated(fstrdynamic%VEC2) )        &
      deallocate( fstrdynamic%VEC2    ,stat=ierror )
    if( ierror /= 0 ) then
      write(idbg,*) 'stop due to deallocation error <fstr_solve_LINEAR_DYNAMIC, VEC2>'
      call flush(idbg)
      call hecmw_abort( hecmw_comm_get_comm())
    end if
    if( associated(fstrdynamic%VEC3) )        &
      deallocate( fstrdynamic%VEC3    ,stat=ierror )
    if( ierror /= 0 ) then
      write(idbg,*) 'stop due to deallocation error <fstr_solve_LINEAR_DYNAMIC, VEC3>'
      call flush(idbg)
      call hecmw_abort( hecmw_comm_get_comm())
    end if
 
  end subroutine
 
  !-----------------------------------------------------------------------------!
  !-----------------------------------------------------------------------------!
 
  subroutine fstr_setup_dynamic( ctrl, P )
    implicit none
    integer(kind=kint) :: ctrl
    type(fstr_param_pack) :: P
    integer(kind=kint) :: rcode
    character(HECMW_NAME_LEN) :: grp_id_name(1)
    integer(kind=kint) :: grp_id(1)
    logical :: is_alpha_given
 
    rcode = fstr_ctrl_get_dynamic( ctrl, &
      p%PARAM%nlgeom,  &
      p%DYN%idx_eqa, &
      p%DYN%idx_resp,&
      p%DYN%n_step,  &
      p%DYN%t_start, &
      p%DYN%t_end,   &
      p%DYN%t_delta, &
      p%DYN%gamma,   &
      p%DYN%beta,    &
      p%DYN%alpha,   &
      p%DYN%idx_mas, &
      p%DYN%idx_dmp, &
      p%DYN%ray_m,   &
      p%DYN%ray_k,   &
      p%DYN%nout,    &
      grp_id_name(1), hecmw_name_len,  &
      p%DYN%nout_monit,  &
      p%DYN%iout_list )
 
    if( rcode /= 0) call fstr_ctrl_err_stop
 
    if( p%DYN%idx_eqa /= 1 .and. p%DYN%idx_eqa /= 2 .and. p%DYN%idx_eqa /= 11 ) then
      write(*,*) '### Error: !DYNAMIC time integration type must be 1 (Newmark-beta), 2 (HHT) or 11 (explicit)'
      write(ilog,*) '### Error: !DYNAMIC time integration type must be 1 (Newmark-beta), 2 (HHT) or 11 (explicit)'
      call flush(6)
      call flush(ilog)
      call fstr_ctrl_err_stop
    endif
    if( p%DYN%idx_eqa == 2 .and. p%DYN%idx_resp /= 1 ) then
      write(*,*) '### Error: HHT method is available only for time history analysis'
      write(ilog,*) '### Error: HHT method is available only for time history analysis'
      call flush(6)
      call flush(ilog)
      call fstr_ctrl_err_stop
    endif
 
    is_alpha_given = ( p%DYN%alpha < huge(0.0d0) )
    if( p%DYN%idx_eqa == 2 ) then
      if( .not. is_alpha_given ) then
        write(*,*) '### Warning: !DYNAMIC HHT parameter alpha is not given; Newmark-beta method is used'
        write(ilog,*) '### Warning: !DYNAMIC HHT parameter alpha is not given; Newmark-beta method is used'
        p%DYN%idx_eqa = 1
        p%DYN%alpha = 0.0d0
      else
        if( p%DYN%alpha < -1.0d0/3.0d0 .or. p%DYN%alpha > 0.0d0 ) then
          write(*,*) '### Warning: !DYNAMIC HHT parameter alpha must be in [-1/3, 0]; alpha = 0 is used'
          write(ilog,*) '### Warning: !DYNAMIC HHT parameter alpha must be in [-1/3, 0]; alpha = 0 is used'
          p%DYN%alpha = 0.0d0
        endif
        ! gamma and beta given with alpha are not used
        p%DYN%beta  = 0.25d0*(1.0d0-p%DYN%alpha)**2
        p%DYN%gamma = 0.5d0-p%DYN%alpha
        if( myrank == 0 ) write(imsg,'(a,3(1x,1pe12.5))') &
          'INFO: HHT method with alpha, beta, gamma =', p%DYN%alpha, p%DYN%beta, p%DYN%gamma
      endif
    else
      if( is_alpha_given ) then
        write(*,*) '### Warning: !DYNAMIC HHT parameter alpha is ignored except for the HHT method (2)'
        write(ilog,*) '### Warning: !DYNAMIC HHT parameter alpha is ignored except for the HHT method (2)'
      endif
      p%DYN%alpha = 0.0d0
    endif
 
    if( p%DYN%idx_mas /= kmasslumped .and. p%DYN%idx_mas /= kmassconsistent ) then
      write(*,*) '### Error: !DYNAMIC mass matrix type must be 1 (lumped) or 2 (consistent)'
      write(ilog,*) '### Error: !DYNAMIC mass matrix type must be 1 (lumped) or 2 (consistent)'
      call flush(6)
      call flush(ilog)
      call fstr_ctrl_err_stop
    endif
    if( p%DYN%idx_mas == kmassconsistent .and. &
        (p%DYN%idx_eqa == 11 .or. p%DYN%idx_resp /= 1) ) then
      write(*,*) '### Error: consistent mass matrix is available only for implicit transient analysis'
      write(ilog,*) '### Error: consistent mass matrix is available only for implicit transient analysis'
      call flush(6)
      call flush(ilog)
      call fstr_ctrl_err_stop
    endif
 
    if (p%DYN%idx_resp == 1) then
      call node_grp_name_to_id_ex( p%MESH, '!DYNAMIC', 1, grp_id_name, grp_id)
      p%DYN%ngrp_monit = grp_id(1)
    else
      read(grp_id_name,*) p%DYN%ngrp_monit
    endif
 
  end subroutine fstr_setup_dynamic
 
end module m_fstr_dynamic_setup