From 01064014fd19087852a4f8e704206ca8254bdc05 Mon Sep 17 00:00:00 2001 From: Joost VandeVondele Date: Tue, 11 Nov 2014 14:26:23 +0000 Subject: [PATCH] fix crashing analytic 1d solver, add some tests. svn-origin-rev: 14602 --- src/common/bessel_lib.F | 4 +-- src/pw/pw_poisson_types.F | 21 +++++++++----- tests/QS/regtest-ot/He_a_x.inp | 48 ++++++++++++++++++++++++++++++++ tests/QS/regtest-ot/He_a_xy.inp | 48 ++++++++++++++++++++++++++++++++ tests/QS/regtest-ot/He_a_xyz.inp | 48 ++++++++++++++++++++++++++++++++ tests/QS/regtest-ot/He_a_xz.inp | 48 ++++++++++++++++++++++++++++++++ tests/QS/regtest-ot/He_a_y.inp | 48 ++++++++++++++++++++++++++++++++ tests/QS/regtest-ot/He_a_yz.inp | 48 ++++++++++++++++++++++++++++++++ tests/QS/regtest-ot/He_a_z.inp | 48 ++++++++++++++++++++++++++++++++ tests/QS/regtest-ot/TEST_FILES | 8 ++++++ 10 files changed, 360 insertions(+), 9 deletions(-) create mode 100644 tests/QS/regtest-ot/He_a_x.inp create mode 100644 tests/QS/regtest-ot/He_a_xy.inp create mode 100644 tests/QS/regtest-ot/He_a_xyz.inp create mode 100644 tests/QS/regtest-ot/He_a_xz.inp create mode 100644 tests/QS/regtest-ot/He_a_y.inp create mode 100644 tests/QS/regtest-ot/He_a_yz.inp create mode 100644 tests/QS/regtest-ot/He_a_z.inp diff --git a/src/common/bessel_lib.F b/src/common/bessel_lib.F index dc81ae8df3..3445c1deaa 100644 --- a/src/common/bessel_lib.F +++ b/src/common/bessel_lib.F @@ -105,7 +105,7 @@ CONTAINS ! ***************************************************************************** !> \brief ... -!> \param x ... +!> \param x must be positive ! ***************************************************************************** FUNCTION bessk0 ( x ) @@ -134,7 +134,7 @@ CONTAINS ! ***************************************************************************** !> \brief ... -!> \param x ... +!> \param x must be positive ! ***************************************************************************** FUNCTION bessk1 ( x ) diff --git a/src/pw/pw_poisson_types.F b/src/pw/pw_poisson_types.F index a9daee46d5..6e2373f2da 100644 --- a/src/pw/pw_poisson_types.F +++ b/src/pw/pw_poisson_types.F @@ -182,8 +182,9 @@ CONTAINS INTEGER :: dim, i, ig, iz, n, nz, stat LOGICAL :: failure - REAL(KIND=dp) :: g2, g3d, gg, gxy, j0g, j1g, & - k0g, k1g, rlength, zlength + REAL(KIND=dp) :: g2, g3d, gg, gxy, gz, j0g, & + j1g, k0g, k1g, rlength, & + zlength REAL(KIND=dp), DIMENSION(3) :: abc TYPE(pw_grid_type), POINTER :: grid TYPE(pw_type), POINTER :: gf @@ -350,7 +351,7 @@ CONTAINS IF ( grid % have_g0 ) gf % cc ( 1 ) = 0.0_dp CASE ( ANALYTIC1D ) - + ! see 'ab initio molecular dynamics' table 3.1 ! iz is the direction of the PBC ( can be 1,2,3 -> x,y,z ) iz = green % special_dimension ! rlength is the radius of the tube @@ -359,12 +360,18 @@ CONTAINS g2 = grid % gsq ( ig ) g3d = fourpi / g2 gxy = SQRT ( g2 - grid % g(iz,ig) * grid % g(iz,ig) ) + gz = ABS(grid % g(iz,ig)) j0g = bessj0 ( rlength * gxy ) j1g = bessj1 ( rlength * gxy ) - k0g = bessk0 ( rlength * grid % g(iz,ig) ) - k1g = bessk1 ( rlength * grid % g(iz,ig) ) - gf % cc ( ig ) = g3d * ( 1.0_dp - rlength * & - ( gxy * j1g * k0g - grid % g(iz,ig) * j0g * k1g ) ) + IF (gz>0) THEN + k0g = bessk0 ( rlength * gz ) + k1g = bessk1 ( rlength * gz ) + ELSE + k0g = 0 + k1g = 0 + ENDIF + gf % cc ( ig ) = g3d * ( 1.0_dp + rlength * & + ( gxy * j1g * k0g - gz * j0g * k1g ) ) END DO IF ( grid % have_g0 ) gf % cc ( 1 ) = 0.0_dp diff --git a/tests/QS/regtest-ot/He_a_x.inp b/tests/QS/regtest-ot/He_a_x.inp new file mode 100644 index 0000000000..cc4316a9cf --- /dev/null +++ b/tests/QS/regtest-ot/He_a_x.inp @@ -0,0 +1,48 @@ +&FORCE_EVAL + METHOD Quickstep + &DFT + BASIS_SET_FILE_NAME ../../../data/BASIS_SET + POTENTIAL_FILE_NAME ../../../data/POTENTIAL + &MGRID + CUTOFF 100 + &END MGRID + &QS + EPS_DEFAULT 1.0E-12 + MAP_CONSISTENT + EXTRAPOLATION PS + EXTRAPOLATION_ORDER 3 + &END QS + &POISSON + PERIODIC X + POISSON_SOLVER ANALYTIC + &END + &SCF + EPS_SCF 1.0E-6 + SCF_GUESS ATOMIC + MAX_SCF 20 + &OT + &END OT + &END SCF + &XC + &XC_FUNCTIONAL Pade + &END XC_FUNCTIONAL + &END XC + &END DFT + &SUBSYS + &CELL + ABC 6.0 6.0 6.0 + &END CELL + &COORD + He 3.0 3.0 3.0 + &END COORD + &KIND He + BASIS_SET DZVP-GTH-PADE + POTENTIAL GTH-PADE-q2 + &END KIND + &END SUBSYS +&END FORCE_EVAL +&GLOBAL + PROJECT He2_none + RUN_TYPE ENERGY + PRINT_LEVEL MEDIUM +&END GLOBAL diff --git a/tests/QS/regtest-ot/He_a_xy.inp b/tests/QS/regtest-ot/He_a_xy.inp new file mode 100644 index 0000000000..3c7f5162a6 --- /dev/null +++ b/tests/QS/regtest-ot/He_a_xy.inp @@ -0,0 +1,48 @@ +&FORCE_EVAL + METHOD Quickstep + &DFT + BASIS_SET_FILE_NAME ../../../data/BASIS_SET + POTENTIAL_FILE_NAME ../../../data/POTENTIAL + &MGRID + CUTOFF 100 + &END MGRID + &QS + EPS_DEFAULT 1.0E-12 + MAP_CONSISTENT + EXTRAPOLATION PS + EXTRAPOLATION_ORDER 3 + &END QS + &POISSON + PERIODIC XY + POISSON_SOLVER ANALYTIC + &END + &SCF + EPS_SCF 1.0E-6 + SCF_GUESS ATOMIC + MAX_SCF 20 + &OT + &END OT + &END SCF + &XC + &XC_FUNCTIONAL Pade + &END XC_FUNCTIONAL + &END XC + &END DFT + &SUBSYS + &CELL + ABC 6.0 6.0 6.0 + &END CELL + &COORD + He 3.0 3.0 3.0 + &END COORD + &KIND He + BASIS_SET DZVP-GTH-PADE + POTENTIAL GTH-PADE-q2 + &END KIND + &END SUBSYS +&END FORCE_EVAL +&GLOBAL + PROJECT He2_none + RUN_TYPE ENERGY + PRINT_LEVEL MEDIUM +&END GLOBAL diff --git a/tests/QS/regtest-ot/He_a_xyz.inp b/tests/QS/regtest-ot/He_a_xyz.inp new file mode 100644 index 0000000000..0dacb6999e --- /dev/null +++ b/tests/QS/regtest-ot/He_a_xyz.inp @@ -0,0 +1,48 @@ +&FORCE_EVAL + METHOD Quickstep + &DFT + BASIS_SET_FILE_NAME ../../../data/BASIS_SET + POTENTIAL_FILE_NAME ../../../data/POTENTIAL + &MGRID + CUTOFF 100 + &END MGRID + &QS + EPS_DEFAULT 1.0E-12 + MAP_CONSISTENT + EXTRAPOLATION PS + EXTRAPOLATION_ORDER 3 + &END QS + &POISSON + PERIODIC XYZ + POISSON_SOLVER ANALYTIC + &END + &SCF + EPS_SCF 1.0E-6 + SCF_GUESS ATOMIC + MAX_SCF 20 + &OT + &END OT + &END SCF + &XC + &XC_FUNCTIONAL Pade + &END XC_FUNCTIONAL + &END XC + &END DFT + &SUBSYS + &CELL + ABC 6.0 6.0 6.0 + &END CELL + &COORD + He 3.0 3.0 3.0 + &END COORD + &KIND He + BASIS_SET DZVP-GTH-PADE + POTENTIAL GTH-PADE-q2 + &END KIND + &END SUBSYS +&END FORCE_EVAL +&GLOBAL + PROJECT He2_none + RUN_TYPE ENERGY + PRINT_LEVEL MEDIUM +&END GLOBAL diff --git a/tests/QS/regtest-ot/He_a_xz.inp b/tests/QS/regtest-ot/He_a_xz.inp new file mode 100644 index 0000000000..c435b34d02 --- /dev/null +++ b/tests/QS/regtest-ot/He_a_xz.inp @@ -0,0 +1,48 @@ +&FORCE_EVAL + METHOD Quickstep + &DFT + BASIS_SET_FILE_NAME ../../../data/BASIS_SET + POTENTIAL_FILE_NAME ../../../data/POTENTIAL + &MGRID + CUTOFF 100 + &END MGRID + &QS + EPS_DEFAULT 1.0E-12 + MAP_CONSISTENT + EXTRAPOLATION PS + EXTRAPOLATION_ORDER 3 + &END QS + &POISSON + PERIODIC XZ + POISSON_SOLVER ANALYTIC + &END + &SCF + EPS_SCF 1.0E-6 + SCF_GUESS ATOMIC + MAX_SCF 20 + &OT + &END OT + &END SCF + &XC + &XC_FUNCTIONAL Pade + &END XC_FUNCTIONAL + &END XC + &END DFT + &SUBSYS + &CELL + ABC 6.0 6.0 6.0 + &END CELL + &COORD + He 3.0 3.0 3.0 + &END COORD + &KIND He + BASIS_SET DZVP-GTH-PADE + POTENTIAL GTH-PADE-q2 + &END KIND + &END SUBSYS +&END FORCE_EVAL +&GLOBAL + PROJECT He2_none + RUN_TYPE ENERGY + PRINT_LEVEL MEDIUM +&END GLOBAL diff --git a/tests/QS/regtest-ot/He_a_y.inp b/tests/QS/regtest-ot/He_a_y.inp new file mode 100644 index 0000000000..fd9bd7cd26 --- /dev/null +++ b/tests/QS/regtest-ot/He_a_y.inp @@ -0,0 +1,48 @@ +&FORCE_EVAL + METHOD Quickstep + &DFT + BASIS_SET_FILE_NAME ../../../data/BASIS_SET + POTENTIAL_FILE_NAME ../../../data/POTENTIAL + &MGRID + CUTOFF 100 + &END MGRID + &QS + EPS_DEFAULT 1.0E-12 + MAP_CONSISTENT + EXTRAPOLATION PS + EXTRAPOLATION_ORDER 3 + &END QS + &POISSON + PERIODIC Y + POISSON_SOLVER ANALYTIC + &END + &SCF + EPS_SCF 1.0E-6 + SCF_GUESS ATOMIC + MAX_SCF 20 + &OT + &END OT + &END SCF + &XC + &XC_FUNCTIONAL Pade + &END XC_FUNCTIONAL + &END XC + &END DFT + &SUBSYS + &CELL + ABC 6.0 6.0 6.0 + &END CELL + &COORD + He 3.0 3.0 3.0 + &END COORD + &KIND He + BASIS_SET DZVP-GTH-PADE + POTENTIAL GTH-PADE-q2 + &END KIND + &END SUBSYS +&END FORCE_EVAL +&GLOBAL + PROJECT He2_none + RUN_TYPE ENERGY + PRINT_LEVEL MEDIUM +&END GLOBAL diff --git a/tests/QS/regtest-ot/He_a_yz.inp b/tests/QS/regtest-ot/He_a_yz.inp new file mode 100644 index 0000000000..0518e21944 --- /dev/null +++ b/tests/QS/regtest-ot/He_a_yz.inp @@ -0,0 +1,48 @@ +&FORCE_EVAL + METHOD Quickstep + &DFT + BASIS_SET_FILE_NAME ../../../data/BASIS_SET + POTENTIAL_FILE_NAME ../../../data/POTENTIAL + &MGRID + CUTOFF 100 + &END MGRID + &QS + EPS_DEFAULT 1.0E-12 + MAP_CONSISTENT + EXTRAPOLATION PS + EXTRAPOLATION_ORDER 3 + &END QS + &POISSON + PERIODIC YZ + POISSON_SOLVER ANALYTIC + &END + &SCF + EPS_SCF 1.0E-6 + SCF_GUESS ATOMIC + MAX_SCF 20 + &OT + &END OT + &END SCF + &XC + &XC_FUNCTIONAL Pade + &END XC_FUNCTIONAL + &END XC + &END DFT + &SUBSYS + &CELL + ABC 6.0 6.0 6.0 + &END CELL + &COORD + He 3.0 3.0 3.0 + &END COORD + &KIND He + BASIS_SET DZVP-GTH-PADE + POTENTIAL GTH-PADE-q2 + &END KIND + &END SUBSYS +&END FORCE_EVAL +&GLOBAL + PROJECT He2_none + RUN_TYPE ENERGY + PRINT_LEVEL MEDIUM +&END GLOBAL diff --git a/tests/QS/regtest-ot/He_a_z.inp b/tests/QS/regtest-ot/He_a_z.inp new file mode 100644 index 0000000000..39d63ed4c8 --- /dev/null +++ b/tests/QS/regtest-ot/He_a_z.inp @@ -0,0 +1,48 @@ +&FORCE_EVAL + METHOD Quickstep + &DFT + BASIS_SET_FILE_NAME ../../../data/BASIS_SET + POTENTIAL_FILE_NAME ../../../data/POTENTIAL + &MGRID + CUTOFF 100 + &END MGRID + &QS + EPS_DEFAULT 1.0E-12 + MAP_CONSISTENT + EXTRAPOLATION PS + EXTRAPOLATION_ORDER 3 + &END QS + &POISSON + PERIODIC Z + POISSON_SOLVER ANALYTIC + &END + &SCF + EPS_SCF 1.0E-6 + SCF_GUESS ATOMIC + MAX_SCF 20 + &OT + &END OT + &END SCF + &XC + &XC_FUNCTIONAL Pade + &END XC_FUNCTIONAL + &END XC + &END DFT + &SUBSYS + &CELL + ABC 6.0 6.0 6.0 + &END CELL + &COORD + He 3.0 3.0 3.0 + &END COORD + &KIND He + BASIS_SET DZVP-GTH-PADE + POTENTIAL GTH-PADE-q2 + &END KIND + &END SUBSYS +&END FORCE_EVAL +&GLOBAL + PROJECT He2_none + RUN_TYPE ENERGY + PRINT_LEVEL MEDIUM +&END GLOBAL diff --git a/tests/QS/regtest-ot/TEST_FILES b/tests/QS/regtest-ot/TEST_FILES index 0a509472e0..b37f1872ba 100644 --- a/tests/QS/regtest-ot/TEST_FILES +++ b/tests/QS/regtest-ot/TEST_FILES @@ -35,3 +35,11 @@ H2-diffBECKE-ET_coupling.inp 1 1e-13 sic_energy.inp 1 # elf C2H4-elf.inp 1 +# analytic poisson solver +He_a_xyz.inp 1 +He_a_xz.inp 1 +He_a_yz.inp 1 +He_a_xy.inp 1 +He_a_x.inp 1 +He_a_y.inp 1 +He_a_z.inp 1