Class hadron_class
In: SingletMesons_Simple_v1.5/hadron_class.F90
lattice_class error_class timer_class comlib field_gauge_class config_profile_class wavefunc_class fftw3_class quark_clover_class quark_wilson_class quark_prop_class hadron_class dot/f_27.png

This contains hadron propagators with up/down(degenerate) and strange quarks.

Methods

Included Modules

lattice_class error_class timer_class comlib field_gauge_class config_profile_class wavefunc_class fftw3_class quark_clover_class quark_wilson_class quark_prop_class

Public Instance methods

HADRON_OUTPUT_FORMAT
Constant :
HADRON_OUTPUT_FORMAT ="(I3,2X,E24.15,2X,E24.15)" :character(len=24), parameter
NUM_BARYON_OP
Constant :
NUM_BARYON_OP =NUM_DELTA_BARYON*2+4*NUM_OCTET_BARYON :integer, parameter
NUM_DELTA_BARYON
Constant :
NUM_DELTA_BARYON =4 :integer, parameter
NUM_MESON_OP
Constant :
NUM_MESON_OP =16+2+9+3 :integer, parameter
NUM_MESON_OP_LOCAL
Constant :
NUM_MESON_OP_LOCAL =18 :integer, parameter
NUM_MESON_OP_MOM
Constant :
NUM_MESON_OP_MOM =4 :integer, parameter
NUM_MOM
Constant :
NUM_MOM =19 :integer, parameter
NUM_OCTET_BARYON
Constant :
NUM_OCTET_BARYON =2 :integer, parameter
baryon_prop_obj
Derived Type :
prop(NTT,0:NUM_BARYON_OP-1) :complex(DP)
quark_mass(3) :type(qprop_params)
fname ="hpropB" :character(len=CHARLEN)
path :character(len=CHARLEN)
run_time :type(timer)
wf_sous(3) :type(src_wavefunc_obj)
wf_sink(3) :type(snk_wavefunc_obj)
config_prof => NULL() :type(config_profile_obj), pointer
smear_sous(3) :integer
smear_sink(3) :integer
iout :integer
itraj :integer

Baryon propagators

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

[Source]

subroutine delete_baryon(this)
  implicit none
  type(baryon_prop_obj), intent(inout) :: this

  if (nodeid==0) close(unit=this%iout)

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

[Source]

subroutine delete_meson_mom(this)
  implicit none
  type(meson_mom_prop_obj), intent(inout) :: this

  if (nodeid==0) close(unit=this%iout)

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

[Source]

subroutine delete_meson(this)
  implicit none
  type(meson_prop_obj), intent(inout) :: this

  if (nodeid==0) close(unit=this%iout)

  return
end subroutine
get_prop( this, u )
Subroutine :
this :type(quark_prop), intent(inout)
u :type(vfield_gluon_wg), intent(in)

compute quark propagator

Original external subprogram is quark_prop_class#get_prop

hadrons_obj
Derived Type :
idummy(4) :integer
meson_mom_prop_obj
Derived Type :
prop(NTT,1:NUM_MESON_OP_MOM,1:NUM_MOM) :complex(DP)
quark_mass(2) :type(qprop_params)
fname ="hpropMP" :character(len=CHARLEN)
path :character(len=CHARLEN)
run_time :type(timer)
wf_sous(2) :type(src_wavefunc_obj)
wf_sink(2) :type(snk_wavefunc_obj)
config_prof => NULL() :type(config_profile_obj), pointer
smear_sous(2) :integer
smear_sink(2) :integer
iout :integer
itraj :integer

Meson propagators with finite momentum

meson_prop_obj
Derived Type :
prop(NTT,0:NUM_MESON_OP-1) :complex(DP)
quark_mass(2) :type(qprop_params)
fname ="hpropM" :character(len=CHARLEN)
path :character(len=CHARLEN)
run_time :type(timer)
wf_sous(2) :type(src_wavefunc_obj)
wf_sink(2) :type(snk_wavefunc_obj)
config_prof => NULL() :type(config_profile_obj), pointer
smear_sous(2) :integer
smear_sink(2) :integer
iout :integer
itraj :integer

Meson propagators

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

Initialize for baryon propagators

[Source]

subroutine new_baryon(this)
!
! Initialize for baryon propagators
!
  implicit none
  type(baryon_prop_obj), intent(inout) :: this
  character(len=CHARLEN) :: fname,str
  integer :: err
  logical :: flag

!  if (nodeid==0) write(*,*)"baryon_mom_prop_obj:",sizeof(this),mod(sizeof(this),16)
  call new(this%run_time)
  this%prop = z0

  if (nodeid==0) then
    do
      inquire(unit=io_unit,opened=flag)
      if (flag) then
        io_unit=io_unit+1
      else
        this%iout=io_unit
        exit
      endif
    enddo
  endif
#ifndef _singlePU
  call comlib_bcast(this%iout,0)
#endif

  write(fname,'(I6.6)')this%itraj
  fname=TRIM(this%path)//'/'//TRIM(this%fname)//'.'//TRIM(fname)
  if (nodeid==0) open(unit=this%iout,file=TRIM(fname),form='formatted',status='new',iostat=err)
#ifndef _singlePU
  call comlib_bcast(err,0)
#endif
  if (err/=0) then
    write(str,'(" Error: file open ",A)')TRIM(fname)
    call error_stop(str)
  endif
  call print_config_prof(this%iout,this%config_prof)

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

Initialize for Hadron measurement

[Source]

subroutine new_hadron(this)
!
! Initialize for Hadron measurement
!
  implicit none
  type(hadrons_obj), intent(inout) :: this

  call set_momentum_table
  this%idummy(:)=1

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

Initialize for meson propagators with finite momentum

[Source]

subroutine new_meson_mom(this)
!
! Initialize for meson propagators with finite momentum
!
  implicit none
  type(meson_mom_prop_obj), intent(inout) :: this
  character(len=CHARLEN) :: fname,str
  integer :: err
  logical :: flag

!  if (nodeid==0) write(*,*)"meson_mom_prop_obj:",sizeof(this),mod(sizeof(this),16)
  call new(this%run_time)
  this%prop = Z0

  if (nodeid==0) then
    do
      inquire(unit=io_unit,opened=flag)
      if (flag) then
        io_unit=io_unit+1
      else
        this%iout=io_unit
        exit
      endif
    enddo
  endif

#ifndef _singlePU
  call comlib_bcast(this%iout,0)
#endif

  write(fname,'(I6.6)')this%itraj
  fname=TRIM(this%path)//'/'//TRIM(this%fname)//'.'//TRIM(fname)
  if (nodeid==0) open(unit=this%iout,file=TRIM(fname),form='formatted',status='new',iostat=err)
#ifndef _singlePU
  call comlib_bcast(err,0)
#endif
  if (err/=0) then
    write(str,'(" Error: file open ",A)')TRIM(fname)
    call error_stop(str)
  endif
  call print_config_prof(this%iout,this%config_prof)

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

Initialize for meson propagators

[Source]

subroutine new_meson(this)
!
! Initialize for meson propagators
!
  implicit none
  type(meson_prop_obj), intent(inout) :: this
  character(len=CHARLEN) :: fname, str
  integer :: err
  logical :: flag

!  if (nodeid==0) write(*,*)"meson_prop_obj:",sizeof(this),mod(sizeof(this),16)
  call new(this%run_time)
  this%prop = Z0

  if (nodeid==0) then
    do
      inquire(unit=io_unit,opened=flag)
      if (flag) then
        io_unit=io_unit+1
      else
        this%iout=io_unit
        exit
      endif
    enddo
  endif

#ifndef _singlePU
  call comlib_bcast(this%iout,0)
#endif

  write(fname,'(I6.6)')this%itraj
  fname=TRIM(this%path)//'/'//TRIM(this%fname)//'.'//TRIM(fname)
  if (nodeid==0) open(unit=this%iout,file=TRIM(fname),form='formatted',status='new',iostat=err)
#ifndef _singlePU
  call comlib_bcast(err,0)
#endif
  if (err/=0) then
    write(str,'(" Error: file open ",A)')TRIM(fname)
    call error_stop(str)
  endif
  call print_config_prof(this%iout,this%config_prof)

  return
end subroutine
Subroutine :
iout :integer, intent(in)
config_prof :type(config_profile_obj), intent(in)

[Source]

subroutine print_config_prof(iout,config_prof)
  implicit none
  integer, intent(in):: iout
  type(config_profile_obj), intent(in) :: config_prof
  integer :: iqk
  if (nodeid==0) then
    write(iout,'("begin config")')
    write(iout,'("  lattice_size(x,y,z,t) = ",4I3)') config_prof%NTX, config_prof%NTY, config_prof%NTZ, config_prof%NTT
    write(iout,'("    PLQ = ",F24.15)')        config_prof%plq
    write(iout,'(" config_number = ",I10)')    config_prof%total_trajectory
    write(iout,'("         CDATE = ",A)') TRIM(config_prof%cdate)
    write(iout,'("         CTIME = ",A)') TRIM(config_prof%ctime)
    write(iout,'("         CZONE = ",A)') TRIM(config_prof%czone)
    write(iout,'("       COMMENT = ",A)') TRIM(config_prof%comment)
    write(iout,'("end config")')
  endif
  return
end subroutine
Subroutine :
this :type(baryon_prop_obj), intent(inout)
iflv :integer, intent(in)

[Source]

subroutine save_baryon_prop(this,iflv)
  implicit none
  type(baryon_prop_obj), intent(inout) :: this
  integer, intent(in) :: iflv
  integer :: itt,iop,iqk

  if (nodeid==0) then
    write(this%iout,'("begin baryon",I3)') iflv

    do iqk=1,3
      write(this%iout,'("begin quark_val ",I1)')iqk
      write(this%iout,'(" kappa_val = ",F11.8)') get_kappa(this%quark_mass(iqk)%mass)
      write(this%iout,'("   csw_val = ",F11.8)')   get_csw(this%quark_mass(iqk)%mass)
      write(this%iout,'(" source_center(x,y,z,t) = ",4I3)') this%wf_sous(iqk)%itx0, this%wf_sous(iqk)%ity0, this%wf_sous(iqk)%itz0, this%wf_sous(iqk)%itt0
      write(this%iout,'(" smear_source = ",I3)') this%smear_sous(iqk)
      write(this%iout,'(" smear_sink   = ",I3)') this%smear_sink(iqk)
      write(this%iout,'("end quark_val ",I1)')iqk
    enddo

    write(this%iout,'(I3,2X,3F15.8,3X,3I3)')iflv, get_kappa(this%quark_mass(1)%mass), get_kappa(this%quark_mass(2)%mass), get_kappa(this%quark_mass(3)%mass), this%smear_sous(1), this%smear_sous(2), this%smear_sous(3)

    do iop=0,NUM_BARYON_OP-1
      write(this%iout,'("op. ",I3," begin")')iop
      do itt=1,NTT
        write(this%iout,HADRON_OUTPUT_FORMAT)itt-1,this%prop(itt,iop)
      enddo
      write(this%iout,'("op. ",I3," end")')iop
    enddo

    write(this%iout,'("end baryon",I3)') iflv
  endif

  return
end subroutine
Subroutine :
this :type(meson_mom_prop_obj), intent(inout)
iflv :integer, intent(in)

[Source]

subroutine save_meson_mom_prop(this,iflv)
  implicit none
  type(meson_mom_prop_obj), intent(inout) :: this
  integer, intent(in) :: iflv
  integer :: itt,iop,imom,iqk

  if (nodeid==0) then
    write(this%iout,'("begin meson",I3)') iflv

    do iqk=1,2
      write(this%iout,'("begin quark_val ",I1)')iqk
      write(this%iout,'(" kappa_val = ",F11.8)') get_kappa(this%quark_mass(iqk)%mass)
      write(this%iout,'("   csw_val = ",F11.8)')   get_csw(this%quark_mass(iqk)%mass)
      write(this%iout,'(" source_center(x,y,z,t) = ",4I3)') this%wf_sous(iqk)%itx0, this%wf_sous(iqk)%ity0, this%wf_sous(iqk)%itz0, this%wf_sous(iqk)%itt0
      write(this%iout,'(" smear_source = ",I3)') this%smear_sous(iqk)
      write(this%iout,'(" smear_sink   = ",I3)') this%smear_sink(iqk)
      write(this%iout,'("end quark_val ",I1)')iqk
    enddo

    write(this%iout,'(I3,2X,2F15.8,2X,2I3)')iflv, get_kappa(this%quark_mass(1)%mass), get_kappa(this%quark_mass(2)%mass), this%smear_sous(1), this%smear_sous(2)

    do iop=1,NUM_MESON_OP_MOM
      write(this%iout,'("op. ",I3," begin")')iop
      do imom=1,NUM_MOM
        write(this%iout,'("%momentum: ",I3," start")')imom
        do itt=1,NTT
          write(this%iout,HADRON_OUTPUT_FORMAT)itt-1,this%prop(itt,iop,imom)
        enddo
        write(this%iout,'("%momentum: ",I3," finish")')imom
      enddo
      write(this%iout,'("op. ",I3," end")')iop
    enddo

    write(this%iout,'("end meson",I3)') iflv
  endif

  return
end subroutine
Subroutine :
this :type(meson_prop_obj), intent(inout)
iflv :integer, intent(in)

[Source]

subroutine save_meson_prop(this,iflv)
  implicit none
  type(meson_prop_obj), intent(inout) :: this
  integer, intent(in) :: iflv
  integer :: itt,iop,iqk

  if (nodeid==0) then
    write(this%iout,'("begin meson",I3)') iflv

    do iqk=1,2
      write(this%iout,'("begin quark_val ",I1)')iqk
      write(this%iout,'(" kappa_val = ",F11.8)') get_kappa(this%quark_mass(iqk)%mass)
      write(this%iout,'("   csw_val = ",F11.8)')   get_csw(this%quark_mass(iqk)%mass)
      write(this%iout,'(" source_center(x,y,z,t) = ",4I3)') this%wf_sous(iqk)%itx0, this%wf_sous(iqk)%ity0, this%wf_sous(iqk)%itz0, this%wf_sous(iqk)%itt0
      write(this%iout,'(" smear_source = ",I3)') this%smear_sous(iqk)
      write(this%iout,'(" smear_sink   = ",I3)') this%smear_sink(iqk)
      write(this%iout,'("end quark_val ",I1)')iqk
    enddo

    write(this%iout,'(I3,2X,2F15.8,2X,2I3)')iflv, get_kappa(this%quark_mass(1)%mass), get_kappa(this%quark_mass(2)%mass), this%smear_sous(1), this%smear_sous(2)

    do iop=0,NUM_MESON_OP-1
      write(this%iout,'("op. ",I3," begin")')iop
      do itt=1,NTT
        write(this%iout,HADRON_OUTPUT_FORMAT)itt-1,this%prop(itt,iop)
      enddo
      write(this%iout,'("op. ",I3," end")')iop
    enddo

    write(this%iout,'("end meson",I3)') iflv
  endif

  return
end subroutine