Class | hmc_class |
In: |
HmcClass/hmc_class.F90
|
type(hmc_algorithm) :: hmc call new(hmc) call read(hmc) call print(hmc) call run(hmc) call save(hmc) call print_statistics(hmc) call delete(hmc)
$Id: hmc_class.F90,v 1.18 2011/10/27 01:12:10 ishikawa Exp $
Subroutine : | |
this : | type(hmc_algorithm), intent(inout) |
Delelete HMC algorithm
subroutine delete_hmc(this) ! ! Delelete HMC algorithm ! use comlib implicit none type(hmc_algorithm), intent(inout) :: this call delete(this%qcd) call delete(this%hmc_status) call delete(this%md_status) #ifndef _singlePU call comlib_finalize #endif return end subroutine
Subroutine : | |
this : | type(hmc_algorithm), intent(inout) |
Initialize and Create HMC algorithm
subroutine new_hmc(this) ! ! Initialize and Create HMC algorithm ! use comlib implicit none type(hmc_algorithm), intent(inout) :: this character(len=8) :: cdate character(len=10) :: ctime integer :: iflg character(CHARLEN) :: str call tic(this%timing_program) call new(this%hmc_status) call new(this%md_status) call new(this%qcd,this%hmc_status) this%file_name=repeat(' ',len(this%file_name)) this%current_directory=repeat(' ',len(this%current_directory)) this%job_id=repeat(' ',len(this%job_id)) call DATE_AND_TIME(date=cdate,time=ctime) iflg=0 if (nodeid == 0) then if(COMMAND_ARGUMENT_COUNT()==1) then call GET_COMMAND_ARGUMENT(1,this%file_name) else iflg=1 endif endif #ifndef _singlePU call comlib_bcast(iflg,0) #endif if (iflg==1) then write(str,'("HMC new error: ",A,I4)')__FILE__,__LINE__ call error_stop(TRIM(str)) endif if (nodeid == 0) then call getcwd(this%current_directory) endif write(this%job_id,'(A,1X,A)')TRIM(ADJUSTL(cdate)),TRIM(ADJUSTL(ctime)) if (nodeid == 0) then write(*,'(" DATE AND TIME :",A)')TRIM(this%job_id) endif #ifndef _singlePU call comlib_bcast(this%file_name,0) call comlib_bcast(this%current_directory,0) call comlib_bcast(this%job_id,0) #endif this%config_path=repeat(' ',len(this%config_path)) this%config_file_name_header=repeat(' ',len(this%config_file_name_header)) return end subroutine
Subroutine : | |
this : | type(hmc_algorithm), intent(inout) |
Print HMC parameters
this should be called after reading HMC parameters by read method.
subroutine print_hmc(this) ! ! Print HMC parameters ! ! this should be called after reading HMC parameters by read method. ! implicit none type(hmc_algorithm), intent(inout) :: this call tic(this%timing_io) if (nodeid == 0) then write(*,'(80("-"))') write(*,'("**** Input file ****")') write(*,'(8X," JOB Input file :",A)')TRIM(this%file_name) write(*,*) write(*,'("**** Configuration files ****")') write(*,'(8X,"Configuration directory :",A)')TRIM(this%config_path) write(*,'(8X,"Configuration file name :",A)')TRIM(this%config_file_name_header) endif call print(this%qcd) if (nodeid == 0) then write(*,*) write(*,'("**** HMC Statistics ****")') endif call print(this%hmc_status) if (nodeid == 0) then write(*,*) write(*,'("**** MD Statistics ****")') endif call print(this%md_status) if (nodeid == 0) then write(*,*) write(*,'(80("-"))') endif call toc(this%timing_io) return end subroutine
Subroutine : | |
this : | type(hmc_algorithm), intent(inout) |
Print HMC statistics at the end of job
subroutine print_stat_hmc(this) ! ! Print HMC statistics at the end of job ! use comlib implicit none type(hmc_algorithm), intent(inout) :: this real(DP) :: hmonte_time real(DP) :: metropolis_time real(DP) :: io_time real(DP) :: program_time hmonte_time = get_elapse(this%timing_hmonte) metropolis_time = get_elapse(this%timing_metropolis) io_time = get_elapse(this%timing_io) call toc(this%timing_program) program_time = get_elapse(this%timing_program) if (nodeid == 0) then write(*,'(80("-"))') write(*,'(" Statistics par 1PU")') write(*,'(" Total time (sec) :",F12.4)') program_time write(*,'(" monte time (sec) :",F12.4," (",F8.4," %)")') hmonte_time, 100*hmonte_time/program_time write(*,'(" metropolis time (sec) :",F12.4," (",F8.4," %)")') metropolis_time, 100*metropolis_time/program_time write(*,'(" IO conf time (sec) :",F12.4," (",F8.4," %)")') io_time, 100*io_time/program_time endif if (nodeid == 0) then write(*,'(80("-"))') write(*,'(" HMC Statistics")') endif call print_statistics(this%hmc_status) call print_statistics(this%qcd,program_time) if (nodeid == 0) then write(*,'(80("-"))') endif return end subroutine
Subroutine : | |
this : | type(hmc_algorithm), intent(inout) |
Read HMC job parameter file
This should be called after calling new method.
subroutine read_hmc(this) ! ! Read HMC job parameter file ! ! This should be called after calling new method. ! use comlib use file_tools_class use print_status_class implicit none type(hmc_algorithm), intent(inout) :: this character(len=CHARLEN) :: char integer :: iout call tic(this%timing_io) iout = search_free_file_unit() if (nodeid==0) then open(iout,file=this%file_name,status='unknown',form='formatted') read(iout,'(a)')this%input_file_comment(1) endif !============================ ! read QCD part !============================ call read(this%qcd,iout) !============================ ! read JOB part !============================ if (nodeid==0) then read(iout,'(a)')this%input_file_comment(2) read(iout,*) this%random_number_seed read(iout,'(a)')this%config_path read(iout,'(a)')this%config_file_name_header read(iout,'(a)')this%input_file_comment(3) endif #ifndef _singlePU call comlib_bcast(this%random_number_seed,0) call comlib_bcast(this%config_path,0) call comlib_bcast(this%config_file_name_header,0) #endif !============================ ! read HMC part !============================ call read(this%hmc_status,iout) if (nodeid==0) then read(iout,'(a)')this%input_file_comment(4) endif !============================ ! read MD part !============================ call read(this%md_status,iout) if (nodeid==0) then close(iout) endif if (nodeid == 0) then write(*,'(80("-"))') endif write(char,'(1X,A," : read OK ",I3)')TRIM(this%file_name),nodeid call print_status(char) if (nodeid == 0) write(*,'(80("-"))') call toc(this%timing_io) !================================== ! initialize/load configuration ! for continuation ! random number generator ! is initialized here !================================== call initialize_config_hmc(this) return end subroutine
Subroutine : | |
this : | type(hmc_algorithm), intent(inout) |
Run Hybrid Monte Carlo (HMC) Algorithm
subroutine run_hmc(this) ! ! Run Hybrid Monte Carlo (HMC) Algorithm ! implicit none type(hmc_algorithm), intent(inout) :: this integer :: ihit,isweep logical :: flag isweep = 0 do !=================================== ! check sweep/iteration count !=================================== flag = is_do_loop_ended(this%hmc_status) if (flag) exit !============================================ ! Hybrid Monte Carlo !============================================ call tic(this%timing_hmonte) ihit=0 !================================ ! P is generated and ths initial ! Hamiltonian is calcd. !================================ call hamil(this%qcd,0) #ifdef _REVCHECK call rev_check(this%qcd,0) #endif !================================ ! Molecular Dynamics !================================ call molecular_dynamics(this) !================================ ! Ths Final Hamiltonian is calcd. !================================ call hamil(this%qcd,1) #ifdef _REVCHECK call rev_check(this%qcd,1) call set_time_step_backward(this%md_status) call molecular_dynamics(this) call hamil(this%qcd,2) call set_time_step_forward(this%md_status) call rev_check(this%qcd,2) #endif call toc(this%timing_hmonte) !================================ ! Do Metropolis test; ! Q are updated by accept/reject ! and P is dropped !================================ call tic(this%timing_metropolis) call metropolis(this%qcd,ihit) call toc(this%timing_metropolis) !================================ ! Do measuremet on the fly !================================ call tic(this%timing_io) call measurement(this%qcd) call toc(this%timing_io) !================================ ! Save configuration for ! measurement !================================ call save_config_hmc(this) !======================= ! increment sweep count !======================= isweep = isweep + 1 enddo return end subroutine
Subroutine : | |
this : | type(hmc_algorithm), intent(inout) |
Save HMC job parameter
This should be called after the method run ended.
subroutine save_hmc(this) ! ! Save HMC job parameter ! ! This should be called after the method run ended. ! use file_tools_class implicit none type(hmc_algorithm), intent(inout) :: this integer :: iout !======================================= ! save configuration for continuation !======================================= call save_contg_hmc(this) call tic(this%timing_io) iout = search_free_file_unit() if (nodeid == 0) then open(iout,file=this%file_name,status='unknown',form='formatted') write(iout,'(a)')TRIM(this%input_file_comment(1)) endif call save(this%qcd,iout) if (nodeid == 0) then write(iout,'(a)')TRIM(this%input_file_comment(2)) write(iout,*) this%random_number_seed write(iout,'(a)')TRIM(this%config_path) write(iout,'(a)')TRIM(this%config_file_name_header) write(iout,'(a)')TRIM(this%input_file_comment(3)) endif call save(this%hmc_status,iout) if (nodeid == 0) then write(iout,'(a)')TRIM(this%input_file_comment(4)) endif call save(this%md_status,iout) if (nodeid == 0) then close(iout) endif call toc(this%timing_io) return end subroutine