Class md_status_class
In: LatticeClass/md_status_class.F90
constants_module comlib lattice_class error_class md_status_class dot/f_246.png

Defines Molecular Dynamics algorithm/controller status class

Version

$Id: md_status_class.F90,v 1.21 2011/09/25 14:04:23 ishikawa Exp $

Methods

Included Modules

constants_module comlib lattice_class error_class

Public Instance methods

MD_PQP_MTS_OMELYAN
Constant :
MD_PQP_MTS_OMELYAN = -3 :integer, parameter, public
MD_PQP_MTS_SIMPLE
Constant :
MD_PQP_MTS_SIMPLE = -1 :integer, parameter, public
MD_PQP_STS_IMPROVED
Constant :
MD_PQP_STS_IMPROVED = -2 :integer, parameter, public
MD_QPQ_MTS_OMELYAN
Constant :
MD_QPQ_MTS_OMELYAN = +3 :integer, parameter, public
MD_QPQ_MTS_SIMPLE
Constant :
MD_QPQ_MTS_SIMPLE = +1 :integer, parameter, public
MD_QPQ_STS_IMPROVED
Constant :
MD_QPQ_STS_IMPROVED = +2 :integer, parameter, public
Subroutine :
this :type(md_status), intent(inout)

delete MD controller

[Source]

subroutine delete_md_status(this)
!
! delete MD controller
!
  implicit none
  type(md_status), intent(inout) :: this
  if (associated(this%state)) deallocate(this%state)
  if (associated(this%nmd))   deallocate(this%nmd)
  if (associated(this%tau_p)) deallocate(this%tau_p)
  return
end subroutine
Function :
depth :integer
this :type(md_status), intent(in)

return MD step depth for current step.

[Source]

function get_depth_md_status(this) result(depth)
!
! return MD step depth for current step.
!
  implicit none
  type(md_status), intent(in) :: this
  integer :: depth
  integer :: i
  i = this%current_step
  depth = this%state(i)%depth
  return
end function
Function :
dt :real(DP)
this :type(md_status), intent(in)

return MD step size dt for current step.

[Source]

function get_dt_md_status(this) result(dt)
!
! return MD step size dt for current step.
!
  implicit none
  type(md_status), intent(in) :: this
  real(DP) :: dt
  integer :: i
  i = this%current_step
  if (1==this%reverse_time) then
    dt = -this%state(i)%dt
  else
    dt =  this%state(i)%dt
  endif
  return
end function
Function :
str :character(CHARLEN)
this :type(md_status), intent(in)

return MD controller parameters in formatted string.

[Source]

function get_parameter_string_md_status(this) result(str)
!
! return MD controller parameters in formatted string.
!
  implicit none
  type(md_status), intent(in) :: this
  character(CHARLEN) :: str
  integer :: i
  write(str,'(A," tau =",E22.15," NMD =",10I3)') TRIM(ADJUSTL(ALGNAME(this%algorithm))),this%tau,(this%nmd(i),i=this%max_depth,1,-1)
  return
end function
Function :
is_end :logical
this :type(md_status), intent(inout)

return MD controller loop status

.true.:MD required loop counts are reached. exit MD step loop.
.false.:MD step loop still remains. go next loop step.

[Source]

function is_do_loop_ended_md_status(this) result(is_end)
!
! return MD controller loop status
!
! .true.:: MD required loop counts are reached. exit MD step loop.
! .false.:: MD step loop still remains.  go next loop step.
!
  implicit none
  type(md_status), intent(inout) :: this
  logical :: is_end
  integer :: i
  i = this%current_step
  if (i == this%total_steps) then
    is_end = .TRUE.
    this%current_step = 0
  else
    is_end = .FALSE.
    this%current_step = this%current_step + 1
  endif
  return
end function
Function :
is_p :logical
this :type(md_status), intent(in)

return MD update type for current step.

.true.:user should update momentum (P)
.false.:user should update corrdinate (Q)

[Source]

function is_p_md_status(this) result(is_p)
!
! return MD update type for current step.
!
! .true.:: user should update momentum (P)
! .false.:: user should update corrdinate (Q)
!
  implicit none
  type(md_status), intent(in) :: this
  logical :: is_p
  integer :: i
  i = this%current_step
  is_p = this%state(i)%is_p
  return
end function
md_status
Derived Type :

MD controller status

Subroutine :
this :type(md_status), intent(inout)

initialize/create MD controller

[Source]

subroutine new_md_status(this)
!
! initialize/create MD controller
!
  implicit none
  type(md_status), intent(inout) :: this
  this%state => NULL()
  this%nmd   => NULL()
  this%tau_p => NULL()
  this%tau   = 0.0_DP
  this%tau_q = 0.0_DP
  this%algorithm    = 0
  this%max_depth    = 1
  this%step_count   = 0
  this%total_steps  = 0
  this%current_step = 0
  this%reverse_time = 0
  return
end subroutine
Subroutine :
this :type(md_status), intent(inout)

Print out MD controller parameters on display

[Source]

subroutine print_md_status(this)
!
! \Print out MD controller parameters on display
!
  use comlib
  implicit none
  type(md_status), intent(inout) :: this
  real(DP) :: dt
  integer :: i,iflg

  iflg=0
  if (nodeid==0) then
    select case(this%algorithm)
    case (MD_QPQ_MTS_SIMPLE)
      write(*,'(8X,"        MD Integ method : Multiple time scale Simple leapflog (QPQ order)")')
    case (MD_PQP_MTS_SIMPLE)
      write(*,'(8X,"        MD Integ method : Multiple time scale Simple leapflog (PQP order)")')
    case (MD_QPQ_STS_IMPROVED)
      write(*,'(8X,"        MD Integ method : Single time scale Improved leapflog (QPQ order)")')
      if (this%max_depth > 1) then
        write(*,'(8X,"        Error in max_detph !")')
        iflg=1
      endif
    case (MD_PQP_STS_IMPROVED)
      write(*,'(8X,"        MD Integ method : Single time scale Improved leapflog (PQP order)")')
      if (this%max_depth > 1) then
        write(*,'(8X,"        Error in max_detph !")')
        iflg=1
      endif
    case (MD_QPQ_MTS_OMELYAN)
      write(*,'(8X,"        MD Integ method : Multiple time scale Omelyan (QPQPQ order)")')
      write(*,'(8X,"                          with Lambda = ",F14.8)')OMELYAN_LAMBDA
    case (MD_PQP_MTS_OMELYAN)
      write(*,'(8X,"        MD Integ method : Multiple time scale Omelyan (PQPQP order)")')
      write(*,'(8X,"                          with Lambda = ",F14.8)')OMELYAN_LAMBDA
    end select
    write(*,*)
    write(*,'(6X,"        Trajectory Length :",E24.16)')this%tau
    write(*,'(6X,"         Time Scale Depth :",I3)')    this%max_depth
    dt = this%tau
    do i=this%max_depth,1,-1
      dt=dt/this%nmd(i)
      write(*,'(8X,"       MD step  dt(",I3,") :",E24.16)')i,dt
      write(*,'(8X,"  # of MD step Nmd(",I3,") :",I4)')    i,this%nmd(i)
    enddo
    write(*,*)
  endif

#ifndef _singlePU
  call comlib_bcast(iflg,0)
#endif

  if (iflg==1) then
    call error_stop("Error in Inputfile.")
  endif

  return
end subroutine
Subroutine :
this :type(md_status), intent(inout)
iout :integer, intent(in)

Read MD controller parameters from formatted file (by unit number)

[Source]

subroutine read_md_status(this,iout)
!
! \Read MD controller parameters from formatted file (by unit number)
!
  implicit none
  type(md_status), intent(inout) :: this
  integer, intent(in) :: iout
  integer :: i
  if (nodeid == 0) then
    read(iout,*)this%algorithm
    read(iout,*)this%tau
    read(iout,*)this%max_depth
  endif
#ifndef _singlePU
  call comlib_bcast(this%algorithm,0)
  call comlib_bcast(this%tau,0)
  call comlib_bcast(this%max_depth,0)
#endif
  allocate(this%nmd(1:this%max_depth))
  if (nodeid == 0) then
    read(iout,*)(this%nmd(i),i=this%max_depth,1,-1)
  endif
#ifndef _singlePU
  do i=1,this%max_depth
    call comlib_bcast(this%nmd(i),0)
  enddo
#endif

  call setup_md_status(this)

  return
end subroutine
Subroutine :
this :type(md_status), intent(inout)
iout :integer,intent(in)

Save MD controller parameters to formatted file (by unit number)

[Source]

subroutine save_md_status(this,iout)
!
! \Save MD controller parameters to formatted file (by unit number)
!
  implicit none
  type(md_status), intent(inout) :: this
  integer,intent(in) :: iout
  integer :: i
  if (nodeid == 0) then
    write(iout,*)         this%algorithm
    write(iout,'(F18.12)')this%tau
    write(iout,*)         this%max_depth
    write(iout,*)(this%nmd(i),i=this%max_depth,1,-1)
  endif
  return
end subroutine
Subroutine :
this :type(md_status), intent(inout)

set MD controller to do backward integration

[Source]

subroutine set_time_step_backward(this)
!
! set MD controller to do backward integration
!
  implicit none
  type(md_status), intent(inout) :: this
  this%reverse_time = 1
  return
end subroutine
Subroutine :
this :type(md_status), intent(inout)

set MD controller to do forward integration

[Source]

subroutine set_time_step_forward(this)
!
! set MD controller to do forward integration
!
  implicit none
  type(md_status), intent(inout) :: this
  this%reverse_time = 0
  return
end subroutine