Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Add cell symmetry keyword hexagonal_gamma_120 #2758

Merged
merged 1 commit into from Apr 29, 2023
Merged

Add cell symmetry keyword hexagonal_gamma_120 #2758

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
18 changes: 11 additions & 7 deletions src/input_cp2k_subsys.F
View file
Open in desktop
Expand Up @@ -19,10 +19,11 @@ MODULE input_cp2k_subsys
VandeVondele2005a,&
VandeVondele2007
USE cell_types, ONLY: &
cell_sym_cubic, cell_sym_hexagonal, cell_sym_monoclinic, cell_sym_monoclinic_gamma_ab, &
cell_sym_none, cell_sym_orthorhombic, cell_sym_rhombohedral, cell_sym_tetragonal_ab, &
cell_sym_tetragonal_ac, cell_sym_tetragonal_bc, cell_sym_triclinic, use_perd_none, &
use_perd_x, use_perd_xy, use_perd_xyz, use_perd_xz, use_perd_y, use_perd_yz, use_perd_z
cell_sym_cubic, cell_sym_hexagonal_gamma_120, cell_sym_hexagonal_gamma_60, &
cell_sym_monoclinic, cell_sym_monoclinic_gamma_ab, cell_sym_none, cell_sym_orthorhombic, &
cell_sym_rhombohedral, cell_sym_tetragonal_ab, cell_sym_tetragonal_ac, &
cell_sym_tetragonal_bc, cell_sym_triclinic, use_perd_none, use_perd_x, use_perd_xy, &
use_perd_xyz, use_perd_xz, use_perd_y, use_perd_yz, use_perd_z
USE cp_output_handling, ONLY: cp_print_key_section_create,&
debug_print_level,&
high_print_level,&
Expand Down Expand Up @@ -211,22 +212,25 @@ SUBROUTINE create_cell_section_low(section, periodic)
usage="SYMMETRY monoclinic", &
enum_desc=s2a("No cell symmetry", &
"Triclinic (a ≠ b ≠ c ≠ a, α ≠ β ≠ γ ≠ α ≠ 90°)", &
"Monoclinic (a ≠ b ≠ c ≠ a, α = γ = 90°, β ≠ 90°)", &
"Monoclinic (a ≠ b ≠ c, α = γ = 90°, β ≠ 90°)", &
"Monoclinic (a = b ≠ c, α = β = 90°, γ ≠ 90°)", &
"Orthorhombic (a ≠ b ≠ c, α = β = γ = 90°)", &
"Tetragonal (a = b ≠ c, α = β = γ = 90°)", &
"Tetragonal (a = c ≠ b, α = β = γ = 90°)", &
"Tetragonal (a ≠ b = c, α = β = γ = 90°)", &
"Tetragonal (alias for TETRAGONAL_AB)", &
"Rhombohedral (a = b = c, α = β = γ ≠ 90°)", &
"Hexagonal (alias for HEXAGONAL_GAMMA_60)", &
"Hexagonal (a = b ≠ c, α = β = 90°, γ = 60°)", &
"Hexagonal (a = b ≠ c, α = β = 90°, γ = 120°)", &
"Cubic (a = b = c, α = β = γ = 90°)"), &
enum_c_vals=s2a("NONE", "TRICLINIC", "MONOCLINIC", "MONOCLINIC_GAMMA_AB", "ORTHORHOMBIC", &
"TETRAGONAL_AB", "TETRAGONAL_AC", "TETRAGONAL_BC", "TETRAGONAL", "RHOMBOHEDRAL", &
"HEXAGONAL", "CUBIC"), &
"HEXAGONAL", "HEXAGONAL_GAMMA_60", "HEXAGONAL_GAMMA_120", "CUBIC"), &
enum_i_vals=(/cell_sym_none, cell_sym_triclinic, cell_sym_monoclinic, cell_sym_monoclinic_gamma_ab, &
cell_sym_orthorhombic, cell_sym_tetragonal_ab, cell_sym_tetragonal_ac, cell_sym_tetragonal_bc, &
cell_sym_tetragonal_ab, cell_sym_rhombohedral, cell_sym_hexagonal, cell_sym_cubic/), &
cell_sym_tetragonal_ab, cell_sym_rhombohedral, cell_sym_hexagonal_gamma_60, &
cell_sym_hexagonal_gamma_60, cell_sym_hexagonal_gamma_120, cell_sym_cubic/), &
default_i_val=cell_sym_none)
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
Expand Down
35 changes: 21 additions & 14 deletions src/motion/cell_opt_utils.F
View file
Open in desktop
Expand Up @@ -13,9 +13,10 @@
! **************************************************************************************************
MODULE cell_opt_utils
USE cell_types, ONLY: &
cell_sym_cubic, cell_sym_hexagonal, cell_sym_monoclinic, cell_sym_monoclinic_gamma_ab, &
cell_sym_orthorhombic, cell_sym_rhombohedral, cell_sym_tetragonal_ab, &
cell_sym_tetragonal_ac, cell_sym_tetragonal_bc, cell_sym_triclinic, cell_type
cell_sym_cubic, cell_sym_hexagonal_gamma_120, cell_sym_hexagonal_gamma_60, &
cell_sym_monoclinic, cell_sym_monoclinic_gamma_ab, cell_sym_orthorhombic, &
cell_sym_rhombohedral, cell_sym_tetragonal_ab, cell_sym_tetragonal_ac, &
cell_sym_tetragonal_bc, cell_sym_triclinic, cell_type
USE cp_files, ONLY: close_file,&
open_file
USE cp_log_handling, ONLY: cp_get_default_logger,&
Expand Down Expand Up @@ -320,9 +321,8 @@ SUBROUTINE apply_cell_constraints(gradient, cell, keep_angles, keep_symmetry, co
INTEGER, INTENT(IN) :: constraint_id

REAL(KIND=dp) :: a, a_length, ab_length, b_length, cosa, &
cosah, cosgamma, deriv_gamma, g, &
gamma, norm, norm_b, norm_c, sina, &
sinah, singamma
cosah, cosg, deriv_gamma, g, gamma, &
norm, norm_b, norm_c, sina, sinah, sing

IF (keep_angles) THEN
! If we want to keep the angles constant we have to project out the
Expand Down Expand Up @@ -378,13 +378,20 @@ SUBROUTINE apply_cell_constraints(gradient, cell, keep_angles, keep_symmetry, co
gradient(2) = 0.0_dp
gradient(4) = 0.0_dp
gradient(5) = 0.0_dp
CASE (cell_sym_hexagonal)
CASE (cell_sym_hexagonal_gamma_60)
g = 0.5_dp*(gradient(1) + 0.5_dp*(gradient(2) + sqrt3*gradient(3)))
gradient(1) = g
gradient(2) = 0.5_dp*g
gradient(3) = sqrt3*gradient(2)
gradient(4) = 0.0_dp
gradient(5) = 0.0_dp
CASE (cell_sym_hexagonal_gamma_120)
g = 0.5_dp*(gradient(1) - 0.5_dp*(gradient(2) - sqrt3*gradient(3)))
gradient(1) = g
gradient(2) = -0.5_dp*g
gradient(3) = -sqrt3*gradient(2)
gradient(4) = 0.0_dp
gradient(5) = 0.0_dp
CASE (cell_sym_rhombohedral)
a = (angle(cell%hmat(:, 3), cell%hmat(:, 2)) + &
angle(cell%hmat(:, 1), cell%hmat(:, 3)) + &
Expand All @@ -407,7 +414,7 @@ SUBROUTINE apply_cell_constraints(gradient, cell, keep_angles, keep_symmetry, co
gradient(2) = 0.0_dp
gradient(5) = 0.0_dp
CASE (cell_sym_monoclinic_gamma_ab)
! Cell symmetry with a=b, alpha=beta=90deg and gammma unequal 90deg
! Cell symmetry with a = b, alpha = beta = 90 degree and gammma not equal 90 degree
a_length = SQRT(cell%hmat(1, 1)*cell%hmat(1, 1) + &
cell%hmat(2, 1)*cell%hmat(2, 1) + &
cell%hmat(3, 1)*cell%hmat(3, 1))
Expand All @@ -416,14 +423,14 @@ SUBROUTINE apply_cell_constraints(gradient, cell, keep_angles, keep_symmetry, co
cell%hmat(3, 2)*cell%hmat(3, 2))
ab_length = 0.5_dp*(a_length + b_length)
gamma = angle(cell%hmat(:, 1), cell%hmat(:, 2))
cosgamma = COS(gamma)
singamma = SIN(gamma)
cosg = COS(gamma)
sing = SIN(gamma)
! Here, g is the average derivative of the cell vector length ab_length, and deriv_gamma is the derivative of the angle gamma
g = 0.5_dp*(gradient(1) + cosgamma*gradient(2) + singamma*gradient(3))
deriv_gamma = (gradient(3)*cosgamma - gradient(2)*singamma)/b_length
g = 0.5_dp*(gradient(1) + cosg*gradient(2) + sing*gradient(3))
deriv_gamma = (gradient(3)*cosg - gradient(2)*sing)/b_length
gradient(1) = g
gradient(2) = g*cosgamma - ab_length*singamma*deriv_gamma
gradient(3) = g*singamma + ab_length*cosgamma*deriv_gamma
gradient(2) = g*cosg - ab_length*sing*deriv_gamma
gradient(3) = g*sing + ab_length*cosg*deriv_gamma
gradient(4) = 0.0_dp
gradient(5) = 0.0_dp
CASE (cell_sym_triclinic)
Expand Down
32 changes: 17 additions & 15 deletions src/subsys/cell_types.F
View file
Open in desktop
Expand Up @@ -39,8 +39,9 @@ MODULE cell_types
cell_sym_tetragonal_ac = 6, &
cell_sym_tetragonal_bc = 7, &
cell_sym_rhombohedral = 8, &
cell_sym_hexagonal = 9, &
cell_sym_cubic = 10
cell_sym_hexagonal_gamma_60 = 9, &
cell_sym_hexagonal_gamma_120 = 10, &
cell_sym_cubic = 11

INTEGER, PARAMETER, PUBLIC :: use_perd_x = 0, &
use_perd_y = 1, &
Expand Down Expand Up @@ -318,9 +319,9 @@ SUBROUTINE init_cell(cell, hmat, periodic)
REAL(KIND=dp), PARAMETER :: eps_hmat = 1.0E-14_dp

INTEGER :: dim
REAL(KIND=dp) :: a, acosa, acosah, acosgamma, alpha, &
asina, asinah, asingamma, beta, gamma, &
norm, norm_c
REAL(KIND=dp) :: a, acosa, acosah, acosg, alpha, asina, &
asinah, asing, beta, gamma, norm, &
norm_c
REAL(KIND=dp), DIMENSION(3) :: abc

IF (PRESENT(hmat)) cell%hmat(:, :) = hmat(:, :)
Expand Down Expand Up @@ -367,13 +368,14 @@ SUBROUTINE init_cell(cell, hmat, periodic)
cell%hmat(1, 1) = abc(1); cell%hmat(1, 2) = 0.0_dp; cell%hmat(1, 3) = 0.0_dp
cell%hmat(2, 1) = 0.0_dp; cell%hmat(2, 2) = abc(2); cell%hmat(2, 3) = 0.0_dp
cell%hmat(3, 1) = 0.0_dp; cell%hmat(3, 2) = 0.0_dp; cell%hmat(3, 3) = abc(3)
CASE (cell_sym_hexagonal)
CASE (cell_sym_hexagonal_gamma_60, cell_sym_hexagonal_gamma_120)
CALL get_cell(cell=cell, abc=abc)
a = 0.5_dp*(abc(1) + abc(2))
acosa = 0.5_dp*a
asina = sqrt3*acosa
cell%hmat(1, 1) = a; cell%hmat(1, 2) = acosa; cell%hmat(1, 3) = 0.0_dp
cell%hmat(2, 1) = 0.0_dp; cell%hmat(2, 2) = asina; cell%hmat(2, 3) = 0.0_dp
acosg = 0.5_dp*a
asing = sqrt3*acosg
IF (cell%symmetry_id == cell_sym_hexagonal_gamma_120) acosg = -acosg
cell%hmat(1, 1) = a; cell%hmat(1, 2) = acosg; cell%hmat(1, 3) = 0.0_dp
cell%hmat(2, 1) = 0.0_dp; cell%hmat(2, 2) = asing; cell%hmat(2, 3) = 0.0_dp
cell%hmat(3, 1) = 0.0_dp; cell%hmat(3, 2) = 0.0_dp; cell%hmat(3, 3) = abc(3)
CASE (cell_sym_rhombohedral)
CALL get_cell(cell=cell, abc=abc)
Expand All @@ -397,14 +399,14 @@ SUBROUTINE init_cell(cell, hmat, periodic)
cell%hmat(2, 1) = 0.0_dp; cell%hmat(2, 2) = abc(2); cell%hmat(2, 3) = 0.0_dp
cell%hmat(3, 1) = 0.0_dp; cell%hmat(3, 2) = 0.0_dp; cell%hmat(3, 3) = abc(3)*SIN(beta)
CASE (cell_sym_monoclinic_gamma_ab)
! Cell symmetry with a=b, alpha=beta=90deg and gammma unequal 90deg
! Cell symmetry with a = b, alpha = beta = 90 degree and gammma not equal 90 degree
CALL get_cell(cell=cell, abc=abc)
a = 0.5_dp*(abc(1) + abc(2))
gamma = angle(cell%hmat(:, 1), cell%hmat(:, 2))
acosgamma = a*COS(gamma)
asingamma = a*SIN(gamma)
cell%hmat(1, 1) = a; cell%hmat(1, 2) = acosgamma; cell%hmat(1, 3) = 0.0_dp
cell%hmat(2, 1) = 0.0_dp; cell%hmat(2, 2) = asingamma; cell%hmat(2, 3) = 0.0_dp
acosg = a*COS(gamma)
asing = a*SIN(gamma)
cell%hmat(1, 1) = a; cell%hmat(1, 2) = acosg; cell%hmat(1, 3) = 0.0_dp
cell%hmat(2, 1) = 0.0_dp; cell%hmat(2, 2) = asing; cell%hmat(2, 3) = 0.0_dp
cell%hmat(3, 1) = 0.0_dp; cell%hmat(3, 2) = 0.0_dp; cell%hmat(3, 3) = abc(3)
CASE (cell_sym_triclinic)
! Nothing to do
Expand Down
103 changes: 103 additions & 0 deletions tests/Fist/regtest-13/Si_nosym_hexagonal.inp
View file
Open in desktop
@@ -0,0 +1,103 @@
@SET Sqrt3 1.7320508075688772

&FORCE_EVAL
METHOD FIST
&MM
&FORCEFIELD
&SPLINE
EPS_SPLINE 1.0E-8
&END
&NONBONDED
&TERSOFF
ATOMS Si Si
A 1.8308E3
B 4.7118E2
lambda1 2.4799
lambda2 1.7322
alpha 0.00
beta 1.0999E-6
n 7.8734E-1
c 1.0039E5
d 1.6218E1
h -5.9826E-1
lambda3 1.7322
bigR 2.85
bigD 0.15
RCUT 10.0
&END TERSOFF
&END NONBONDED
&END FORCEFIELD
&POISSON
&EWALD
EWALD_TYPE none
&END EWALD
&END POISSON
&END MM
&SUBSYS
&CELL
A 3.850/2 -3.850*${Sqrt3}/2 0
B 3.850/2 3.850*${Sqrt3}/2 0
C 0 0 6.366
MULTIPLE_UNIT_CELL 3 3 3
PERIODIC XYZ
&END CELL
&TOPOLOGY
MULTIPLE_UNIT_CELL 3 3 3
CONNECTIVITY OFF
&END TOPOLOGY
&COORD
SCALED
Si 1/3 2/3 0.063
Si 1/3 2/3 0.437
Si 2/3 1/3 0.563
Si 2/3 1/3 0.937
&END COORD
&PRINT
&ATOMIC_COORDINATES
&END
&CELL
&END
&END
&END SUBSYS
STRESS_TENSOR ANALYTICAL
&PRINT
&STRESS_TENSOR
&END STRESS_TENSOR
&END PRINT
&END FORCE_EVAL
&GLOBAL
PREFERRED_DIAG_LIBRARY ScaLAPACK
PRINT_LEVEL medium
PROJECT Si_direct_cell_opt
RUN_TYPE CELL_OPT
&END GLOBAL
&MOTION
&CELL_OPT
EXTERNAL_PRESSURE [bar] 1.0
CONSTRAINT none
KEEP_ANGLES no
KEEP_SYMMETRY no
MAX_DR 0.00001
MAX_FORCE 0.000001
MAX_ITER 400
OPTIMIZER BFGS
PRESSURE_TOLERANCE [bar] 10.0
RMS_DR 0.000003
RMS_FORCE 0.0000003
TYPE direct_cell_opt
&BFGS
TRUST_RADIUS 0.1
USE_MODEL_HESSIAN no
USE_RAT_FUN_OPT on
&END BFGS
&END CELL_OPT
&PRINT
&CELL on
&END CELL
&RESTART on
BACKUP_COPIES 0
&END RESTART
&RESTART_HISTORY off
&END RESTART_HISTORY
&END PRINT
&END MOTION
100 changes: 100 additions & 0 deletions tests/Fist/regtest-13/Si_nosym_hexagonal_abc.inp
View file
Open in desktop
@@ -0,0 +1,100 @@
&FORCE_EVAL
METHOD FIST
&MM
&FORCEFIELD
&SPLINE
EPS_SPLINE 1.0E-8
&END
&NONBONDED
&TERSOFF
ATOMS Si Si
A 1.8308E3
B 4.7118E2
lambda1 2.4799
lambda2 1.7322
alpha 0.00
beta 1.0999E-6
n 7.8734E-1
c 1.0039E5
d 1.6218E1
h -5.9826E-1
lambda3 1.7322
bigR 2.85
bigD 0.15
RCUT 10.0
&END TERSOFF
&END NONBONDED
&END FORCEFIELD
&POISSON
&EWALD
EWALD_TYPE none
&END EWALD
&END POISSON
&END MM
&SUBSYS
&CELL
ABC 3.850 3.850 6.366
ANGLES 90.0 90.0 120.0
MULTIPLE_UNIT_CELL 3 3 3
PERIODIC XYZ
&END CELL
&TOPOLOGY
MULTIPLE_UNIT_CELL 3 3 3
CONNECTIVITY OFF
&END TOPOLOGY
&COORD
SCALED
Si 1/3 2/3 0.063
Si 1/3 2/3 0.437
Si 2/3 1/3 0.563
Si 2/3 1/3 0.937
&END COORD
&PRINT
&ATOMIC_COORDINATES
&END
&CELL
&END
&END
&END SUBSYS
STRESS_TENSOR ANALYTICAL
&PRINT
&STRESS_TENSOR
&END STRESS_TENSOR
&END PRINT
&END FORCE_EVAL
&GLOBAL
PREFERRED_DIAG_LIBRARY ScaLAPACK
PRINT_LEVEL medium
PROJECT Si_direct_cell_opt
RUN_TYPE CELL_OPT
&END GLOBAL
&MOTION
&CELL_OPT
EXTERNAL_PRESSURE [bar] 1.0
CONSTRAINT none
KEEP_ANGLES no
KEEP_SYMMETRY no
MAX_DR 0.00001
MAX_FORCE 0.000001
MAX_ITER 400
OPTIMIZER BFGS
PRESSURE_TOLERANCE [bar] 10.0
RMS_DR 0.000003
RMS_FORCE 0.0000003
TYPE direct_cell_opt
&BFGS
TRUST_RADIUS 0.1
USE_MODEL_HESSIAN no
USE_RAT_FUN_OPT on
&END BFGS
&END CELL_OPT
&PRINT
&CELL on
&END CELL
&RESTART on
BACKUP_COPIES 0
&END RESTART
&RESTART_HISTORY off
&END RESTART_HISTORY
&END PRINT
&END MOTION