Class sf_measure_class
In: MeasureClass/sf_measure_class.F90
lattice_class logfile_class quark_prop_class wavefunc_param_class wavefunc_class field_gauge_class quark_clover_class quark_wilson_class hmc_status_class comlib file_tools_class hmc_gluon_class error_class field_fermion_class sf_measure_class dot/f_81.png

This module contains the following SF measurements.

 - SF scheme coupling measurement.
 - PCAC correlation function measurement with SF boundary. (f_A,f_P,f_A',f_P')

Methods

delete   measure   new   sf_measurement  

Included Modules

lattice_class logfile_class quark_prop_class wavefunc_param_class wavefunc_class field_gauge_class quark_clover_class quark_wilson_class hmc_status_class comlib file_tools_class hmc_gluon_class error_class field_fermion_class

Public Instance methods

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

Delete SF measurement

[Source]

subroutine delete_sf_measurement(this)
!
! Delete SF measurement
!
  implicit none
  type(sf_measurement),      intent(inout) :: this
  this%status => NULL()
  call delete(this%logfile_gsf)
  call delete(this%logfile_pcac)
  call delete(this%logfile_solver)
  deallocate(this%logfile_gsf)
  deallocate(this%logfile_pcac)
  deallocate(this%logfile_solver)
  this%logfile_gsf => NULL()
  this%logfile_pcac => NULL()
  this%logfile_solver => NULL()
  if (0==nodeid) then
    write(*,'(80("="))')
  endif
  return
end subroutine
Subroutine :
this :type(sf_measurement), intent(inout)
gauge :type(gluon_rg_sf), intent(in)
u :type(vfield_gluon_wg), intent(in)

Measure SF scheme coupling and PCAC correlation functions

  this : SF measurement
 gauge : gauge action couplings (SF boundary couplings are used)
     u : gauge link

[Source]

subroutine measure_sf(this,gauge,u)
!
! Measure SF scheme coupling and PCAC correlation functions
!
!   this : SF measurement
!  gauge : gauge action couplings (SF boundary couplings are used)
!      u : gauge link
!
  use hmc_gluon_class, only : gluon_rg_sf
  implicit none
  type(sf_measurement),  intent(inout) :: this
  type(gluon_rg_sf),     intent(in) :: gauge
  type(vfield_gluon_wg), intent(in) :: u

  type(quark_prop), allocatable :: qp
  type(qprop_params), target :: params
  type(src_wavefunc_obj), target, allocatable :: wf



  !========================================
  ! compute SF scheme couplings
  !========================================
  call get_sf_scheme_coupling(this,gauge,u)

  !========================================
  ! compute SF PCAC corelation functions
  !========================================
  call new(params)
  call new(params%mass,id=0)
  call set_kappa(params%mass,this%kappa)
#ifdef _CLOVER
  call set_csw(params%mass,this%csw)
#endif
  params%tol = this%tolerance
  params%solver_log => this%logfile_solver
  params%smear_sous = 1
  params%smear_sink = 0

  allocate(qp,wf)
  call new(qp)
  call new(wf)
  qp%parameter     => params
  qp%wavefunc_sous => wf
  qp%wavefunc_sink => NULL()

  call set_sf_source(qp,u)
  call get_prop(qp,u)
  call get_pcac(this,qp)

  call delete(qp)
  call delete(wf)
  call delete(params)
  deallocate(qp,wf)

  return
end subroutine
Subroutine :
this :type(sf_measurement), intent(inout)
quark_fname :character(len=*), intent(in)
status :type(hmc_status), target, intent(in)

Initialize SF measurement

       this : SF measurement

quark_fname : SF PCAC quark mass parameter file name

     status : HMC job status

[Source]

subroutine new_sf_measurement(this,quark_fname,status)
!
! Initialize SF measurement
!
!        this : SF measurement
! quark_fname : SF PCAC quark mass parameter file name
!      status : HMC job status
!
  use comlib, only : comlib_bcast
  use file_tools_class
  implicit none
  type(sf_measurement),     intent(inout) :: this
  character(len=*),         intent(in)    :: quark_fname
  type(hmc_status), target, intent(in)    :: status
  character(len=CHARLEN) :: str
  integer :: iout,itraj

  this%fname_quark = TRIM(quark_fname)
  this%status => status

  if (0==nodeid) then
    iout = search_free_file_unit()
    open(iout,file=TRIM(this%fname_quark),status='old',form='formatted')
    read(iout,*)str
#ifdef _CLOVER
    read(iout,*)this%kappa,this%csw
#else
    this%csw = ZERO
    read(iout,*)this%kappa
#endif
    read(iout,*)this%nflavor
    read(iout,*)str
    read(iout,*)this%tolerance
    close(iout)
  endif
#ifndef _singlePU  
  call comlib_bcast(this%kappa,0)
  call comlib_bcast(this%csw,0)
  call comlib_bcast(this%nflavor,0)
  call comlib_bcast(this%tolerance,0)
#endif

  allocate(this%logfile_gsf)
  allocate(this%logfile_pcac)
  allocate(this%logfile_solver)

  !===============================================
  ! open SF scheme coupling log file
  !===============================================
  itraj = get_trajectory_number(this%status)
  write(str,'(A,I6.6)')TRIM(ADJUSTL(this%fname_gsf)),itraj
  call new(this%logfile_gsf,TRIM(str))

  call print(this%logfile_gsf,"# SF scheme coupling")
  write(str,'("# traj",8X,"1/(g_SF^2)",16X,"v-term",16X,"k/(g_SF^2)",15X,"k*v-term",20X,"k")')
  call print(this%logfile_gsf,TRIM(str))


  !===============================================
  ! open SF PCAC log file
  !===============================================
  write(str,'(A,I6.6)')TRIM(ADJUSTL(this%fname_pcac)),itraj
  call new(this%logfile_pcac,TRIM(str))

  call print(this%logfile_pcac,"# SF PCAC correlation functions")
  write(str,'("# traj itt",9X,"fA(itt)",17X,"fP(itt)",17X,"fA''(itt)",16X,"fP''(itt)")')
  call print(this%logfile_pcac,TRIM(str))


  !===============================================
  ! open SF PCAC quark solver log file
  !===============================================
  write(str,'(A,I6.6)')TRIM(ADJUSTL(this%fname_solver)),itraj
  call new(this%logfile_solver,TRIM(str))
  call print(this%logfile_solver,"# SF PCAC correlation solver log")


  if (0==nodeid) then
    write(*,'(80("="))')
    write(*,'("  SF measuement : SF coupling and PCAC")')
    write(*,'("        Trajectory :",I10)')  itraj
    write(*,'("       PCAC: Kappa :",F16.9)')this%kappa
#ifdef _CLOVER
    write(*,'("       PCAC:   Csw :",F16.9)')this%csw
#endif
    write(*,'("   SF coupling: Nf :",I3)')   this%nflavor
    write(*,'("  Solver tolerance :",E24.15)')this%tolerance
    write(*,'(80("-"))')
  endif

  return
end subroutine
sf_measurement
Derived Type :
fname_quark = "" :character(len=CHARLEN)
fname_gsf = "Sf_coupling." :character(len=CHARLEN)
fname_pcac = "Sf_PCAC." :character(len=CHARLEN)
fname_solver = "Sf_PCAC_solver_log." :character(len=CHARLEN)
logfile_gsf => NULL() :type(logfile), pointer
logfile_pcac => NULL() :type(logfile), pointer
logfile_solver => NULL() :type(logfile), pointer
status => NULL() :type(hmc_status), pointer
kappa :real(DP)
csw :real(DP)
nflavor :integer
tolerance :real(DP)
fA0(2:NTT) :complex(DP)
fP0(2:NTT) :complex(DP)
fA1(2:NTT) :complex(DP)
fP1(2:NTT) :complex(DP)
gSF2inv :real(DP)
vSF :real(DP)

SF measurement data (Wilson/Clover)