convert_ildg_config.F90

Path: ILDGClass/convert_ildg_config.F90
Last Update: Wed Jan 18 16:14:48 +0900 2012
dot/f_2.png

Required files

Methods

Included Modules

lattice_class file_tools_class ildg_class timer_class

Public Instance methods

Main Program :

This converts the ILG configuration binary data to the configuration binary data

 used in this HMC program for measurements.

Version

 $Id: convert_ildg_config.F90,v 1.5 2011/05/30 13:38:30 ishikawa Exp $

[Source]

program convert_ildg_config
!
!= This converts the ILG configuration binary data to the configuration binary data
!  used in this HMC program for measurements.
!
!== Version
!
!  $Id: convert_ildg_config.F90,v 1.5 2011/05/30 13:38:30 ishikawa Exp $
!
  use lattice_class
  use file_tools_class
  use ildg_class
  use timer_class
  implicit none

  integer, parameter :: NARGS = 10
  type(gauge_link_ildg) :: u

  character(len=1) :: cbyte(1)
  character(len=CHARLEN) :: fname,ftraj,fhead,fpath,comment,arg(NARGS)

  real(DP) :: plq,Kud,Ks,CSW,beta,cRG0
  integer :: traj, i
  type(timer) :: etime

  if ( COMMAND_ARGUMENT_COUNT() /= NARGS ) then
    write(*,'("Usage:")',advance='no')
    write(*,'(" convert_ildg_config <input file> <comment> <beta>")',advance='no')
    write(*,'(" <crg0> <Kud> <Ks> <Csw> <traj> <output dir> <output file name header>")')
    write(*,'("This converts the ILDG format data ordering to the parallel/even-odd data ordering.")')
    write(*,*)
    stop
  endif
  call new(etime,master_only=.true.)
  call tic(etime)
  do i=1,NARGS
    call GET_COMMAND_ARGUMENT(i,arg(i))
  enddo
  read(arg(1),'(A)')fname
  read(arg(2),'(A)')comment
  read(arg(3),*)beta
  read(arg(4),*)cRG0
  read(arg(5),*)Kud
  read(arg(6),*)Ks
  read(arg(7),*)Csw
  read(arg(8),*)traj
  read(arg(9), '(A)')fpath
  read(arg(10),'(A)')fhead

  write(*,'("INPUTFILE=",A)')TRIM(fname)
  write(*,'("  COMMENT=",A)')TRIM(comment)
  write(*,'("     Beta=",F16.7)')beta
  write(*,'("     cRG0=",F16.7)')cRG0
  write(*,'("      Kud=",F16.7)')Kud
  write(*,'("      Ks_=",F16.7)')Ks
  write(*,'("      Csw=",F16.7)')Csw
  write(*,'("     Traj=",I6)')traj
  write(*,'("    FPath=",A)')TRIM(fpath)
  write(*,'("    FNAME=",A)')TRIM(fhead)
  write(*,'("SIZEOF(U)=",I12)')u%sizeof

  call new(u)
  call read_ildg_config(TRIM(ADJUSTL(fname)),u,endian_conversion=.true.)

  plq = plaquette(u)

  write(*,'("PLQ=",F20.15)')plq

  write(ftraj,'(I6.6)')traj
  fpath = TRIM(fpath)//"/"//TRIM(ftraj)

  call convert_and_save_config(TRIM(ADJUSTL(fpath)),TRIM(ADJUSTL(fhead)),u,plq)

  call delete(u)
  call toc(etime)
  write(*,'("ETIME:",F10.4)')get_elapse(etime)

  stop
end program