1!--------------------------------------------------------------------------------------------------! 2! CP2K: A general program to perform molecular dynamics simulations ! 3! Copyright (C) 2000 - 2020 CP2K developers group ! 4!--------------------------------------------------------------------------------------------------! 5 6! ************************************************************************************************** 7!> \brief routines for handling splines_types 8!> \par History 9!> 2001-09-21-HAF added this doc entry and changed formatting 10!> \author various 11! ************************************************************************************************** 12MODULE splines_types 13 14 USE kinds, ONLY: dp 15#include "./base/base_uses.f90" 16 17 IMPLICIT NONE 18 19 PRIVATE 20 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'splines_types' 21 INTEGER, PRIVATE, SAVE :: last_spline_env_id_nr = 0 22 INTEGER, PRIVATE, SAVE :: last_spline_data_id_nr = 0 23 24 PUBLIC :: spline_env_release, spline_environment_type 25 PUBLIC :: spline_env_create, spline_data_p_type 26 PUBLIC :: spline_data_create, spline_data_p_copy 27 PUBLIC :: spline_data_retain, spline_data_p_retain 28 PUBLIC :: spline_data_release, spline_data_p_release 29 PUBLIC :: spline_factor_copy, spline_factor_create, spline_factor_release 30 PUBLIC :: spline_data_type ! the data structure for spline table 31 PUBLIC :: spline_factor_type ! the multiplicative factors for splines 32 33! ************************************************************************************************** 34!> \brief Data-structure that holds all needed information about 35!> a specific spline interpolation. 36!> \par History 37!> 2001-09-19-HAF added this doc entry and changed formatting 38!> \author unknown 39! ************************************************************************************************** 40 TYPE spline_data_type 41 INTEGER :: ref_count, id_nr 42 REAL(KIND=dp), POINTER :: y(:) ! the function values y(x) 43 REAL(KIND=dp), POINTER :: y2(:) ! the 2nd derivative via interpolation 44 INTEGER :: n ! dimension of above arrays 45 ! not used if uniform increments 46 REAL(KIND=dp) :: h ! uniform increment of x if applicable 47 REAL(KIND=dp) :: invh ! inverse of h 48 REAL(KIND=dp) :: h26 ! 1/6 * h**2 if uniform increments 49 ! 1/6 otherwise 50 REAL(KIND=dp) :: x1 ! starting x value if uniform incr. 51 REAL(KIND=dp) :: xn ! end x value if uniform incr. 52 END TYPE spline_data_type 53 54! ************************************************************************************************** 55 TYPE spline_data_p_type 56 TYPE(spline_data_type), POINTER :: spline_data 57 END TYPE spline_data_p_type 58 59! ************************************************************************************************** 60 TYPE spline_data_pp_type 61 TYPE(spline_data_p_type), POINTER, DIMENSION(:) :: spl_p 62 END TYPE spline_data_pp_type 63 64! ************************************************************************************************** 65 TYPE spline_environment_type 66 INTEGER :: ref_count, id_nr 67 TYPE(spline_data_pp_type), POINTER, DIMENSION(:) :: spl_pp 68 INTEGER, POINTER, DIMENSION(:, :) :: spltab 69 END TYPE spline_environment_type 70 71! ************************************************************************************************** 72 TYPE spline_factor_type 73 REAL(KIND=dp) :: rcutsq_f, cutoff 74 REAL(KIND=dp), DIMENSION(:), POINTER :: rscale 75 REAL(KIND=dp), DIMENSION(:), POINTER :: fscale 76 REAL(KIND=dp), DIMENSION(:), POINTER :: dscale 77 END TYPE spline_factor_type 78 79CONTAINS 80 81! ************************************************************************************************** 82!> \brief releases spline_env 83!> \param spline_env ... 84!> \author unknown 85! ************************************************************************************************** 86 SUBROUTINE spline_env_release(spline_env) 87 TYPE(spline_environment_type), POINTER :: spline_env 88 89 CHARACTER(len=*), PARAMETER :: routineN = 'spline_env_release', & 90 routineP = moduleN//':'//routineN 91 92 INTEGER :: i 93 TYPE(spline_data_p_type), DIMENSION(:), POINTER :: spl_p 94 95 IF (ASSOCIATED(spline_env)) THEN 96 CPASSERT(spline_env%ref_count > 0) 97 spline_env%ref_count = spline_env%ref_count - 1 98 IF (spline_env%ref_count < 1) THEN 99 DEALLOCATE (spline_env%spltab) 100 DO i = 1, SIZE(spline_env%spl_pp) 101 spl_p => spline_env%spl_pp(i)%spl_p 102 CALL spline_data_p_release(spl_p) 103 END DO 104 DEALLOCATE (spline_env%spl_pp) 105 DEALLOCATE (spline_env) 106 END IF 107 END IF 108 109 END SUBROUTINE spline_env_release 110 111! ************************************************************************************************** 112!> \brief releases spline_data 113!> \param spline_data ... 114!> \author CJM 115! ************************************************************************************************** 116 SUBROUTINE spline_data_release(spline_data) 117 TYPE(spline_data_type), POINTER :: spline_data 118 119 CHARACTER(len=*), PARAMETER :: routineN = 'spline_data_release', & 120 routineP = moduleN//':'//routineN 121 122 IF (ASSOCIATED(spline_data)) THEN 123 CPASSERT(spline_data%ref_count > 0) 124 spline_data%ref_count = spline_data%ref_count - 1 125 IF (spline_data%ref_count < 1) THEN 126 IF (ASSOCIATED(spline_data%y)) THEN 127 DEALLOCATE (spline_data%y) 128 END IF 129 IF (ASSOCIATED(spline_data%y2)) THEN 130 DEALLOCATE (spline_data%y2) 131 END IF 132 DEALLOCATE (spline_data) 133 END IF 134 END IF 135 END SUBROUTINE spline_data_release 136 137! ************************************************************************************************** 138!> \brief releases spline_data_p 139!> \param spl_p ... 140!> \author CJM 141! ************************************************************************************************** 142 SUBROUTINE spline_data_p_release(spl_p) 143 TYPE(spline_data_p_type), DIMENSION(:), POINTER :: spl_p 144 145 CHARACTER(len=*), PARAMETER :: routineN = 'spline_data_p_release', & 146 routineP = moduleN//':'//routineN 147 148 INTEGER :: i 149 LOGICAL :: release_kind 150 151 IF (ASSOCIATED(spl_p)) THEN 152 release_kind = .TRUE. 153 DO i = 1, SIZE(spl_p) 154 CALL spline_data_release(spl_p(i)%spline_data) 155 release_kind = release_kind .AND. (.NOT. ASSOCIATED(spl_p(i)%spline_data)) 156 END DO 157 IF (release_kind) THEN 158 DEALLOCATE (spl_p) 159 END IF 160 END IF 161 162 END SUBROUTINE spline_data_p_release 163 164! ************************************************************************************************** 165!> \brief retains spline_env 166!> \param spline_data ... 167!> \author CJM 168! ************************************************************************************************** 169 SUBROUTINE spline_data_retain(spline_data) 170 TYPE(spline_data_type), POINTER :: spline_data 171 172 CHARACTER(len=*), PARAMETER :: routineN = 'spline_data_retain', & 173 routineP = moduleN//':'//routineN 174 175 CPASSERT(ASSOCIATED(spline_data)) 176 CPASSERT(spline_data%ref_count > 0) 177 spline_data%ref_count = spline_data%ref_count + 1 178 END SUBROUTINE spline_data_retain 179 180! ************************************************************************************************** 181!> \brief retains spline_data_p_type 182!> \param spl_p ... 183!> \author CJM 184! ************************************************************************************************** 185 SUBROUTINE spline_data_p_retain(spl_p) 186 TYPE(spline_data_p_type), DIMENSION(:), POINTER :: spl_p 187 188 CHARACTER(len=*), PARAMETER :: routineN = 'spline_data_p_retain', & 189 routineP = moduleN//':'//routineN 190 191 INTEGER :: i 192 193 CPASSERT(ASSOCIATED(spl_p)) 194 DO i = 1, SIZE(spl_p) 195 CALL spline_data_retain(spl_p(i)%spline_data) 196 END DO 197 END SUBROUTINE spline_data_p_retain 198 199! ************************************************************************************************** 200!> \brief retains spline_env 201!> \param spline_env ... 202!> \par History 203!> 2001-09-19-HAF added this doc entry and changed formatting 204!> \author unknown 205! ************************************************************************************************** 206 SUBROUTINE spline_env_retain(spline_env) 207 TYPE(spline_environment_type), POINTER :: spline_env 208 209 CHARACTER(len=*), PARAMETER :: routineN = 'spline_env_retain', & 210 routineP = moduleN//':'//routineN 211 212 CPASSERT(ASSOCIATED(spline_env)) 213 CPASSERT(spline_env%ref_count > 0) 214 spline_env%ref_count = spline_env%ref_count + 1 215 END SUBROUTINE spline_env_retain 216 217! ************************************************************************************************** 218!> \brief Data-structure that holds all needed information about 219!> a specific spline interpolation. 220!> \param spline_env ... 221!> \param ntype ... 222!> \param ntab_in ... 223!> \par History 224!> 2001-09-19-HAF added this doc entry and changed formatting 225!> \author unknown 226! ************************************************************************************************** 227 SUBROUTINE spline_env_create(spline_env, ntype, ntab_in) 228 TYPE(spline_environment_type), POINTER :: spline_env 229 INTEGER, INTENT(IN) :: ntype 230 INTEGER, INTENT(IN), OPTIONAL :: ntab_in 231 232 CHARACTER(len=*), PARAMETER :: routineN = 'spline_env_create', & 233 routineP = moduleN//':'//routineN 234 235 INTEGER :: handle, i, isize, j, ntab 236 237 CALL timeset(routineN, handle) 238 239 ALLOCATE (spline_env) 240 NULLIFY (spline_env%spl_pp) 241 NULLIFY (spline_env%spltab) 242 spline_env%ref_count = 1 243 last_spline_env_id_nr = last_spline_env_id_nr + 1 244 spline_env%id_nr = last_spline_env_id_nr 245 ! Allocate the number of spline data tables (upper triangular) 246 IF (PRESENT(ntab_in)) THEN 247 ntab = ntab_in 248 ELSE 249 ntab = (ntype*ntype + ntype)/2 250 END IF 251 ALLOCATE (spline_env%spl_pp(ntab)) 252 253 ALLOCATE (spline_env%spltab(ntype, ntype)) 254 255 DO i = 1, ntab 256 NULLIFY (spline_env%spl_pp(i)%spl_p) 257 isize = 1 258 ALLOCATE (spline_env%spl_pp(i)%spl_p(isize)) 259 DO j = 1, SIZE(spline_env%spl_pp(i)%spl_p) 260 CALL spline_data_create(spline_env%spl_pp(i)%spl_p(j)%spline_data) 261 END DO 262 END DO 263 264 CALL timestop(handle) 265 266 END SUBROUTINE spline_env_create 267 268! ************************************************************************************************** 269!> \brief Copy Data-structure of spline_data_p_type 270!> \param spl_p_source ... 271!> \param spl_p_dest ... 272!> \author teo 06.2007 273! ************************************************************************************************** 274 SUBROUTINE spline_data_p_copy(spl_p_source, spl_p_dest) 275 TYPE(spline_data_p_type), DIMENSION(:), POINTER :: spl_p_source, spl_p_dest 276 277 CHARACTER(len=*), PARAMETER :: routineN = 'spline_data_p_copy', & 278 routineP = moduleN//':'//routineN 279 280 INTEGER :: i, nsized, nsizes 281 282 CPASSERT(ASSOCIATED(spl_p_source)) 283 nsizes = SIZE(spl_p_source) 284 IF (.NOT. ASSOCIATED(spl_p_dest)) THEN 285 ALLOCATE (spl_p_dest(nsizes)) 286 DO i = 1, nsizes 287 NULLIFY (spl_p_dest(i)%spline_data) 288 END DO 289 ELSE 290 nsized = SIZE(spl_p_dest) 291 CPASSERT(nsizes == nsized) 292 DO i = 1, nsizes 293 CALL spline_data_release(spl_p_dest(i)%spline_data) 294 END DO 295 END IF 296 DO i = 1, nsizes 297 CALL spline_data_copy(spl_p_source(i)%spline_data, spl_p_dest(i)%spline_data) 298 END DO 299 END SUBROUTINE spline_data_p_copy 300 301! ************************************************************************************************** 302!> \brief Copy Data-structure that constains spline table 303!> \param spline_data_source ... 304!> \param spline_data_dest ... 305!> \author teo 11.2005 306! ************************************************************************************************** 307 SUBROUTINE spline_data_copy(spline_data_source, spline_data_dest) 308 TYPE(spline_data_type), POINTER :: spline_data_source, spline_data_dest 309 310 CHARACTER(len=*), PARAMETER :: routineN = 'spline_data_copy', & 311 routineP = moduleN//':'//routineN 312 313 CPASSERT(ASSOCIATED(spline_data_source)) 314 IF (.NOT. ASSOCIATED(spline_data_dest)) CALL spline_data_create(spline_data_dest) 315 316 spline_data_dest%ref_count = spline_data_source%ref_count 317 spline_data_dest%id_nr = spline_data_source%id_nr 318 spline_data_dest%n = spline_data_source%n 319 spline_data_dest%h = spline_data_source%h 320 spline_data_dest%invh = spline_data_source%invh 321 spline_data_dest%h26 = spline_data_source%h26 322 spline_data_dest%x1 = spline_data_source%x1 323 spline_data_dest%xn = spline_data_source%xn 324 IF (ASSOCIATED(spline_data_source%y)) THEN 325 ALLOCATE (spline_data_dest%y(SIZE(spline_data_source%y))) 326 spline_data_dest%y = spline_data_source%y 327 END IF 328 IF (ASSOCIATED(spline_data_source%y2)) THEN 329 ALLOCATE (spline_data_dest%y2(SIZE(spline_data_source%y2))) 330 spline_data_dest%y2 = spline_data_source%y2 331 END IF 332 END SUBROUTINE spline_data_copy 333 334! ************************************************************************************************** 335!> \brief Data-structure that constains spline table 336!> \param spline_data ... 337!> \author unknown 338! ************************************************************************************************** 339 SUBROUTINE spline_data_create(spline_data) 340 TYPE(spline_data_type), POINTER :: spline_data 341 342 CHARACTER(len=*), PARAMETER :: routineN = 'spline_data_create', & 343 routineP = moduleN//':'//routineN 344 345 ALLOCATE (spline_data) 346 spline_data%ref_count = 1 347 last_spline_data_id_nr = last_spline_data_id_nr + 1 348 spline_data%id_nr = last_spline_data_id_nr 349 NULLIFY (spline_data%y) 350 NULLIFY (spline_data%y2) 351 END SUBROUTINE spline_data_create 352 353! ************************************************************************************************** 354!> \brief releases spline_factor 355!> \param spline_factor ... 356!> \author teo 357! ************************************************************************************************** 358 SUBROUTINE spline_factor_release(spline_factor) 359 TYPE(spline_factor_type), POINTER :: spline_factor 360 361 CHARACTER(len=*), PARAMETER :: routineN = 'spline_factor_release', & 362 routineP = moduleN//':'//routineN 363 364 IF (ASSOCIATED(spline_factor)) THEN 365 IF (ASSOCIATED(spline_factor%rscale)) THEN 366 DEALLOCATE (spline_factor%rscale) 367 END IF 368 IF (ASSOCIATED(spline_factor%fscale)) THEN 369 DEALLOCATE (spline_factor%fscale) 370 END IF 371 IF (ASSOCIATED(spline_factor%dscale)) THEN 372 DEALLOCATE (spline_factor%dscale) 373 END IF 374 DEALLOCATE (spline_factor) 375 END IF 376 END SUBROUTINE spline_factor_release 377 378! ************************************************************************************************** 379!> \brief releases spline_factor 380!> \param spline_factor ... 381!> \author teo 382! ************************************************************************************************** 383 SUBROUTINE spline_factor_create(spline_factor) 384 TYPE(spline_factor_type), POINTER :: spline_factor 385 386 CHARACTER(len=*), PARAMETER :: routineN = 'spline_factor_create', & 387 routineP = moduleN//':'//routineN 388 389 CPASSERT(.NOT. ASSOCIATED(spline_factor)) 390 ALLOCATE (spline_factor) 391 ALLOCATE (spline_factor%rscale(1)) 392 ALLOCATE (spline_factor%fscale(1)) 393 ALLOCATE (spline_factor%dscale(1)) 394 spline_factor%rscale = 1.0_dp 395 spline_factor%fscale = 1.0_dp 396 spline_factor%dscale = 1.0_dp 397 spline_factor%rcutsq_f = 1.0_dp 398 spline_factor%cutoff = 0.0_dp 399 END SUBROUTINE spline_factor_create 400 401! ************************************************************************************************** 402!> \brief releases spline_factor 403!> \param spline_factor_source ... 404!> \param spline_factor_dest ... 405!> \author teo 406! ************************************************************************************************** 407 SUBROUTINE spline_factor_copy(spline_factor_source, spline_factor_dest) 408 TYPE(spline_factor_type), POINTER :: spline_factor_source, spline_factor_dest 409 410 CHARACTER(len=*), PARAMETER :: routineN = 'spline_factor_copy', & 411 routineP = moduleN//':'//routineN 412 413 INTEGER :: isize, jsize, ksize 414 415 IF (ASSOCIATED(spline_factor_dest)) CALL spline_factor_release(spline_factor_dest) 416 IF (ASSOCIATED(spline_factor_source)) THEN 417 isize = SIZE(spline_factor_source%rscale) 418 jsize = SIZE(spline_factor_source%fscale) 419 ksize = SIZE(spline_factor_source%dscale) 420 CPASSERT(isize == jsize) 421 CPASSERT(isize == ksize) 422 CALL spline_factor_create(spline_factor_dest) 423 spline_factor_dest%rscale = spline_factor_source%rscale 424 spline_factor_dest%fscale = spline_factor_source%fscale 425 spline_factor_dest%dscale = spline_factor_source%dscale 426 spline_factor_dest%rcutsq_f = spline_factor_source%rcutsq_f 427 spline_factor_dest%cutoff = spline_factor_source%cutoff 428 END IF 429 END SUBROUTINE spline_factor_copy 430 431END MODULE splines_types 432