Class quark_prop_class
In: MeasureClass/quark_prop_class.F90
comlib lattice_class timer_class logfile_class field_gauge_class field_fermion_class quark_clover_class quark_wilson_class wavefunc_class error_class print_status_class hmc_std_quark_wilson_class hmc_blk_quark_wilson_class quark_prop_class dot/f_78.png

This defines quark propagtor for measurement program

 The quark parameter is defined as a pointer.
 The wavefunction for source and sink are defined as pointers.
 These pointers should be assigned or initialized before calling get_qprop.
 User should set the source (numerical) data on quark propagator before calling get_qprop.

 set_source is defined for typical wavefunction setting,
 with whith the numerical data are generated on a quark propagator according
 to the wavefunction source data.

Methods

assign   assign   assign_mult   assign_mult   clear_prop   delete   delete   get_prop   new   new   print   print   print_statistics   qprop_params   quark_prop   read   read   set_source  

Included Modules

comlib lattice_class timer_class logfile_class field_gauge_class field_fermion_class quark_clover_class quark_wilson_class wavefunc_class error_class print_status_class hmc_std_quark_wilson_class hmc_blk_quark_wilson_class

Public Instance methods

Subroutine :
param1 :type(qprop_params), intent(inout)
param2 :type(qprop_params), intent(in)

assignment parameters

param1 := param2

[Source]

subroutine assign_qprop_params(param1,param2)
!
! assignment parameters
!
! param1 := param2
!
  type(qprop_params), intent(inout) :: param1
  type(qprop_params), intent(in)    :: param2
  call assign(param1%mass,param2%mass)
  param1%solver_log  => param2%solver_log
  param1%fname       =  param2%fname
  param1%path        =  param2%path
  param1%smear_sous  =  param2%smear_sous
  param1%smear_sink  =  param2%smear_sink
  return
end subroutine
Subroutine :
qp1 :type(quark_prop), intent(inout)
qp2 :type(quark_prop), intent(in)

assign(=duplicate) quark propagator (incuding parameters and data)

q1 := q2

[Source]

subroutine assign_qprop(qp1,qp2)
!
! assign(=duplicate) quark propagator (incuding parameters and data)
!
! q1 := q2
!
  implicit none
  type(quark_prop), intent(inout) :: qp1
  type(quark_prop), intent(in)    :: qp2
  if (allocated(qp1%prop)) deallocate(qp1%prop)
  allocate(qp1%prop(COL,COL,SPIN,SPIN,NX*NY*NZ,NTT))
  qp1%prop = qp2%prop
  qp1%psum = qp2%psum
  qp1%parameter     => qp2%parameter
  qp1%wavefunc_sous => qp2%wavefunc_sous
  qp1%wavefunc_sink => qp2%wavefunc_sink
  return
end subroutine
Subroutine :
uq :type(quark_prop), intent(inout)
q :type(quark_prop), intent(in)
u :type(vfield_gluon_wg), intent(in)
mu :integer, intent(in)

assign a quark propagtor shifted and multiplied by gauge link

 uq(n) := u_mu(n) q(n+mu)

 uq : quark prop * link
  q : quark prop
  u : gauge link
 mu : direction

[Source]

subroutine assign_mult_u_q(uq,q,u,mu)
!  
! assign a quark propagtor shifted and multiplied by gauge link
!  
!  uq(n) := u_mu(n) q(n+mu) 
!  
!  uq : quark prop * link
!   q : quark prop
!   u : gauge link 
!  mu : direction
!  
  implicit none
  integer, intent(in) :: mu
  type(vfield_gluon_wg), intent(in)    :: u
  type(quark_prop),      intent(in)    :: q
  type(quark_prop),      intent(inout) :: uq
  type(field_quark_wg), allocatable :: y,uy
  integer :: ic,is

  call assign(uq,q)
  allocate(y,uy)
  call new(y)
  call new(uy)

  do is=1,SPIN
  do ic=1,COL
    call conv_prop_to_y(uq, y,ic,is)
    call assign_mult_u_y(uy,y,u,mu)
    call conv_y_to_prop(uy,uq,ic,is)
  enddo
  enddo

  deallocate(y,uy)

  return
end subroutine
Subroutine :
uy :type(field_quark_wg), intent(inout)
y :type(field_quark_wg), intent(inout)
u :type(vfield_gluon_wg), intent(in)
mu :integer, intent(in)

assign a fermion field shifted and multiplied by gauge link

 uy(n) := u_mu(n) y(n+mu)

 uy : quark field * link
  y : quark field
  u : gauge link
 mu : direction

[Source]

subroutine assign_mult_u_y(uy,y,u,mu)
!  
! assign a fermion field shifted and multiplied by gauge link
!  
!  uy(n) := u_mu(n) y(n+mu) 
!  
!  uy : quark field * link
!   y : quark field
!   u : gauge link 
!  mu : direction
!  
  implicit none
  type(field_quark_wg),  intent(inout) :: uy, y
  type(vfield_gluon_wg), intent(in)    :: u
  integer, intent(in) :: mu

  call assign_mult_u_eo(uy%eo(0),y%eo(1),u%eo(0),mu)
  call assign_mult_u_eo(uy%eo(1),y%eo(0),u%eo(1),mu)

  return
end subroutine
Subroutine :
this :type(quark_prop), intent(inout)
: $OMP PARALLEL WORKSHARE

clear quark propagator data

[Source]

subroutine clear_qprop(this)
!
! clear quark propagator data
!
  implicit none
  type(quark_prop), intent(inout) :: this
!$OMP PARALLEL WORKSHARE
  this%prop(:,:,:,:,:,:) = Z0
!$OMP END PARALLEL WORKSHARE
  return
end subroutine
Subroutine :
this :type(qprop_params), intent(inout)

delete quark parameter parameters

[Source]

subroutine delete_qprop_params(this)
!
! delete quark parameter parameters
!
  implicit none
  type(qprop_params), intent(inout) :: this
  this%solver_log => NULL()
  call delete(this%mass)
  return
end subroutine
Subroutine :
this :type(quark_prop), intent(inout)

set default values for a quark propagator

[Source]

subroutine delete_qprop(this)
!
! set default values for a quark propagator
!
  implicit none
  type(quark_prop), intent(inout) :: this
  if (allocated(this%prop)) deallocate(this%prop)
  this%parameter     => NULL()
  this%wavefunc_sous => NULL()
  this%wavefunc_sink => NULL()
  return
end subroutine
Subroutine :
this :type(quark_prop), intent(inout)
u :type(vfield_gluon_wg), intent(in)

compute quark propagator

[Source]

subroutine get_qprop(this,u)
!
! compute quark propagator
!
  use error_class
  implicit none
  type(quark_prop),      intent(inout) :: this
  type(vfield_gluon_wg), intent(in)    :: u
  type(timer) :: solver_time
#ifdef _USE_BLOCK_SOLVER
  type(field_quark_wg), allocatable :: y(:),Dy(:)
  integer :: ic,is,ics
  real(DP), allocatable :: rtmp0(:),rtmp1(:)
#else
  type(field_quark_wg), allocatable :: y,Dy
  integer :: ic,is
  real(DP) :: rtmp0,rtmp1
#endif
  character(len=CHARLEN) :: str
  integer :: maxiter,mode,iter
  integer :: iout
  real(DP) :: tol,etime
  real(DP) :: copy_y_time0,copy_y_time1,copy_y_etime
#ifdef _CLOVER
  type(tfield_gluon_wg), allocatable :: wl34,wl98
  real(DP) :: dummy
#endif

  call new(this%parameter%etime)
  call tic(this%parameter%etime)


#ifdef _CLOVER
!=======================
! compute clover term
!=======================
  allocate(wl34,wl98)
  call make_two_links(u,wl34,wl98)
  call make_clover_leaf(u,wl34,wl98)
  deallocate(wl34,wl98)
  call make_inverse_clover_term(this%parameter%mass,0,dummy)
  call delete_clover_leaf()
#endif

  if (.not.allocated(this%prop)) then
    call error_stop("get_qprop: Error, source is not setted/allocated.")
  endif

#ifdef _USE_BLOCK_SOLVER
!===============================
!  Blocked solver 
!===============================
  call new(solver_time)
  call tic(solver_time)

  copy_y_time0 = get_elapse(copy_fq_time)

  write(str,'("# BLOCK_SOLVER FOR COL=1,2,3 SPIN=1,2,3,4")')
  call print(this%parameter%solver_log,TRIM(str))
  iout = get_file_unit(this%parameter%solver_log)

  allocate(y(CLSP),Dy(CLSP))
  allocate(rtmp0(CLSP),rtmp1(CLSP))
  do ics=1,CLSP
    call new( y(ics))
    call new(dy(ics))
  enddo

  !------------------------
  ! Dy = D \ y
  !------------------------
  tol    = this%parameter%tol
  do is=1,SPIN
  do ic=1,COL
    ics = ic + (is-1)*COL
    call conv_prop_to_y(this, y(ics),ic,is)
    rtmp0(ics) = abs2(y(ics))
  enddo
  enddo

  call assign_inv_blk_mult_wd(this%parameter%mass,iout,tol,iter,Dy,y,u)

  do is=1,SPIN
  do ic=1,COL
    ics = ic + (is-1)*COL
    rtmp1(ics) = abs2(Dy(ics))
    call conv_y_to_prop(Dy(ics),this,ic,is)
  enddo
  enddo

  copy_y_time1 = get_elapse(copy_fq_time)
  copy_y_etime = copy_y_time1-copy_y_time0

  call toc(solver_time)
  etime = get_elapse(solver_time)
  if (nodeid==0) then
    write(*,'(" block solver done  ETime:",E24.16," COPY_Y_ETime:",E24.16)')etime,copy_y_etime
    do is=1,SPIN
    do ic=1,COL
      ics = ic + (is-1)*COL
      write(*,'(" (col,spin)=",2I3," |y|^2:",E24.16," |Dy|^2:",E24.16)')ic,is,rtmp0(ics),rtmp1(ics)
    enddo
    enddo
  endif
  deallocate(rtmp0,rtmp1)

#else

!===============================
!  Non-Blocked solver 
!===============================
  allocate(y,Dy)
  do is=1,SPIN
  do ic=1,COL

    call new(solver_time)
    call tic(solver_time)

    copy_y_time0 = get_elapse(copy_fq_time)

    write(str,'("# SOLVER FOR COL=",I2," SPIN=",I2)')ic,is
    call print(this%parameter%solver_log,TRIM(str))
    iout = get_file_unit(this%parameter%solver_log)

    call new(y)
    call new(Dy)

    !------------------------
    ! Dy = D \ y
    !------------------------
    tol    = this%parameter%tol
    call conv_prop_to_y(this, y,ic,is)
    rtmp0 = abs2(y)
    call assign_inv_mult_wd(this%parameter%mass,iout,tol,iter,Dy,y,u)
    rtmp1 = abs2(Dy)
    call conv_y_to_prop(Dy,this,ic,is)

    copy_y_time1 = get_elapse(copy_fq_time)
    copy_y_etime = copy_y_time1-copy_y_time0

    call toc(solver_time)
    etime = get_elapse(solver_time)
    if (nodeid==0) then
      write(*,'(" solver done ",2I3," ETime:",E24.16," COPY_Y_ETime:",E24.16," |y|^2:",E24.16," |Dy|^2:",E24.16)') ic,is,etime,copy_y_etime,rtmp0,rtmp1
    endif
  enddo
  enddo
#endif

  deallocate(y,Dy)
#ifdef _CLOVER
  call delete_inverse_clover_term(this%parameter%mass)
#endif

  call toc(this%parameter%etime)

  return
end subroutine
Subroutine :
this :type(qprop_params), intent(inout)

set default values for a quark propagator parameter

[Source]

subroutine new_qprop_params(this)
!
! set default values for a quark propagator parameter
!
  implicit none
  type(qprop_params), intent(inout) :: this

  this%solver_log  => NULL()
  this%tol = 10*EPSILON(1.0_DP)
  call new(this%etime)

  return
end subroutine
Subroutine :
this :type(quark_prop), intent(inout)

set default values for a quark propagator

[Source]

subroutine new_qprop(this)
!
! set default values for a quark propagator
!
  implicit none
  type(quark_prop), intent(inout) :: this
  if (allocated(this%prop)) deallocate(this%prop)
  this%parameter     => NULL()
  this%wavefunc_sous => NULL()
  this%wavefunc_sink => NULL()
  return
end subroutine
Subroutine :
this :type(qprop_params), intent(inout)

print parameters for a quark propagator

[Source]

subroutine print_qprop_params(this)
!
! print parameters for a quark propagator
!
  implicit none
  type(qprop_params), intent(inout) :: this

  if (nodeid==0) then
    write(*,'("==== Mass Params ====")')
    call print(this%mass)
  endif

  if (nodeid==0) then
    write(*,'("==== Solver Params ====")')
    write(*,'(3X," Solver tol      :",E24.16)')this%tol


    write(*,'("==== Smearing type ====")')
    write(*,'(3X,"          Source :",I3)')this%smear_sous
    write(*,'(3X,"            Sink :",I3)')this%smear_sink
  endif

  return
end subroutine
Subroutine :
this :type(quark_prop), intent(inout)

print parameters for a quark propagator

[Source]

subroutine print_qprop(this)
!
! print parameters for a quark propagator
!
  implicit none
  type(quark_prop), intent(inout) :: this

  call print(this%parameter)

  if (nodeid==0) then
    write(*,'("==== Source Wavefunction type ====")')
  endif
  call print(this%wavefunc_sous)

  if (0 < this%parameter%smear_sink) then
    if (nodeid==0) then
      write(*,'("==== Sink Wavefunction type ====")')
    endif
    call print(this%wavefunc_sink)
  endif

  if (nodeid==0) then
    write(*,'(80("="))')
  endif

  return
end subroutine
Subroutine :
this :type(quark_prop), intent(inout)

print out quark propagator statistics

[Source]

subroutine print_stat_qprop(this)
!
! print out quark propagator statistics
!
  implicit none
  type(quark_prop), intent(inout) :: this
  integer :: itt
  real(DP) :: etime

  call psum(this)

  etime = get_elapse(this%parameter%etime)
  if (nodeid==0) then
    write(*,'("==== Solver Statistics ====")')
    write(*,'(3X,"             Time :",F12.4)')etime
  endif

  if (nodeid==0) then
    write(*,'("==== Quark Prop |S(itt)|^2 ====")')
    do itt=1,NTT
      write(*,'(I3,E24.16)')itt,this%psum(itt)
    enddo
    write(*,'("")')
  endif

  return
end subroutine
qprop_params
Derived Type :
mass :type(quark_clover)
mass :type(quark_wilson)
mass :type(quark_clover)
mass :type(quark_wilson)
etime :type(timer)
solver_log => NULL() :type(logfile), pointer
fname :character(len=CHARLEN)
path :character(len=CHARLEN)
smear_sous :integer
: 0 for local, 1 for smeared
smear_sink :integer
: 0 for local, 1 for smeared
tol :real(DP)

quark propagator parameter for physics measurement

ifdef _CLOVER

quark_prop
Derived Type :
prop(:,:,:,:,:,:) :complex(DP), allocatable
: prop(sink,sous,sink,sous,sink,sink)
psum(NTT) :real(DP)
: check sum
parameter => NULL() :type(qprop_params), pointer
wavefunc_sous => NULL() :type(src_wavefunc_obj), pointer
wavefunc_sink => NULL() :type(snk_wavefunc_obj), pointer

quark propagator for physics measurement

Subroutine :
this :type(qprop_params), intent(inout)
iout :integer, intent(in)

read parameters for a quark propagator parameter

[Source]

subroutine read_qprop_params(this,iout)
!
! read parameters for a quark propagator parameter
!
  implicit none
  type(qprop_params), intent(inout) :: this
  integer,            intent(in)    :: iout
  real(DP) :: kappa,csw
  integer, parameter :: Nfone=1

  !
  ! read wilson quark mass parameters
  !

  if (nodeid==0) then
    !
    ! read Solver parameters
    !
#ifdef _CLOVER
    read(iout,*)kappa,csw
#else
    read(iout,*)kappa
#endif
    read(iout,*)this%tol
  endif

  if (NPU > 1) then
    call comlib_bcast(kappa,0)
#ifdef _CLOVER
    call comlib_bcast(csw,0)
#endif
    call comlib_bcast(this%tol,0)
  endif

  call set_nflavor(this%mass,Nfone)
  call set_kappa(this%mass,kappa) 
#ifdef _CLOVER
  call   set_csw(this%mass,csw)
#endif

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

read parameters for a quark propagator

[Source]

subroutine read_qprop(this,iout)
!
! read parameters for a quark propagator
!
  use error_class
  implicit none
  type(quark_prop), intent(inout) :: this
  integer,     intent(in)    :: iout
  call read(this%parameter,iout)
  return
end subroutine
Subroutine :
this :type(quark_prop), intent(inout)

set source on a quark propagator via wavefunction

[Source]

subroutine set_source(this)
!
! set source on a quark propagator via wavefunction
!
  use error_class
  use print_status_class
  implicit none
  type(quark_prop), intent(inout) :: this
  integer :: it,ic,is,isp,jc,js
  integer :: ix,iy,iz,itx,ity,itz,itt
  character(len=CHARLEN) :: str
  complex(DP) :: ctmp

  if (allocated(this%prop)) deallocate(this%prop)
  allocate(this%prop(COL,COL,SPIN,SPIN,NX*NY*NZ,NTT))

  if (.not.associated(this%wavefunc_sous)) then
    call error_stop("set_source: Error, assign wavefunc_sous for qprop.")
  endif

!$OMP PARALLEL DO PRIVATE(it,isp,ic,is,jc,js)
  do it=1,NT
  do isp=1,NSPACE
    do js=1,SPIN
    do is=1,SPIN
    do jc=1,COL
    do ic=1,COL
      this%prop(ic,jc,is,js,isp,it) = Z0
    enddo
    enddo
    enddo
    enddo
    do is=1,SPIN
    do ic=1,COL
      this%prop(ic,ic,is,is,isp,it) = this%wavefunc_sous%wavefunc(ic,is,isp,it)
    enddo
    enddo
  enddo
  enddo

!#define _DDD_
#ifdef _DDD_
!==============
  do ix=1,NX
  do iy=1,NY
  do iz=1,NZ
  do it=1,NT
    itx = ix + ipsite(1)*NX
    ity = iy + ipsite(2)*NY
    itz = iz + ipsite(3)*NZ
    itt = it
    isp = ispace(ix,iy,iz)
    do js=1,SPIN
    do jc=1,COL
    do is=1,SPIN
    do ic=1,COL
      ctmp = this%prop(ic,jc,is,js,isp,it)
      call comlib_sumcast(ctmp)
      if (ABS(ctmp) > 0.0_DP) then
        write(str,'(4I4,2I2,2I2,2E24.15," SETSRC")') itx,ity,itz,itt,ic,is,jc,js, this%prop(ic,jc,is,js,isp,it)
         call print_status(str)
      endif
    enddo
    enddo
    enddo
    enddo
  enddo
  enddo
  enddo
  enddo
!==============
#endif

  return
end subroutine