Class quark_wilson_class
In: QuarkWilsonClass/quark_wilson_class.F90
lattice_class hmc_status_class field_gauge_class field_fermion_class action_base_class comlib counter_class quark_wilson_class dot/f_95.png

Definition of Wilson quark related parameters and routines

version

$Id: quark_wilson_class.F90,v 1.20 2011/09/25 06:13:16 ishikawa Exp $

Methods

Included Modules

lattice_class hmc_status_class field_gauge_class field_fermion_class action_base_class comlib counter_class

Public Instance methods

Subroutine :
qp1 :type(quark_wilson), intent(inout)
qp2 :type(quark_wilson), intent(in)

[Source]

subroutine assign_quark_wilson(qp1,qp2)
  implicit none
  type(quark_wilson), intent(inout) :: qp1
  type(quark_wilson), intent(in)    :: qp2
  qp1 = qp2
  return
end subroutine
Subroutine :
this :type(quark_wilson), intent(in)
: wilson quark parameter
yde :type(field_quark_eo_wg), intent(inout)
: even/odd fermion field
ye :type(field_quark_eo_wg), intent(inout)
: even/odd fermion field / boundary copied
u :type(vfield_gluon_wg), intent(in)
: gauge field

Multiply even/odd-site preconditioned Dirac operator :

       v
 yde = Dee ye,

 where  v
        Dee = 1ee - kappa^2 Meo Moe

[Source]

subroutine assign_mult_eoprec_wd_eo(this,yde,ye,u)
!
!= Multiply even/odd-site preconditioned Dirac operator :
!        v
!  yde = Dee ye,
!
!  where  v
!         Dee = 1ee - kappa^2 Meo Moe
! 
  implicit none
  type(quark_wilson),      intent(in)    :: this  ! wilson quark parameter
  type(field_quark_eo_wg), intent(inout) :: yde   ! even/odd fermion field
  type(field_quark_eo_wg), intent(inout) ::  ye   ! even/odd fermion field / boundary copied
  type(vfield_gluon_wg),   intent(in)    ::   u   ! gauge field

  type(field_quark_eo_wg), allocatable :: Myo
  real(DP) :: kappa2

  allocate(Myo)
  kappa2=-(this%kappa)**2

!***********************************************************
!* Calc:
!*        v
!*  yde = D ye
!*      = ( 1 - kappa**2 Meo Moe ) ye
!*      
!***********************************************************

!***********************************
!*  Myo = F^-1oo Moe ye
!*  yde = F^-1ee Meo Myo = F^-1ee Meo F^-1oo Moe ye
!***********************************
  call assign_mult_hop(this,Myo, ye,u)
  call assign_mult_hop(this,yde,Myo,u)

!***********************************
!*  yde = ye - kappa^2 yde
!*      = ye - kappa^2 Meo Moe ye
!***********************************
  call accum_mult_add(yde,kappa2,ye)

  deallocate(Myo)

  return
end subroutine
Subroutine :
this :type(quark_wilson), intent(in)
: Wilson quark parameters
My :type(field_quark_wg), intent(out)
: even/odd site fermion field
y :type(field_quark_wg), intent(inout)
: odd/even site fermion field
u :type(vfield_gluon_wg), intent(in)
: gauge field

Multiply Hopping operator

  My = M[u] y

[Source]

subroutine assign_mult_hopping(this,My,y,u)
!
!= Multiply Hopping operator
!
!   My = M[u] y
!
  use field_gauge_class
  implicit none
  type(quark_wilson),    intent(in)    :: this ! Wilson quark parameters
  type(field_quark_wg),  intent(out)   ::  My  ! even/odd site fermion field
  type(field_quark_wg),  intent(inout) ::   y  ! odd/even site fermion field
  type(vfield_gluon_wg), intent(in)    ::   u  ! gauge field
  call assign_mult_hop(this,My%eo(0),y%eo(1),u)
  call assign_mult_hop(this,My%eo(1),y%eo(0),u)
  return
end subroutine
Subroutine :
this :type(quark_wilson), intent(in)
: Wilson quark parameters
Mye :type(field_quark_eo_wg), intent(out)
: even/odd site fermion field
yo :type(field_quark_eo_wg), intent(inout)
: odd/even site fermion field
u :type(vfield_gluon_wg), intent(in)
: gauge field

Multiply Hopping operator

  Mye = Meo[u] yo (hoppng odd->even/even->odd)

[Source]

subroutine assign_mult_hopping_eo(this,Mye,yo,u)
!
!= Multiply Hopping operator
!
!   Mye = Meo[u] yo (hoppng odd->even/even->odd)
!
  use field_gauge_class
  use counter_class
  implicit none
  type(quark_wilson),      intent(in)    :: this ! Wilson quark parameters
  type(field_quark_eo_wg), intent(out)   :: Mye  ! even/odd site fermion field
  type(field_quark_eo_wg), intent(inout) ::  yo  ! odd/even site fermion field
  type(vfield_gluon_wg),   intent(in)    ::   u  ! gauge field

  call copy_boundary(yo)
  call mult_hopping_tzyx_eo(Mye,yo,u)
  call inc(mult_iter)

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

[Source]

subroutine delete_quark_wilson(this)
  implicit none
  type(quark_wilson), intent(inout) :: this
  this%action_name = REPEAT(' ',LEN(this%action_name))
  return
end subroutine
get_id( this ) result(id)
Function :
id :integer
this :class(action_parameters), intent(in)

Original external subprogram is action_base_class#get_id

Function :
kappa :real(DP)
this :type(quark_wilson), intent(in)

[Source]

function get_kappa(this) result(kappa)
  type(quark_wilson), intent(in) :: this
  real(DP) :: kappa
  kappa = this%kappa
  return
end function
Function :
nf :integer
this :type(quark_wilson), intent(in)

[Source]

function get_nflavor(this) result(nf)
  type(quark_wilson), intent(in) :: this
  integer :: nf
  nf = this%nflavor
  return
end function
Subroutine :
this :type(quark_wilson), intent(inout)
id :integer, intent(in)

[Source]

subroutine new_quark_wilson(this,id)
  implicit none
  type(quark_wilson), intent(inout) :: this
  integer,            intent(in)    :: id
  this%kappa   = 0.0_DP
  this%nflavor = 0
  this%action_name = ACTION_NAME
  return
end subroutine
Subroutine :
this :type(quark_wilson), intent(inout)

[Source]

subroutine print_quark_wilson(this)
  implicit none
  type(quark_wilson), intent(inout) :: this
  write(*,'(14X,"    Num of Flavor :",I3)')   this%nflavor
  write(*,'(14X,"            Kappa :",F12.8)')this%kappa
  return
end subroutine
quark_wilson
Derived Type :

Wilson quark type extended from the above type (extends from action_parameters)

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

[Source]

subroutine read_quark_wilson(this,iout)
  implicit none
  type(quark_wilson), intent(inout) :: this
  integer, intent(in) :: iout
  if (nodeid==0) then
    read(iout,*)this%nflavor
    read(iout,*)this%kappa
  endif
#ifndef _singlePU
  call comlib_bcast(this%kappa,0)
  call comlib_bcast(this%nflavor,0)
#endif
  return
end subroutine
Subroutine :
this :type(quark_wilson), intent(in)
iout :integer, intent(in)

[Source]

subroutine save_config_quark_wilson(this,iout)
  implicit none
  type(quark_wilson), intent(in) :: this
  integer, intent(in) :: iout
  write(iout)this%nflavor
  write(iout)this%kappa
  return
end subroutine
set_gaussian_noise( p )
Subroutine :
p :type(vfield_gluon_wog), intent(inout)
: gauge momentum field

set Gaussian noise on canonical momentum (su(3) Lie algebra) of gauge field (SU(3) Lie group)

Original external subprogram is field_gauge_class#set_gaussian_noise

set_gaussian_noise( y )
Subroutine :
y :type(field_quark_wg), intent(inout)

set Gaussian noise on y

Original external subprogram is field_fermion_class#set_gaussian_noise

set_gaussian_noise( ye )
Subroutine :
ye :type(field_quark_eo_wg), intent(inout)

set Gaussian noise on even/odd sites only

Original external subprogram is field_fermion_class#set_gaussian_noise

set_id( this, id )
Subroutine :
this :class(action_parameters), intent(inout)
id :integer, intent(in)

Original external subprogram is action_base_class#set_id

Subroutine :
this :type(quark_wilson), intent(inout)
kappa :real(DP), intent(in)

[Source]

subroutine set_kappa(this,kappa)
  type(quark_wilson), intent(inout) :: this
  real(DP),           intent(in)    :: kappa
  this%kappa = kappa
  return
end subroutine
Subroutine :
this :type(quark_wilson), intent(inout)
nflavor :integer, intent(in)

[Source]

subroutine set_nflavor(this,nflavor)
  type(quark_wilson), intent(inout) :: this
  integer,            intent(in)    :: nflavor
  this%nflavor = nflavor
  return
end subroutine