Class metropolis_test_class
In: LatticeClass/metropolis_test_class.F90
comlib lattice_class hmc_logfile_class counter_class metropolis_test_class dot/f_240.png

Define Metropolis test class

Version

$Id: metropolis_test_class.F90,v 1.7 2011/06/28 14:00:11 ishikawa Exp $

Methods

Included Modules

comlib lattice_class hmc_logfile_class counter_class

Public Instance methods

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

delete Metropolis test status

[Source]

subroutine delete_mp(this)
!
! delete Metropolis test status
!
  implicit none
  type(metropolis_status), intent(inout) :: this
  call delete(this%mp_logfile)
  this%mp_logfname    = REPEAT(' ',CHARLEN)
  this%mp_test_name   = REPEAT(' ',CHARLEN)
  this%mp_weight_name = REPEAT(' ',CHARLEN)
  return
end subroutine
Function :
flag :logical
: .true. for accept, .false. for reject.
this :type(metropolis_status), intent(inout)
id :integer, intent(in)
: test ID (traj number or run number...)
dH :real(DP), intent(in)
: Metropolis test weight, dH = H1 - H0
mp_onoff :logical, intent(in)
: switch, .false. for always accept, .true. for do test.

return accept/reject result of Metropolis test

 Acceptance rate P is

 P = MIN[1, exp(-dH)]
.true.:accept
.false.:reject
  id : identifier in log (traj num or run num...)
  dH : weight
      mp_onoff : test_switch,  if true do mp test, if false always accept.

[Source]

function is_accepted_mp(this,id,dH,mp_onoff) result(flag)
!
! return accept/reject result of Metropolis test
!
!  Acceptance rate P is
!
!  P = MIN[1, exp(-dH)]
!
! .true.::  accept
! .false.:: reject
!
!             id : identifier in log (traj num or run num...)
!             dH : weight
!       mp_onoff : test_switch,  if true do mp test, if false always accept.
!
  implicit none
  type(metropolis_status), intent(inout) :: this
  integer,  intent(in) :: id        ! test ID (traj number or run number...)
  real(DP), intent(in) :: dH        ! Metropolis test weight, dH = H1 - H0
  logical,  intent(in) :: mp_onoff  ! switch, .false. for always accept, .true. for do test.
  logical :: flag                   ! .true. for accept, .false. for reject.
  integer :: ihit
  character(len=CHARLEN) :: str(5)
  real(DP) :: yr,th

  call get(g_rand,yr)

  !=========================================
  ! do Metropolis test only on Master node
  !=========================================
  if (nodeid == 0) then

    select case (mp_onoff)
    case (.FALSE.)
      yr = 0.0_DP
      th = 1.0_DP
    case (.TRUE.)
      th = EXP(-dH)
    end select

    str(:) = REPEAT(' ',LEN(str))
    if ( yr <= th ) then
      call inc(this%mp_counter,1)
      ihit = 1
      write(str(1),'("   Accepted : ",A," =",E24.16)')TRIM(this%mp_weight_name),dH
      write(str(2),'(" yr =< exp(-",A,") = ",E24.16," <=",E24.16)')TRIM(this%mp_weight_name),yr,th
    else
      call inc(this%mp_counter,0)
      ihit = 0
      write(str(1),'("   Rejected : ",A," =",E24.16)')TRIM(this%mp_weight_name),dH
      write(str(2),'(" yr  > exp(-",A,") = ",E24.16," > ",E24.16)')TRIM(this%mp_weight_name),yr,th
    endif

    write(str(3),'("@",I8,I2,"  ",A,":",E24.16)') id,ihit,TRIM(this%mp_weight_name),dH
    write(str(4),'("  TRY_=",I6," HIT_=",I6," Acceptance(%)_=",F8.4)') get_called(this%mp_counter), get_count(this%mp_counter), 100*get_average(this%mp_counter)
    write(str(5),*)
  endif

  call print(this%mp_logfile,str)

  !=========================================
  ! Bcast the MP test result to all nodes.
  !=========================================
#ifndef _singlePU
  call comlib_bcast(ihit,0)
#endif

  select case(ihit)
  case (1)
    flag = .TRUE.
  case (0)
    flag = .FALSE.
  end select

  return
end function
metropolis_status
Derived Type :

Metropolis test parameters

Subroutine :
this :type(metropolis_status), intent(inout)
fname :character(CHARLEN), intent(in)
mp_test_name :character(len=*), intent(in)
mp_weight_name :character(len=*), intent(in)

Initialize/create Metropolis test status

[Source]

subroutine new_mp(this,fname,mp_test_name,mp_weight_name)
!
! Initialize/create Metropolis test status
!
  implicit none
  type(metropolis_status), intent(inout) :: this
  character(CHARLEN), intent(in) :: fname
  character(len=*),   intent(in) :: mp_test_name
  character(len=*),   intent(in) :: mp_weight_name
  this%mp_logfname    = TRIM(ADJUSTL(fname))
  this%mp_test_name   = TRIM(ADJUSTL(mp_test_name))
  this%mp_weight_name = TRIM(ADJUSTL(mp_weight_name))
  call new(this%mp_counter)
  call new(this%mp_logfile,this%mp_logfname)
  return
end subroutine
Subroutine :
this :type(metropolis_status), intent(inout)
str :character(CHARLEN), intent(in)
: formatted string to be written to log file.

Print out Metropolis test status in formatted string to log file

[Source]

subroutine print_log_mp(this,str)
!
! \Print out Metropolis test status in formatted string to log file
!
  implicit none
  type(metropolis_status), intent(inout) :: this
  character(CHARLEN),      intent(in)    :: str  ! formatted string to be written to log file.
  if (nodeid==0) call print(this%mp_logfile,str)
  return
end subroutine
Subroutine :
this :type(metropolis_status), intent(inout)

Print out Metropolis test acceptance rate on display.

User should supply key name for this Metropolis test.

[Source]

subroutine print_stat_mp(this)
!
! \Print out Metropolis test acceptance rate on display.
!
! User should supply key name for this Metropolis test.
!
  implicit none
  type(metropolis_status), intent(inout) :: this
  if (nodeid==0) then
    write(*,'(1X,A," Acceptance (%) :",F12.4," (",I8,"/",I8,")")') TRIM(this%mp_test_name), 100*get_average(this%mp_counter), get_count(this%mp_counter), get_called(this%mp_counter)
  endif
  return
end subroutine