Class | hmc_quark_class |
In: |
HmcQuarkClass/hmc_quark_class.F90
|
$Id: hmc_quark_class.F90,v 1.18 2011/10/18 07:18:00 ishikawa Exp $
Variable : | |||
copy_fq_time : | type(timer), save
|
Original external subprogram is field_fermion_class#copy_fq_time
Subroutine : | |
this(:) : | type(hmc_quark), intent(inout) |
delete HMC quark actions (array)
subroutine delete_quark_array(this) ! ! delete HMC quark actions (array) ! implicit none type(hmc_quark), intent(inout) :: this(:) integer :: iqk,nqk nqk = SIZE(this) do iqk=1,nqk call delete_quark(this(iqk)) enddo return end subroutine
Subroutine : | |||
this(:) : | type(hmc_quark),
intent(inout)
| ||
BB : | type(vfield_gluon_wg), intent(inout)
| ||
idepth : | integer, intent(in)
| ||
u : | type(vfield_gluon_wg), intent(in)
|
compute MD Force from HMC QCD quark actions.
subroutine force_quark_array(this,BB,idepth,u) ! ! compute MD Force from HMC QCD quark actions. ! implicit none type(hmc_quark), intent(inout) :: this(:) ! quark actions (array) type(vfield_gluon_wg), intent(inout) :: BB ! MD pre-force (accumulated) integer, intent(in) :: idepth ! MD depth for recursive MD integrator type(vfield_gluon_wg), intent(in) :: u ! gauge field #ifdef _CLOVER type(tfield_gluon_wg), allocatable, target :: XX ! clover term contribution for MD force (need only for clover quarks) #endif integer :: iqk,nqk nqk = SIZE(this) #ifdef _CLOVER !================================= ! Clear clover term contribution !================================= allocate(XX) call new(XX) call clear(XX) !================================= ! Calc clover leaf F_{mu,nu} !================================= call make_clover_leaf(u,m_wl34,m_wl98) #endif do iqk=1,nqk select case (this(iqk)%algorithm) case (PHMC_ALGORITHM) ! PHMC #ifdef _CLOVER this(iqk)%poly_hmc%XX => XX #endif call force(this(iqk)%poly_hmc,BB,idepth,u) case (HMC_ALGORITHM) ! HMC #ifdef _CLOVER this(iqk)%std_hmc%XX => XX #endif call force(this(iqk)%std_hmc,BB,idepth,u) case (MPHMC_ALGORITHM) ! MPHMC #ifdef _CLOVER this(iqk)%mp_hmc%XX => XX #endif call force(this(iqk)%mp_hmc,BB,idepth,u) case (HMC_BLK_ALGORITHM) ! HMCBLK #ifdef _CLOVER this(iqk)%blk_hmc%XX => XX #endif call force(this(iqk)%blk_hmc,BB,idepth,u) case(HMC_ALGORITHM_TRUNCD_DW_OVERLAP) call force(this(iqk)%tdwov_hmc,BB,idepth,u) end select enddo ! end of do iqk #ifdef _CLOVER call delete_clover_leaf() call force_clover_terms(BB,XX,u,m_wl34,m_wl98) deallocate(XX) #endif return end subroutine
Subroutine : | |||
this(:) : | type(hmc_quark),
intent(inout)
| ||
ifi : | integer, intent(in)
| ||
u : | type(vfield_gluon_wg), intent(in)
| ||
SQ(:) : | real(DP), intent(out)
|
compute HMC QCD quark Hamiltonian/Action value
subroutine hamil_quark_array(this,ifi,u,SQ) ! ! compute HMC QCD quark Hamiltonian/Action value ! implicit none type(hmc_quark), intent(inout) :: this(:) ! quark action (array) integer, intent(in) :: ifi ! target switch: 0 for initial Hamiltonian, 1 for final, 2 initial type(vfield_gluon_wg), intent(in) :: u ! gauge field real(DP), intent(out) :: SQ(:) ! Hamiltonian/Action values (array) integer :: iqk,nqk SQ(:) = 0.0_DP nqk = SIZE(this) #ifdef _CLOVER call make_clover_leaf(u,m_wl34,m_wl98) #endif do iqk=1,nqk select case (this(iqk)%algorithm) case (PHMC_ALGORITHM) ! PHMC call hamil(this(iqk)%poly_hmc,ifi,u,SQ(iqk)) case (HMC_ALGORITHM) ! HMC call hamil(this(iqk)%std_hmc,ifi,u,SQ(iqk)) case (MPHMC_ALGORITHM) ! MPHMC call hamil(this(iqk)%mp_hmc,ifi,u,SQ(iqk)) case (HMC_BLK_ALGORITHM) ! HMCBLK call hamil(this(iqk)%blk_hmc,ifi,u,SQ(iqk)) case(HMC_ALGORITHM_TRUNCD_DW_OVERLAP) call hamil(this(iqk)%tdwov_hmc,ifi,u,SQ(iqk)) end select ! end of quark algorithm enddo #ifdef _CLOVER call delete_clover_leaf() #endif return end subroutine
Derived Type : | |
file_name : | character(CHARLEN), public |
HMC QCD quark action
This can have several actions.
Subroutine : | |||
this(:) : | type(hmc_quark),
intent(inout)
| ||
u0 : | type(vfield_gluon_wg), intent(in)
| ||
u1 : | type(vfield_gluon_wg), intent(in)
| ||
ihit : | integer, intent(inout)
|
Do Metropolis test for HMC QCD quark actions (array)
This routine do Metropoils test when the quark action array contains the actions with Global Metropolis test.
subroutine metropolis_quark_array(this,u0,u1,ihit) ! ! Do Metropolis test for HMC QCD quark actions (array) ! ! This routine do Metropoils test when the quark action array contains ! the actions with Global Metropolis test. ! implicit none type(hmc_quark), intent(inout) :: this(:) ! HMC QCD quark actions type(vfield_gluon_wg), intent(in) :: u0 ! MD Initial configuration type(vfield_gluon_wg), intent(in) :: u1 ! MD final configuration integer, intent(inout) :: ihit ! 1 for accept , 0 for reject integer :: iqk,nqk if (ihit == 0) return nqk = SIZE(this) do iqk=1,nqk select case(this(iqk)%algorithm) case(PHMC_ALGORITHM) if (has_global_metropolis_test(this(iqk)%poly_hmc)) then call metropolis(this(iqk)%poly_hmc,u1,u0,ihit) endif case(HMC_ALGORITHM) if (has_global_metropolis_test(this(iqk)%std_hmc)) then call metropolis(this(iqk)%std_hmc,u1,u0,ihit) endif case(MPHMC_ALGORITHM) if (has_global_metropolis_test(this(iqk)%mp_hmc)) then call metropolis(this(iqk)%mp_hmc,u1,u0,ihit) endif case(HMC_BLK_ALGORITHM) if (has_global_metropolis_test(this(iqk)%blk_hmc)) then call metropolis(this(iqk)%blk_hmc,u1,u0,ihit) endif case(HMC_ALGORITHM_TRUNCD_DW_OVERLAP) if (has_global_metropolis_test(this(iqk)%tdwov_hmc)) then call metropolis(this(iqk)%tdwov_hmc,u1,u0,ihit) endif case default call error_stop(" Error in quark param file (algorithm)") end select if (ihit == 0) return enddo return end subroutine
Variable : | |||
mult_iter : | type(counter), save
|
Original external subprogram is field_fermion_class#mult_iter
Subroutine : | |
this(:) : | type(hmc_quark), intent(inout) |
Initialize/create HMC quark action array
subroutine new_quark_array(this) ! ! Initialize/create HMC quark action array ! implicit none type(hmc_quark), intent(inout) :: this(:) integer :: iqk,nqk nqk = SIZE(this) do iqk=1,nqk call new_quark(this(iqk),iqk) enddo return end subroutine
Subroutine : | |
this(:) : | type(hmc_quark), intent(inout) |
print out HMC QCD quark action parameters on display
subroutine print_quark_array(this) ! ! print out HMC QCD quark action parameters on display ! implicit none type(hmc_quark), intent(inout) :: this(:) integer :: iqk,nqk nqk = SIZE(this) do iqk=1,nqk call print_quark(this(iqk)) enddo return end subroutine
Subroutine : |
Original external subprogram is field_fermion_class#print_copy_fq_statistics
Subroutine : | |
this(:) : | type(hmc_quark), intent(inout) |
print out HMC QCD quark action job statistics on display
subroutine print_stat_quark_array(this) ! ! print out HMC QCD quark action job statistics on display ! implicit none type(hmc_quark), intent(inout) :: this(:) integer :: iqk,nqk nqk = SIZE(this) do iqk=1,nqk call print_stat_quark(this(iqk)) enddo return end subroutine
Subroutine : | |
this(:) : | type(hmc_quark), intent(inout) |
Read HMC quark action parameter from formatted file (array)
subroutine read_quark_array(this) ! ! Read HMC quark action parameter from formatted file (array) ! use comlib implicit none type(hmc_quark), intent(inout) :: this(:) integer :: iqk,nqk nqk = SIZE(this) do iqk=1,nqk call read_quark(this(iqk)) enddo return end subroutine
Subroutine : | |||
this(:) : | type(hmc_quark),
intent(inout)
|
subroutine remove_two_links(this) implicit none type(hmc_quark), intent(inout) :: this(:) ! quark action (array) if (associated(m_wl34)) m_wl34 => NULL() if (associated(m_wl98)) m_wl98 => NULL() return end subroutine
Subroutine : | |||
this(:) : | type(hmc_quark),
intent(inout)
| ||
u : | type(vfield_gluon_wg), intent(in) |
Compute Reweighting factor for quark determinant
This routine do Metropoils test when the quark action array contains the actions with Global Metropolis test.
subroutine reweighting_quark_array(this,u) ! ! Compute Reweighting factor for quark determinant ! ! This routine do Metropoils test when the quark action array contains ! the actions with Global Metropolis test. ! implicit none type(hmc_quark), intent(inout) :: this(:) ! HMC QCD quark actions type(vfield_gluon_wg), intent(in) :: u integer :: iqk,nqk nqk = SIZE(this) do iqk=1,nqk select case(this(iqk)%algorithm) case(PHMC_ALGORITHM) if (has_reweighting(this(iqk)%poly_hmc)) then call reweighting(this(iqk)%poly_hmc,u) endif case(HMC_ALGORITHM) if (has_reweighting(this(iqk)%std_hmc)) then call reweighting(this(iqk)%std_hmc,u) endif case(MPHMC_ALGORITHM) if (has_reweighting(this(iqk)%mp_hmc)) then call reweighting(this(iqk)%mp_hmc,u) endif case(HMC_BLK_ALGORITHM) if (has_reweighting(this(iqk)%blk_hmc)) then call reweighting(this(iqk)%blk_hmc,u) endif case(HMC_ALGORITHM_TRUNCD_DW_OVERLAP) if (has_reweighting(this(iqk)%tdwov_hmc)) then call reweighting(this(iqk)%tdwov_hmc,u) endif case default call error_stop(" Error in quark param file (algorithm)") end select enddo return end subroutine
Subroutine : | |
this(:) : | type(hmc_quark), intent(in) |
iout : | integer, intent(in) |
Save HMC quark action parameter on the measurement config file in unformatted form.
subroutine save_config_quark_array(this,iout) ! ! Save HMC quark action parameter on the measurement config file in unformatted form. ! implicit none type(hmc_quark), intent(in) :: this(:) integer, intent(in) :: iout integer :: iqk,nqk nqk = SIZE(this) do iqk=1,nqk call save_config_quark(this(iqk),iout) enddo return end subroutine
Subroutine : | |||
this(:) : | type(hmc_quark),
intent(inout)
| ||
wl34 : | type(tfield_gluon_wg), intent(in), target
| ||
wl98 : | type(tfield_gluon_wg), intent(in), target
|
subroutine set_two_links(this,wl34,wl98) implicit none type(hmc_quark), intent(inout) :: this(:) ! quark action (array) type(tfield_gluon_wg), intent(in), target :: wl34,wl98 ! 2-links fields prepared by make_two_links (need only for clover quarks) if (.not.associated(m_wl34)) m_wl34 => wl34 if (.not.associated(m_wl98)) m_wl98 => wl98 return end subroutine