Class | main_singlet_class |
In: |
SingletMesons_Simple_v1.5/main_singlet_class.F90
|
This contains singlet meson measurement controll routines
Subroutine : | |
this : | type(main_singlet_obj), intent(inout) |
subroutine delete_main_singlet(this) implicit none type(main_singlet_obj), intent(inout) :: this call delete(m_rand) call delete(this%hadrons) return end subroutine
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
Original external subprogram is main_hadron_class#main_hadron_obj
Derived Type : | |||
hadrons : | type(main_hadron_obj) | ||
singlet_parameter_file : | character(len=CHARLEN) | ||
seed_fpath : | character(len=CHARLEN) | ||
seed_fname : | character(len=CHARLEN) | ||
run_number =0 : | integer
| ||
random_number_seed =1 : | integer
| ||
nnoise =1 : | integer
| ||
idummy1 : | integer |
singlet meson measurement parameters
Subroutine : | |
this : | type(main_singlet_obj), intent(inout) |
subroutine new_main_singlet(this) implicit none type(main_singlet_obj), intent(inout) :: this type(gamma_matrix_obj), save :: gamma_matrix call new(this%hadrons) this%singlet_parameter_file="" this%seed_fpath="" this%seed_fname="singlet_seed." call new(gamma_matrix) return end subroutine
Subroutine : | |
this : | type(main_singlet_obj), intent(inout) |
subroutine print_main_singlet(this) use file_tools_class use hmc_status_class, only : fname_cont implicit none type(main_singlet_obj), intent(inout) :: this character(len=CHARLEN) :: fname,char call print(this%hadrons) if (nodeid==0) then write(*,'(80("="))') write(*,'(" Singlet Parameter file :",A)')TRIM(this%singlet_parameter_file) write(*,'(" Run # :",I10)')this%run_number write(*,'(" Random Seed :",I10)')this%random_number_seed write(*,'(" Seed Save path :",A)')TRIM(this%seed_fpath) write(*,'(" Seed Save name :",A)')TRIM(this%seed_fname) write(*,'(" # of Noise Estimate :",I10)')this%nnoise write(*,'(80("="))') endif !===================================== ! Initialize Random Number generator !===================================== fname=fname_cont(this%run_number,this%seed_fpath,this%seed_fname) select case (this%run_number) case (0) call new(m_rand,0,this%random_number_seed,NPU,nodeid) if (nodeid == 0) then write(*,'(1X,"Random Number Generator Initialized")') write(*,'(80("-"))') endif case (1:) call new(m_rand,1,this%random_number_seed,NPU,nodeid) m_seed_out = search_free_file_unit() open(m_seed_out,file=fname,status='unknown',form='unformatted',action='read') call read(m_rand,m_seed_out) close(m_seed_out) if (nodeid == 0) then write(*,'(1X,"Random Number Generator Recovered")') write(*,'(80("-"))') endif write(char,'(1X,A," : read OK ",I3)')TRIM(fname),nodeid call print_status(char) end select this%run_number=this%run_number+1 m_is_initialized = .TRUE. return end subroutine
Subroutine : | |
iout : | integer, intent(in) |
this : | type(main_singlet_obj), intent(inout) |
subroutine read_main_singlet(iout,this) use file_tools_class implicit none integer, intent(in) :: iout type(main_singlet_obj), intent(inout) :: this character(len=CHARLEN) :: str integer :: istat call read(this%hadrons,iout) if (nodeid == 0) then read(iout,*) read(iout,'(A)')this%singlet_parameter_file endif #ifndef _singlePU call comlib_bcast(this%singlet_parameter_file,0) #endif istat=0 if (nodeid == 0) then m_js_out = search_free_file_unit() open(m_js_out,file=TRIM(this%singlet_parameter_file),form='formatted',status='old',iostat=istat) endif #ifndef _singlePU call comlib_bcast(istat,0) #endif if (istat /= 0) goto 100 if (nodeid == 0) then read(m_js_out,*)this%run_number read(m_js_out,*)this%random_number_seed read(m_js_out,'(A)')this%seed_fpath read(m_js_out,'(A)')this%seed_fname read(m_js_out,*)this%nnoise close(m_js_out) endif #ifndef _singlePU call comlib_bcast(this%run_number,0) call comlib_bcast(this%random_number_seed,0) call comlib_bcast(this%seed_fpath,0) call comlib_bcast(this%seed_fname,0) call comlib_bcast(this%nnoise,0) #endif return 100 continue write(str,'("Error Reading Singlet Parameter file.",A)')TRIM(this%singlet_parameter_file) call error_stop(str) end subroutine
Subroutine : | |
this : | type(main_singlet_obj), intent(inout) |
subroutine save_main_singlet(this) use file_tools_class use hmc_status_class, only : fname_cont implicit none type(main_singlet_obj), intent(inout) :: this character(len=CHARLEN) :: fname !===================================== ! Save Random Seed !===================================== fname = fname_cont(this%run_number,this%seed_fpath,this%seed_fname) m_seed_out = search_free_file_unit() open(m_seed_out,file=fname,status='unknown',form='unformatted',action='write') call save(m_rand,m_seed_out) close(m_seed_out) !===================================== ! Save Input Parameters !===================================== if (nodeid == 0) then m_js_out = search_free_file_unit() open(m_js_out,file=TRIM(this%singlet_parameter_file),form='formatted',status='old') write(m_js_out,*)this%run_number write(m_js_out,*)this%random_number_seed write(m_js_out,'(A)')TRIM(this%seed_fpath) write(m_js_out,'(A)')TRIM(this%seed_fname) write(m_js_out,*)this%nnoise close(m_js_out) endif return end subroutine
Subroutine : | |
u : | type(vfield_gluon_wg), intent(inout) |
itx : | integer, intent(in) |
ity : | integer, intent(in) |
itz : | integer, intent(in) |
This is a test routine for lattice translational symmetry. This shift the corrdinate oringin of the gauge link field. Do not use this production run.
subroutine shift(u,itx,ity,itz) ! ! This is a test routine for lattice translational symmetry. ! This shift the corrdinate oringin of the gauge link field. ! Do not use this production run. ! implicit none type(vfield_gluon_wg), intent(inout) :: u integer, intent(in) :: itx,ity,itz type(sfield_gluon_wg), allocatable :: w integer :: idir,ix,iy,iz,mu,Nshift(1:NDIM-1) allocate(w) call new(w) call copy_boundary(u) Nshift(1) = itx Nshift(2) = ity Nshift(3) = itz do idir = 1,NDIM-1 do ix=1,Nshift(idir)-1 do mu=1,NDIM call assign_shift(w%eo(0),u%eo(1)%mu(mu),-idir) call assign_shift(w%eo(1),u%eo(0)%mu(mu),-idir) call assign(u%eo(0)%mu(mu),w%eo(0)) call assign(u%eo(1)%mu(mu),w%eo(1)) enddo call copy_boundary(u) enddo enddo deallocate(w) return end subroutine