Class main_hadron_class
In: SingletMesons_Simple_v1.5/main_hadron_class.F90
config_profile_class timer_class comlib lattice_class logfile_class field_gauge_class field_fermion_class quark_clover_class quark_wilson_class quark_prop_class wavefunc_class main_hadron_class dot/f_23.png

This contains Hadron measurement controll routines.

Methods

delete   main_hadron_obj   new   print   read  

Included Modules

config_profile_class timer_class comlib lattice_class logfile_class field_gauge_class field_fermion_class quark_clover_class quark_wilson_class quark_prop_class wavefunc_class

Public Instance methods

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

delete hadron measurement

[Source]

subroutine delete_main(this)
!
! delete hadron measurement
!
  implicit none
  type(main_hadron_obj), intent(inout) :: this
  type(lattice_world)   :: lattice
  integer :: iqk
  real(DP) :: etime

  do iqk=1,this%num_wave
    call delete(this%wavefunc(iqk))
  enddo

  do iqk=1,this%num_quark
    call delete(this%quark_parameter(iqk))
  enddo

  if (associated(this%solver_log)) then
     call delete(this%solver_log)
     deallocate(this%solver_log)
     this%solver_log => NULL()
  endif

  deallocate(this%config_prof)
  deallocate(this%wavefunc)
  deallocate(this%quark_parameter)

  call toc(this%run_time)
  etime = get_elapse(this%run_time)
  if (nodeid==0) then
    write(*,'(80("="))')
    write(*,'(" Hadron ETime:",E24.16)')etime
  endif
  call delete(lattice)

  return
end subroutine
main_hadron_obj
Derived Type :
u => NULL() :type(vfield_gluon_wg), pointer
quark_parameter(:) => NULL() :type(qprop_params), pointer
wavefunc(:) => NULL() :type(src_wavefunc_obj), pointer
config_prof => NULL() :type(config_profile_obj), pointer
config_name :character(len=CHARLEN)
config_path :character(len=CHARLEN)
fname :character(len=CHARLEN)
input_fname :character(len=CHARLEN)
cwd :character(len=CHARLEN)
job_id :character(len=CHARLEN)
cdate :character(len=CHARLEN)
ctime :character(len=CHARLEN)
czone :character(len=CHARLEN)
comment :character(len=CHARLEN)
str :character(len=CHARLEN)
hadron_path :character(len=CHARLEN)
num_quark :integer
num_wave :integer
ini_traj :integer
end_traj :integer
traj_skip :integer
traj :integer
solver_log => NULL() :type(logfile), pointer
run_time :type(timer)

Hadron measurement parameters

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

initialize hadron measurement

[Source]

subroutine new_main(this)
!
! initialize hadron measurement
!
  implicit none
  type(main_hadron_obj), intent(inout) :: this
  type(lattice_world)   :: lattice

  character(len=8)  :: cdate
  character(len=10) :: ctime
  character(len=CHARLEN) :: fname

  call new(lattice)

  call new(this%run_time)
  call tic(this%run_time)

  this%input_fname=repeat(' ',len(this%input_fname))
  this%cwd        =repeat(' ',len(this%cwd))
  this%job_id     =repeat(' ',len(this%job_id))

  call DATE_AND_TIME(date=cdate,time=ctime)
  write(this%job_id,'(A,1X,A)')TRIM(ADJUSTL(cdate)),TRIM(ADJUSTL(ctime))

  if (nodeid == 0) then
    call GET_COMMAND_ARGUMENT(1,this%input_fname)
    call getcwd(this%cwd)
  endif

  if (nodeid == 0) then
    write(*,'(" DATE AND TIME :",A)')TRIM(this%job_id)
  endif

  if (NPU > 1) then
    call comlib_bcast(this%input_fname,0)
    call comlib_bcast(this%cwd,0)
    call comlib_bcast(this%job_id,0)
  endif

  allocate(this%u)
  call new(this%u)

  allocate(this%config_prof)

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

print out hadron measurement parameters

[Source]

subroutine print_main(this)
!
! print out hadron measurement parameters
!
  implicit none
  type(main_hadron_obj), intent(inout) :: this
  type(lattice_world) :: lattice
  integer :: iqk

  call print(lattice)

  if (nodeid==0) then
    write(*,'(80("="))')
    write(*,'("   Hadron Save path :",A)')TRIM(this%hadron_path)
    write(*,'(80("="))')
  endif
  
  return
end subroutine
Subroutine :
this :type(main_hadron_obj), intent(inout)
iout :integer, intent(in)

read Hadron measurement parameters

[Source]

subroutine read_main(this,iout)
!
! read Hadron measurement parameters
!
  implicit none
  type(main_hadron_obj), intent(inout) :: this
  integer,               intent(in) :: iout
  integer :: iqk,iwf
  character(len=CHARLEN) :: char,solver_log_fname

  if (nodeid == 0) then
    open(iout,file=this%input_fname,form='formatted',action='read',status='old')
    read(iout,*)
    read(iout,*)this%ini_traj,this%end_traj,this%traj_skip
    read(iout,'(A)')this%config_path
    read(iout,'(A)')this%config_name
    read(iout,*)
    read(iout,'(A)')this%hadron_path
    read(iout,*)
    read(iout,*)this%num_wave
  endif

  if (NPU > 1) then
    call comlib_bcast(this%ini_traj,0)
    call comlib_bcast(this%end_traj,0)
    call comlib_bcast(this%traj_skip,0)
    call comlib_bcast(this%config_path,0)
    call comlib_bcast(this%config_name,0)
    call comlib_bcast(this%hadron_path,0)
    call comlib_bcast(this%num_wave,0)
  endif

  allocate(this%solver_log)
  write(char,'(I6.6,"-",I6.6)')this%ini_traj,this%end_traj
  solver_log_fname = "quark_solver_log."//TRIM(ADJUSTL(char))
  call new(this%solver_log,solver_log_fname)

  allocate(this%wavefunc(this%num_wave))
  do iqk=1,this%num_wave
    call new(this%wavefunc(iqk))
  enddo

  do iwf=1,this%num_wave
    if (nodeid == 0) then
      read(iout,*)
    endif
    call read(this%wavefunc(iwf),iout)
  enddo

  if (nodeid == 0) then
    read(iout,*)
    read(iout,*)this%num_quark
  endif

  if (NPU > 1) then
    call comlib_bcast(this%num_quark,0)
  endif

  allocate(this%quark_parameter(this%num_quark))
  do iqk=1,this%num_quark
    call new(this%quark_parameter(iqk))
    this%quark_parameter(iqk)%solver_log => this%solver_log
  enddo

  do iqk=1,this%num_quark
    if (nodeid == 0) then
      read(iout,*)
    endif
    call read(this%quark_parameter(iqk),iout)
  enddo

  return
end subroutine