Skip to content

Commit

Permalink
HELIUM: Add keywords for worm move control
Browse files Browse the repository at this point in the history
  • Loading branch information
cschran authored and hforbert committed Sep 12, 2019
1 parent e7023df commit 1cbdb02
Show file tree
Hide file tree
Showing 4 changed files with 135 additions and 43 deletions.
45 changes: 44 additions & 1 deletion src/motion/helium_methods.F
Original file line number Diff line number Diff line change
Expand Up @@ -477,10 +477,53 @@ SUBROUTINE helium_create(helium_env, input, solute)
END IF
helium_env(k)%helium%worm_ln_openclose_scale = LOG(rtmp)

! deal with accptance statistics without changing the ceperley stuff
! deal with acceptance statistics without changing the ceperley stuff
helium_env(k)%helium%maxcycle = 1
helium_env(k)%helium%bisctlog2 = 0

! get the absolute weights of the individual moves
helium_env(k)%helium%worm_all_limit = 0
CALL section_vals_val_get(helium_section, "WORM%CENTROID_WEIGHT", &
i_val=itmp)
helium_env(k)%helium%worm_centroid_min = 1
helium_env(k)%helium%worm_centroid_max = itmp
helium_env(k)%helium%worm_all_limit = helium_env(k)%helium%worm_all_limit+itmp
CALL section_vals_val_get(helium_section, "WORM%STAGING_WEIGHT", &
i_val=itmp)
helium_env(k)%helium%worm_staging_min = helium_env(k)%helium%worm_centroid_max+1
helium_env(k)%helium%worm_staging_max = helium_env(k)%helium%worm_centroid_max+itmp
helium_env(k)%helium%worm_all_limit = helium_env(k)%helium%worm_all_limit+itmp
IF (helium_env(k)%helium%worm_allow_open) THEN
CALL section_vals_val_get(helium_section, "WORM%CRAWL_WEIGHT", &
i_val=itmp)
helium_env(k)%helium%worm_fcrawl_min = helium_env(k)%helium%worm_staging_max+1
helium_env(k)%helium%worm_fcrawl_max = helium_env(k)%helium%worm_staging_max+itmp
helium_env(k)%helium%worm_all_limit = helium_env(k)%helium%worm_all_limit+itmp
helium_env(k)%helium%worm_bcrawl_min = helium_env(k)%helium%worm_fcrawl_max+1
helium_env(k)%helium%worm_bcrawl_max = helium_env(k)%helium%worm_fcrawl_max+itmp
helium_env(k)%helium%worm_all_limit = helium_env(k)%helium%worm_all_limit+itmp
CALL section_vals_val_get(helium_section, "WORM%HEAD_TAIL_WEIGHT", &
i_val=itmp)
helium_env(k)%helium%worm_head_min = helium_env(k)%helium%worm_bcrawl_max+1
helium_env(k)%helium%worm_head_max = helium_env(k)%helium%worm_bcrawl_max+itmp
helium_env(k)%helium%worm_all_limit = helium_env(k)%helium%worm_all_limit+itmp
helium_env(k)%helium%worm_tail_min = helium_env(k)%helium%worm_head_max+1
helium_env(k)%helium%worm_tail_max = helium_env(k)%helium%worm_head_max+itmp
helium_env(k)%helium%worm_all_limit = helium_env(k)%helium%worm_all_limit+itmp
CALL section_vals_val_get(helium_section, "WORM%SWAP_WEIGHT", &
i_val=itmp)
helium_env(k)%helium%worm_swap_min = helium_env(k)%helium%worm_tail_max+1
helium_env(k)%helium%worm_swap_max = helium_env(k)%helium%worm_tail_max+itmp
helium_env(k)%helium%worm_all_limit = helium_env(k)%helium%worm_all_limit+itmp
CALL section_vals_val_get(helium_section, "WORM%OPEN_CLOSE_WEIGHT", &
i_val=itmp)
helium_env(k)%helium%worm_open_close_min = helium_env(k)%helium%worm_swap_max+1
helium_env(k)%helium%worm_open_close_max = helium_env(k)%helium%worm_swap_max+itmp
helium_env(k)%helium%worm_all_limit = helium_env(k)%helium%worm_all_limit+itmp
CALL section_vals_val_get(helium_section, "WORM%CRAWL_REPETION", &
i_val=helium_env(k)%helium%worm_repeat_crawl)
END IF

!CPPostcondition(i<helium_env(k)%helium%beads,cp_failure_level,routineP,failure)
! end of worm
CASE DEFAULT
Expand Down
14 changes: 12 additions & 2 deletions src/motion/helium_types.F
Original file line number Diff line number Diff line change
Expand Up @@ -125,10 +125,20 @@ MODULE helium_types
INTEGER :: sampling_method
! worm sampling parameters
REAL(KIND=dp) :: worm_centroid_drmax
INTEGER :: worm_staging_l
INTEGER :: worm_staging_l
INTEGER :: worm_repeat_crawl
INTEGER :: worm_all_limit
INTEGER :: worm_centroid_min, worm_centroid_max
INTEGER :: worm_staging_min, worm_staging_max
INTEGER :: worm_fcrawl_min, worm_fcrawl_max
INTEGER :: worm_bcrawl_min, worm_bcrawl_max
INTEGER :: worm_head_min, worm_head_max
INTEGER :: worm_tail_min, worm_tail_max
INTEGER :: worm_swap_min, worm_swap_max
INTEGER :: worm_open_close_min, worm_open_close_max
REAL(KIND=dp) :: worm_open_close_scale
REAL(KIND=dp) :: worm_ln_openclose_scale
LOGICAL :: worm_allow_open, worm_show_statistics
LOGICAL :: worm_allow_open, worm_show_statistics
! worm specific variables
REAL(KIND=dp), DIMENSION(3) :: worm_xtra_bead, worm_xtra_bead_work
Expand Down
69 changes: 33 additions & 36 deletions src/motion/helium_worm.F
Original file line number Diff line number Diff line change
Expand Up @@ -91,17 +91,18 @@ SUBROUTINE helium_sample_worm(helium, pint_env)
IF (helium%worm_allow_open) THEN
DO ! Exit criterion at the end of the loop
DO iMC = 1, nMC
imove = next_random_number(helium%rng_stream_uniform, 1, helium%worm_all_limit)
IF (helium%worm_is_closed) THEN
imove = next_random_number(helium%rng_stream_uniform, 1, 100)
SELECT CASE (imove)
CASE (1:10)
IF ((imove >= helium%worm_centroid_min) .AND. (imove <= helium%worm_centroid_max)) THEN
! centroid move
iatom = next_random_number(helium%rng_stream_uniform, 1, helium%atoms)
CALL worm_centroid_move(pint_env, helium, &
iatom, helium%worm_centroid_drmax, ac)
ncentratt = ncentratt+1
ncentracc = ncentracc+ac
CASE (11:90)
! Note: weights for open and centroid move are taken from open sampling
! staging is adjusted to conserve these weights
ELSE IF ((imove >= helium%worm_centroid_max+1) .AND. (imove <= helium%worm_open_close_min-1)) THEN
! staging move
iatom = next_random_number(helium%rng_stream_uniform, 1, helium%atoms)
ibead = next_random_number(helium%rng_stream_uniform, 1, helium%beads)
Expand All @@ -110,7 +111,7 @@ SUBROUTINE helium_sample_worm(helium, pint_env)
iatom, ibead, staging_l, ac)
nstagatt = nstagatt+1
nstagacc = nstagacc+ac
CASE (91:100)
ELSE IF ((imove >= helium%worm_open_close_min) .AND. (imove <= helium%worm_open_close_max)) THEN
! attempt opening of worm
iatom = next_random_number(helium%rng_stream_uniform, 1, helium%atoms)
ibead = next_random_number(helium%rng_stream_uniform, 1, helium%beads)
Expand All @@ -119,22 +120,19 @@ SUBROUTINE helium_sample_worm(helium, pint_env)
iatom, ibead, staging_l, ac)
nopenatt = nopenatt+1
nopenacc = nopenacc+ac
CASE DEFAULT
! this must not occour
ELSE
! this must not occur
CPABORT("Undefined move selected in helium worm sampling!")
END SELECT
END IF
ELSE ! worm is open
imove = next_random_number(helium%rng_stream_uniform, 1, 100)
SELECT CASE (imove)
CASE (1:10)
IF ((imove >= helium%worm_centroid_min) .AND. (imove <= helium%worm_centroid_max)) THEN
! centroid move
iatom = next_random_number(helium%rng_stream_uniform, 1, helium%atoms)
CALL worm_centroid_move(pint_env, helium, &
iatom, helium%worm_centroid_drmax, ac)
ncentratt = ncentratt+1
ncentracc = ncentracc+ac
!CASE (11:60)
CASE (11:40)
ELSE IF ((imove >= helium%worm_staging_min) .AND. (imove <= helium%worm_staging_max)) THEN
! staging move
iatom = next_random_number(helium%rng_stream_uniform, 1, helium%atoms)
ibead = next_random_number(helium%rng_stream_uniform, 1, helium%beads)
Expand All @@ -143,56 +141,56 @@ SUBROUTINE helium_sample_worm(helium, pint_env)
iatom, ibead, staging_l, ac)
nstagatt = nstagatt+1
nstagacc = nstagacc+ac
CASE (41:50)
ELSE IF ((imove >= helium%worm_fcrawl_min) .AND. (imove <= helium%worm_fcrawl_max)) THEN
! crawl forward
DO icrawl = 1, 5
DO icrawl = 1, helium%worm_repeat_crawl
staging_l = next_random_number(helium%rng_stream_uniform, 2, helium%worm_staging_l)
CALL worm_crawl_move_forward(pint_env, helium, &
staging_l, ac)
ncrawlfwdatt = ncrawlfwdatt+1
ncrawlfwdacc = ncrawlfwdacc+ac
END DO
CASE (51:60)
ELSE IF ((imove >= helium%worm_bcrawl_min) .AND. (imove <= helium%worm_bcrawl_max)) THEN
! crawl backward
DO icrawl = 1, 5
DO icrawl = 1, helium%worm_repeat_crawl
staging_l = next_random_number(helium%rng_stream_uniform, 2, helium%worm_staging_l)
CALL worm_crawl_move_backward(pint_env, helium, &
staging_l, ac)
ncrawlbwdatt = ncrawlbwdatt+1
ncrawlbwdacc = ncrawlbwdacc+ac
END DO
CASE (61:70)
ELSE IF ((imove >= helium%worm_head_min) .AND. (imove <= helium%worm_head_max)) THEN
! move head
staging_l = next_random_number(helium%rng_stream_uniform, 2, helium%worm_staging_l)
CALL worm_head_move(pint_env, helium, &
staging_l, ac)
nmoveheadatt = nmoveheadatt+1
nmoveheadacc = nmoveheadacc+ac
CASE (71:80)
ELSE IF ((imove >= helium%worm_tail_min) .AND. (imove <= helium%worm_tail_max)) THEN
! move tail
staging_l = next_random_number(helium%rng_stream_uniform, 2, helium%worm_staging_l)
CALL worm_tail_move(pint_env, helium, &
staging_l, ac)
nmovetailatt = nmovetailatt+1
nmovetailacc = nmovetailacc+ac
CASE (81:90)
! attempt closing of worm
staging_l = next_random_number(helium%rng_stream_uniform, 2, helium%worm_staging_l)
CALL worm_close_move(pint_env, helium, &
staging_l, ac)
ncloseatt = ncloseatt+1
ncloseacc = ncloseacc+ac
CASE (91:100)
ELSE IF ((imove >= helium%worm_swap_min) .AND. (imove <= helium%worm_swap_max)) THEN
staging_l = next_random_number(helium%rng_stream_uniform, 2, helium%worm_staging_l)
CALL worm_swap_move(pint_env, helium, &
helium%atoms, staging_l, ac)
npswapacc = npswapacc+ac
nswapacc = nswapacc+ac
nswapatt = nswapatt+1
CASE DEFAULT
! this must not occour
ELSE IF ((imove >= helium%worm_open_close_min) .AND. (imove <= helium%worm_open_close_max)) THEN
! attempt closing of worm
staging_l = next_random_number(helium%rng_stream_uniform, 2, helium%worm_staging_l)
CALL worm_close_move(pint_env, helium, &
staging_l, ac)
ncloseatt = ncloseatt+1
ncloseacc = ncloseacc+ac
ELSE
! this must not occur
CPABORT("Undefined move selected in helium worm sampling!")
END SELECT
END IF
END IF

! Accumulate statistics if we are in the Z-sector:
Expand Down Expand Up @@ -222,28 +220,27 @@ SUBROUTINE helium_sample_worm(helium, pint_env)
END DO !attempts loop
ELSE ! only closed configurations allowed
DO iMC = 1, nMC
imove = next_random_number(helium%rng_stream_uniform, 1, 2)
imove = next_random_number(helium%rng_stream_uniform, 1, helium%worm_all_limit)

SELECT CASE (imove)
CASE (1)
IF ((imove >= helium%worm_centroid_min) .AND. (imove <= helium%worm_centroid_max)) THEN
! centroid move
iatom = next_random_number(helium%rng_stream_uniform, 1, helium%atoms)
CALL worm_centroid_move(pint_env, helium, &
iatom, helium%worm_centroid_drmax, ac)
ncentratt = ncentratt+1
ncentracc = ncentracc+ac
CASE (2)
ELSE IF ((imove >= helium%worm_staging_min) .AND. (imove <= helium%worm_staging_max)) THEN
! staging move
iatom = next_random_number(helium%rng_stream_uniform, 1, helium%atoms)
ibead = next_random_number(helium%rng_stream_uniform, 1, helium%beads)
CALL worm_staging_move(pint_env, helium, &
iatom, ibead, helium%worm_staging_l, ac)
nstagatt = nstagatt+1
nstagacc = nstagacc+ac
CASE DEFAULT
ELSE
! this must not occour
CPABORT("Undefined move selected in helium worm sampling!")
END SELECT
END IF

! Accumulate statistics if we are in closed configurations (which we always are)
nstat = nstat+1
Expand Down
50 changes: 46 additions & 4 deletions src/start/input_cp2k_motion.F
Original file line number Diff line number Diff line change
Expand Up @@ -2279,8 +2279,8 @@ SUBROUTINE create_helium_section(section)
! worm algorithm parameters:
NULLIFY (subsection)
CALL section_create(subsection, __LOCATION__, name="WORM", &
description="Enables sampling with the worm algorithm by Bonisegni", &
n_keywords=5, n_subsections=0, repeats=.FALSE.)
description="Enables sampling via the canonical worm algorithm adapted from Bonisegni", &
n_keywords=11, n_subsections=0, repeats=.FALSE.)

CALL keyword_create(keyword, __LOCATION__, name="CENTROID_DRMAX", &
description="Maximum displacement allowed for the centroid moves", &
Expand All @@ -2289,7 +2289,7 @@ SUBROUTINE create_helium_section(section)
CALL keyword_release(keyword)

CALL keyword_create(keyword, __LOCATION__, name="STAGING_L", &
description="L-1 beads will be moved", &
description="From 2 up to max. L-1 beads will be moved", &
repeats=.FALSE., default_i_val=5)
CALL section_add_keyword(subsection, keyword)
CALL keyword_release(keyword)
Expand All @@ -2301,7 +2301,7 @@ SUBROUTINE create_helium_section(section)
CALL keyword_release(keyword)

CALL keyword_create(keyword, __LOCATION__, name="ALLOW_OPEN", &
description="Enable boltzmanonic sampling", &
description="Enable bosonic exchange sampling", &
repeats=.FALSE., default_l_val=.TRUE.)
CALL section_add_keyword(subsection, keyword)
CALL keyword_release(keyword)
Expand All @@ -2312,6 +2312,48 @@ SUBROUTINE create_helium_section(section)
CALL section_add_keyword(subsection, keyword)
CALL keyword_release(keyword)

CALL keyword_create(keyword, __LOCATION__, name="CENTROID_WEIGHT", &
description="Absolute weight of the centroid move", &
repeats=.FALSE., default_i_val=10)
CALL section_add_keyword(subsection, keyword)
CALL keyword_release(keyword)

CALL keyword_create(keyword, __LOCATION__, name="STAGING_WEIGHT", &
description="Absolute weight of the staging move", &
repeats=.FALSE., default_i_val=30)
CALL section_add_keyword(subsection, keyword)
CALL keyword_release(keyword)

CALL keyword_create(keyword, __LOCATION__, name="OPEN_CLOSE_WEIGHT", &
description="Absolute weight of the open/close move", &
repeats=.FALSE., default_i_val=10)
CALL section_add_keyword(subsection, keyword)
CALL keyword_release(keyword)

CALL keyword_create(keyword, __LOCATION__, name="HEAD_TAIL_WEIGHT", &
description="Absolute weight of the head/tail moves (both)", &
repeats=.FALSE., default_i_val=10)
CALL section_add_keyword(subsection, keyword)
CALL keyword_release(keyword)

CALL keyword_create(keyword, __LOCATION__, name="CRAWL_WEIGHT", &
description="Absolute weight of the crawl bwd/fwd moves (both)", &
repeats=.FALSE., default_i_val=10)
CALL section_add_keyword(subsection, keyword)
CALL keyword_release(keyword)

CALL keyword_create(keyword, __LOCATION__, name="CRAWL_REPETION", &
description="Number of repeated crawl moves", &
repeats=.FALSE., default_i_val=4)
CALL section_add_keyword(subsection, keyword)
CALL keyword_release(keyword)

CALL keyword_create(keyword, __LOCATION__, name="SWAP_WEIGHT", &
description="Absolute weight of the crawl move", &
repeats=.FALSE., default_i_val=10)
CALL section_add_keyword(subsection, keyword)
CALL keyword_release(keyword)

CALL section_add_subsection(section, subsection)
CALL section_release(subsection) ! release WORM subsection

Expand Down

0 comments on commit 1cbdb02

Please sign in to comment.