1!--------------------------------------------------------------------------------------------------! 2! CP2K: A general program to perform molecular dynamics simulations ! 3! Copyright (C) 2000 - 2019 CP2K developers group ! 4!--------------------------------------------------------------------------------------------------! 5MODULE optimize_basis_utils 6 USE cp_files, ONLY: close_file,& 7 open_file 8 USE cp_log_handling, ONLY: cp_get_default_logger,& 9 cp_logger_get_default_unit_nr,& 10 cp_logger_type,& 11 cp_to_string 12 USE cp_para_types, ONLY: cp_para_env_type 13 USE cp_parser_methods, ONLY: parser_get_object,& 14 parser_search_string 15 USE cp_parser_types, ONLY: cp_parser_type,& 16 parser_create,& 17 parser_release 18 USE input_constants, ONLY: do_opt_all,& 19 do_opt_coeff,& 20 do_opt_exps,& 21 do_opt_none 22 USE input_section_types, ONLY: section_vals_get,& 23 section_vals_get_subs_vals,& 24 section_vals_type,& 25 section_vals_val_get 26 USE kinds, ONLY: default_path_length,& 27 default_string_length,& 28 dp 29 USE machine, ONLY: default_output_unit,& 30 m_getcwd 31 USE optimize_basis_types, ONLY: basis_optimization_type,& 32 derived_basis_info,& 33 flex_basis_type,& 34 subset_type 35 USE powell, ONLY: opt_state_type 36 USE string_utilities, ONLY: uppercase 37#include "./base/base_uses.f90" 38 39 IMPLICIT NONE 40 PRIVATE 41 42 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'optimize_basis_utils' 43 44 PUBLIC :: optimize_basis_init_read_input, get_set_and_basis_id, & 45 update_derived_basis_sets 46 47CONTAINS 48 49! ************************************************************************************************** 50!> \brief initialize all parts of the optimization type and read input settings 51!> \param opt_bas ... 52!> \param root_section ... 53!> \param para_env ... 54!> \author Florian Schiffmann 55! ************************************************************************************************** 56 57 SUBROUTINE optimize_basis_init_read_input(opt_bas, root_section, para_env) 58 TYPE(basis_optimization_type) :: opt_bas 59 TYPE(section_vals_type), POINTER :: root_section 60 TYPE(cp_para_env_type), POINTER :: para_env 61 62 CHARACTER(len=*), PARAMETER :: routineN = 'optimize_basis_init_read_input', & 63 routineP = moduleN//':'//routineN 64 65 CHARACTER(LEN=default_path_length) :: main_dir 66 INTEGER :: iset, iweight, nrep 67 TYPE(section_vals_type), POINTER :: kind_section, optbas_section, & 68 powell_section, train_section 69 70 optbas_section => section_vals_get_subs_vals(root_section, "OPTIMIZE_BASIS") 71 powell_section => section_vals_get_subs_vals(optbas_section, "OPTIMIZATION") 72 train_section => section_vals_get_subs_vals(optbas_section, "TRAINING_FILES") 73 kind_section => section_vals_get_subs_vals(optbas_section, "FIT_KIND") 74 75 CALL section_vals_val_get(optbas_section, "BASIS_TEMPLATE_FILE", c_val=opt_bas%template_basis_file) 76 CALL section_vals_val_get(optbas_section, "BASIS_WORK_FILE", c_val=opt_bas%work_basis_file) 77 CALL section_vals_val_get(optbas_section, "BASIS_OUTPUT_FILE", c_val=opt_bas%output_basis_file) 78 CALL m_getcwd(main_dir) 79 opt_bas%work_basis_file = TRIM(ADJUSTL(main_dir))//"/"//TRIM(ADJUSTL(opt_bas%work_basis_file)) 80 81 CALL section_vals_val_get(optbas_section, "WRITE_FREQUENCY", i_val=opt_bas%write_frequency) 82 CALL section_vals_val_get(optbas_section, "USE_CONDITION_NUMBER", l_val=opt_bas%use_condition_number) 83 84 CALL generate_initial_basis(kind_section, opt_bas, para_env) 85 86 CALL section_vals_get(train_section, n_repetition=opt_bas%ntraining_sets) 87 IF (opt_bas%ntraining_sets == 0) & 88 CPABORT("No training set was specified in the Input") 89 90 ALLOCATE (opt_bas%training_input(opt_bas%ntraining_sets)) 91 ALLOCATE (opt_bas%training_dir(opt_bas%ntraining_sets)) 92 DO iset = 1, opt_bas%ntraining_sets 93 CALL section_vals_val_get(train_section, "DIRECTORY", c_val=opt_bas%training_dir(iset), & 94 i_rep_section=iset) 95 CALL section_vals_val_get(train_section, "INPUT_FILE_NAME", c_val=opt_bas%training_input(iset), & 96 i_rep_section=iset) 97 END DO 98 99 CALL init_powell_var(opt_bas%powell_param, powell_section) 100 opt_bas%powell_param%nvar = SIZE(opt_bas%x_opt) 101 102 CALL generate_derived_basis_sets(opt_bas, para_env) 103 104 CALL generate_basis_combinations(opt_bas, optbas_section) 105 106 CALL section_vals_val_get(optbas_section, "RESIDUUM_WEIGHT", n_rep_val=nrep) 107 ALLOCATE (opt_bas%fval_weight(0:opt_bas%ncombinations)) 108 opt_bas%fval_weight = 1.0_dp 109 DO iweight = 1, nrep 110 CALL section_vals_val_get(optbas_section, "RESIDUUM_WEIGHT", r_val=opt_bas%fval_weight(iweight - 1), & 111 i_rep_val=iweight) 112 END DO 113 114 CALL section_vals_val_get(optbas_section, "CONDITION_WEIGHT", n_rep_val=nrep) 115 ALLOCATE (opt_bas%condition_weight(0:opt_bas%ncombinations)) 116 opt_bas%condition_weight = 1.0_dp 117 DO iweight = 1, nrep 118 CALL section_vals_val_get(optbas_section, "CONDITION_WEIGHT", r_val=opt_bas%condition_weight(iweight - 1), & 119 i_rep_val=iweight) 120 END DO 121 122 CALL generate_computation_groups(opt_bas, optbas_section, para_env) 123 124 CALL print_opt_info(opt_bas) 125 126 END SUBROUTINE optimize_basis_init_read_input 127 128! ************************************************************************************************** 129!> \brief ... 130!> \param opt_bas ... 131! ************************************************************************************************** 132 SUBROUTINE print_opt_info(opt_bas) 133 TYPE(basis_optimization_type) :: opt_bas 134 135 CHARACTER(len=*), PARAMETER :: routineN = 'print_opt_info', routineP = moduleN//':'//routineN 136 137 INTEGER :: icomb, ikind, unit_nr 138 TYPE(cp_logger_type), POINTER :: logger 139 140 logger => cp_get_default_logger() 141 unit_nr = -1 142 IF (logger%para_env%ionode) & 143 unit_nr = cp_logger_get_default_unit_nr(logger) 144 145 IF (unit_nr > 0) THEN 146 WRITE (unit_nr, '(1X,A,A)') "BASOPT| Total number of calculations ", & 147 TRIM(cp_to_string(opt_bas%ncombinations*opt_bas%ntraining_sets)) 148 WRITE (unit_nr, '(A)') "" 149 DO icomb = 1, opt_bas%ncombinations 150 WRITE (unit_nr, '(1X,A,A)') "BASOPT| Content of basis combination ", TRIM(cp_to_string(icomb)) 151 DO ikind = 1, opt_bas%nkind 152 WRITE (unit_nr, '(1X,A,A4,4X,A,3X,A)') "BASOPT| Element: ", TRIM(opt_bas%kind_basis(ikind)%element), & 153 "Basis set: ", TRIM(opt_bas%kind_basis(ikind)%flex_basis(opt_bas%combination(icomb, ikind))%basis_name) 154 END DO 155 WRITE (unit_nr, '(A)') "" 156 END DO 157 END IF 158 END SUBROUTINE print_opt_info 159 160! ************************************************************************************************** 161!> \brief Generation of the requested basis set combinations if multiple kinds 162!> are fitted at the same time (if not specified create all possible) 163!> \param opt_bas ... 164!> \param optbas_section ... 165!> \author Florian Schiffmann 166! ************************************************************************************************** 167 SUBROUTINE generate_basis_combinations(opt_bas, optbas_section) 168 TYPE(basis_optimization_type) :: opt_bas 169 TYPE(section_vals_type), POINTER :: optbas_section 170 171 CHARACTER(len=*), PARAMETER :: routineN = 'generate_basis_combinations', & 172 routineP = moduleN//':'//routineN 173 174 INTEGER :: i, ikind, j, n_rep 175 INTEGER, DIMENSION(:), POINTER :: i_vals, tmp_i, tmp_i2 176 LOGICAL :: explicit, raise 177 178!setup the basis combinations to optimize 179 180 CALL section_vals_val_get(optbas_section, "BASIS_COMBINATIONS", explicit=explicit, n_rep_val=n_rep) 181 IF (.NOT. explicit) THEN 182 opt_bas%ncombinations = 1 183 ALLOCATE (tmp_i(opt_bas%nkind)) 184 ALLOCATE (tmp_i2(opt_bas%nkind)) 185 DO ikind = 1, opt_bas%nkind 186 opt_bas%ncombinations = opt_bas%ncombinations*SIZE(opt_bas%kind_basis(ikind)%flex_basis) 187 tmp_i(ikind) = opt_bas%kind_basis(ikind)%nbasis_deriv 188 END DO 189 ALLOCATE (opt_bas%combination(opt_bas%ncombinations, opt_bas%nkind)) 190 tmp_i2 = 0 191 DO i = 1, opt_bas%ncombinations 192 DO j = 1, opt_bas%nkind 193 opt_bas%combination(i, j) = tmp_i2(j) 194 END DO 195 tmp_i2(opt_bas%nkind) = tmp_i2(opt_bas%nkind) + 1 196 raise = .FALSE. 197 DO j = opt_bas%nkind, 1, -1 198 IF (raise) tmp_i2(j) = tmp_i2(j) + 1 199 IF (tmp_i2(j) .GT. tmp_i(j)) THEN 200 tmp_i2(j) = 0 201 raise = .TRUE. 202 END IF 203 END DO 204 END DO 205 DEALLOCATE (tmp_i) 206 DEALLOCATE (tmp_i2) 207 ELSE 208 opt_bas%ncombinations = n_rep 209 ALLOCATE (opt_bas%combination(opt_bas%ncombinations, opt_bas%nkind)) 210 DO i = 1, n_rep 211 CALL section_vals_val_get(optbas_section, "BASIS_COMBINATIONS", i_vals=i_vals, i_rep_val=i) 212 opt_bas%combination(i, :) = i_vals(:) 213 END DO 214 END IF 215 216 END SUBROUTINE generate_basis_combinations 217 218! ************************************************************************************************** 219!> \brief returns a mapping from the calculation id to the trainings set id and 220!> basis combination id 221!> \param calc_id ... 222!> \param opt_bas ... 223!> \param set_id ... 224!> \param bas_id ... 225!> \author Florian Schiffmann 226! ************************************************************************************************** 227 228 SUBROUTINE get_set_and_basis_id(calc_id, opt_bas, set_id, bas_id) 229 230 INTEGER :: calc_id 231 TYPE(basis_optimization_type) :: opt_bas 232 INTEGER :: set_id, bas_id 233 234 CHARACTER(len=*), PARAMETER :: routineN = 'get_set_and_basis_id', & 235 routineP = moduleN//':'//routineN 236 237 INTEGER :: ncom, nset 238 239 ncom = opt_bas%ncombinations 240 nset = opt_bas%ntraining_sets 241 242 set_id = (calc_id)/ncom + 1 243 bas_id = MOD(calc_id, ncom) + 1 244 245 END SUBROUTINE 246 247! ************************************************************************************************** 248!> \brief Pack calculations in groups. If less MPI tasks than systems are used 249!> multiple systems will be assigned to a single MPI task 250!> \param opt_bas ... 251!> \param optbas_section ... 252!> \param para_env ... 253!> \author Florian Schiffmann 254! ************************************************************************************************** 255 256 SUBROUTINE generate_computation_groups(opt_bas, optbas_section, para_env) 257 TYPE(basis_optimization_type) :: opt_bas 258 TYPE(section_vals_type), POINTER :: optbas_section 259 TYPE(cp_para_env_type), POINTER :: para_env 260 261 CHARACTER(len=*), PARAMETER :: routineN = 'generate_computation_groups', & 262 routineP = moduleN//':'//routineN 263 264 INTEGER :: iadd1, iadd2, icount, igroup, isize, j, & 265 ncalc, nproc, nptot 266 INTEGER, DIMENSION(:), POINTER :: i_vals 267 LOGICAL :: explicit 268 269 nproc = para_env%num_pe 270 ncalc = opt_bas%ncombinations*opt_bas%ntraining_sets 271 CALL section_vals_val_get(optbas_section, "GROUP_PARTITION", explicit=explicit) 272 273 ! No input information available, try to equally distribute 274 IF (.NOT. explicit) THEN 275 IF (nproc .GE. ncalc) THEN 276 iadd1 = nproc/ncalc 277 iadd2 = MOD(nproc, ncalc) 278 ALLOCATE (opt_bas%comp_group(ncalc)) 279 ALLOCATE (opt_bas%group_partition(0:ncalc - 1)) 280 DO igroup = 0, ncalc - 1 281 ALLOCATE (opt_bas%comp_group(igroup + 1)%member_list(1)) 282 opt_bas%comp_group(igroup + 1)%member_list(1) = igroup 283 opt_bas%group_partition(igroup) = iadd1 284 IF (igroup .LT. iadd2) opt_bas%group_partition(igroup) = opt_bas%group_partition(igroup) + 1 285 END DO 286 ELSE 287 iadd1 = ncalc/nproc 288 iadd2 = MOD(ncalc, nproc) 289 ALLOCATE (opt_bas%comp_group(nproc)) 290 ALLOCATE (opt_bas%group_partition(0:nproc - 1)) 291 icount = 0 292 DO igroup = 0, nproc - 1 293 opt_bas%group_partition(igroup) = 1 294 isize = iadd1 295 IF (igroup .LT. iadd2) isize = isize + 1 296 ALLOCATE (opt_bas%comp_group(igroup + 1)%member_list(isize)) 297 DO j = 1, isize 298 opt_bas%comp_group(igroup + 1)%member_list(j) = icount 299 icount = icount + 1 300 END DO 301 END DO 302 END IF 303 ELSE 304 305 ! Group partition from input. see if all systems can be assigned. If not add to existing group 306 CALL section_vals_val_get(optbas_section, "GROUP_PARTITION", i_vals=i_vals) 307 isize = SIZE(i_vals) 308 nptot = SUM(i_vals) 309 IF (nptot /= nproc) & 310 CALL cp_abort(__LOCATION__, & 311 "Number of processors in group distribution does not match number of MPI tasks."// & 312 " Please change input.") 313 IF (.NOT. isize .LE. ncalc) & 314 CALL cp_abort(__LOCATION__, & 315 "Number of Groups larger than number of calculations"// & 316 " Please change input.") 317 CPASSERT(nptot == nproc) 318 ALLOCATE (opt_bas%comp_group(isize)) 319 ALLOCATE (opt_bas%group_partition(0:isize - 1)) 320 IF (isize .LT. ncalc) THEN 321 iadd1 = ncalc/isize 322 iadd2 = MOD(ncalc, isize) 323 icount = 0 324 DO igroup = 0, isize - 1 325 opt_bas%group_partition(igroup) = i_vals(igroup + 1) 326 isize = iadd1 327 IF (igroup .LT. iadd2) isize = isize + 1 328 ALLOCATE (opt_bas%comp_group(igroup + 1)%member_list(isize)) 329 DO j = 1, isize 330 opt_bas%comp_group(igroup + 1)%member_list(j) = icount 331 icount = icount + 1 332 END DO 333 END DO 334 ELSE 335 DO igroup = 0, isize - 1 336 opt_bas%group_partition(igroup) = i_vals(igroup + 1) 337 ALLOCATE (opt_bas%comp_group(igroup + 1)%member_list(1)) 338 opt_bas%comp_group(igroup + 1)%member_list(1) = igroup 339 END DO 340 END IF 341 END IF 342 END SUBROUTINE generate_computation_groups 343 344! ************************************************************************************************** 345!> \brief Regenerate the basis sets from reference 0 after an update from the 346!> optimizer to reference was performed, and print basis file if required 347!> \param opt_bas ... 348!> \param write_it ... 349!> \param output_file ... 350!> \param para_env ... 351!> \author Florian Schiffmann 352! ************************************************************************************************** 353 SUBROUTINE update_derived_basis_sets(opt_bas, write_it, output_file, para_env) 354 TYPE(basis_optimization_type) :: opt_bas 355 LOGICAL :: write_it 356 CHARACTER(LEN=default_path_length) :: output_file 357 TYPE(cp_para_env_type), POINTER :: para_env 358 359 CHARACTER(len=*), PARAMETER :: routineN = 'update_derived_basis_sets', & 360 routineP = moduleN//':'//routineN 361 362 INTEGER :: ibasis, ikind, unit_nr 363 364 DO ikind = 1, opt_bas%nkind 365 DO ibasis = 1, opt_bas%kind_basis(ikind)%nbasis_deriv 366 CALL update_used_parts(opt_bas%kind_basis(ikind)%deriv_info(ibasis), & 367 opt_bas%kind_basis(ikind)%flex_basis(0), & 368 opt_bas%kind_basis(ikind)%flex_basis(ibasis)) 369 END DO 370 END DO 371 372 IF (write_it) THEN 373 IF (para_env%ionode) THEN 374 CALL open_file(file_name=output_file, file_status="UNKNOWN", & 375 file_action="WRITE", unit_number=unit_nr) 376 ELSE 377 unit_nr = -999 378 END IF 379 DO ikind = 1, opt_bas%nkind 380 DO ibasis = 0, opt_bas%kind_basis(ikind)%nbasis_deriv 381 CALL write_basis(opt_bas%kind_basis(ikind)%flex_basis(ibasis), opt_bas%kind_basis(ikind)%element, & 382 unit_nr) 383 END DO 384 END DO 385 IF (para_env%ionode) CALL close_file(unit_number=unit_nr) 386 END IF 387 388 END SUBROUTINE update_derived_basis_sets 389 390! ************************************************************************************************** 391!> \brief Update the the information in a given basis set 392!> \param info_new ... 393!> \param basis ... 394!> \param basis_new ... 395!> \author Florian Schiffmann 396! ************************************************************************************************** 397 398 SUBROUTINE update_used_parts(info_new, basis, basis_new) 399 TYPE(derived_basis_info) :: info_new 400 TYPE(flex_basis_type) :: basis, basis_new 401 402 CHARACTER(len=*), PARAMETER :: routineN = 'update_used_parts', & 403 routineP = moduleN//':'//routineN 404 405 INTEGER :: icont, iset, jcont, jset 406 407 jset = 0 408 DO iset = 1, basis%nsets 409 IF (info_new%in_use_set(iset)) THEN 410 jset = jset + 1 411 basis_new%subset(jset)%exps(:) = basis%subset(iset)%exps 412 jcont = 0 413 DO icont = 1, basis%subset(iset)%ncon_tot 414 IF (info_new%use_contr(iset)%in_use(icont)) THEN 415 jcont = jcont + 1 416 basis_new%subset(jset)%coeff(:, jcont) = basis%subset(iset)%coeff(:, icont) 417 END IF 418 END DO 419 END IF 420 END DO 421 422 END SUBROUTINE update_used_parts 423 424! ************************************************************************************************** 425!> \brief Initial generation of the basis set from the file and DERIVED_SET 426!> \param opt_bas ... 427!> \param para_env ... 428!> \author Florian Schiffmann 429! ************************************************************************************************** 430 431 SUBROUTINE generate_derived_basis_sets(opt_bas, para_env) 432 TYPE(basis_optimization_type) :: opt_bas 433 TYPE(cp_para_env_type), POINTER :: para_env 434 435 CHARACTER(len=*), PARAMETER :: routineN = 'generate_derived_basis_sets', & 436 routineP = moduleN//':'//routineN 437 438 INTEGER :: ibasis, ikind, iref, jbasis, unit_nr 439 440 DO ikind = 1, opt_bas%nkind 441 CALL init_deriv_info_ref(opt_bas%kind_basis(ikind)%deriv_info(0), opt_bas%kind_basis(ikind)%flex_basis(0)) 442 opt_bas%kind_basis(ikind)%deriv_info(0)%basis_name = TRIM(ADJUSTL(opt_bas%kind_basis(ikind)%basis_name)) 443 ! initialize the reference set used as template for the rest 444 DO ibasis = 1, opt_bas%kind_basis(ikind)%nbasis_deriv 445 iref = opt_bas%kind_basis(ikind)%deriv_info(ibasis)%reference_set 446 DO jbasis = 0, opt_bas%kind_basis(ikind)%nbasis_deriv 447 IF (iref == jbasis) CALL setup_used_parts_init_basis(opt_bas%kind_basis(ikind)%deriv_info(ibasis), & 448 opt_bas%kind_basis(ikind)%deriv_info(iref), & 449 opt_bas%kind_basis(ikind)%flex_basis(0), & 450 opt_bas%kind_basis(ikind)%flex_basis(ibasis)) 451 END DO 452 END DO 453 END DO 454 455 IF (para_env%ionode) THEN 456 CALL open_file(file_name=opt_bas%work_basis_file, file_status="UNKNOWN", & 457 file_action="WRITE", unit_number=unit_nr) 458 ELSE 459 unit_nr = -999 460 END IF 461 DO ikind = 1, opt_bas%nkind 462 DO ibasis = 0, opt_bas%kind_basis(ikind)%nbasis_deriv 463 IF (LEN_TRIM(opt_bas%kind_basis(ikind)%deriv_info(ibasis)%basis_name) > 0) THEN 464 opt_bas%kind_basis(ikind)%flex_basis(ibasis)%basis_name = & 465 TRIM(ADJUSTL(opt_bas%kind_basis(ikind)%deriv_info(ibasis)%basis_name)) 466 ELSE 467 opt_bas%kind_basis(ikind)%flex_basis(ibasis)%basis_name = & 468 TRIM(ADJUSTL(opt_bas%kind_basis(ikind)%basis_name))//"-DERIVED_SET-"//ADJUSTL(cp_to_string(ibasis)) 469 END IF 470 CALL write_basis(opt_bas%kind_basis(ikind)%flex_basis(ibasis), opt_bas%kind_basis(ikind)%element, & 471 unit_nr) 472 END DO 473 END DO 474 IF (para_env%ionode) CALL close_file(unit_number=unit_nr) 475 476 END SUBROUTINE generate_derived_basis_sets 477 478! ************************************************************************************************** 479!> \brief Write a basis set file which can be used from CP2K 480!> \param basis ... 481!> \param element ... 482!> \param unit_nr ... 483!> \author Florian Schiffmann 484! ************************************************************************************************** 485 486 SUBROUTINE write_basis(basis, element, unit_nr) 487 TYPE(flex_basis_type) :: basis 488 CHARACTER(LEN=default_string_length) :: element 489 INTEGER :: unit_nr 490 491 CHARACTER(len=*), PARAMETER :: routineN = 'write_basis', routineP = moduleN//':'//routineN 492 493 INTEGER :: iexp, iset 494 495 IF (unit_nr > 0) THEN 496 WRITE (UNIT=unit_nr, FMT="(A)") TRIM(ADJUSTL(element))//" "//TRIM(ADJUSTL(basis%basis_name)) 497 WRITE (UNIT=unit_nr, FMT="(1X,I0)") basis%nsets 498 DO iset = 1, basis%nsets 499 WRITE (UNIT=unit_nr, FMT="(30(1X,I0))") basis%subset(iset)%n, basis%subset(iset)%lmin, basis%subset(iset)%lmax, & 500 basis%subset(iset)%nexp, basis%subset(iset)%l 501 DO iexp = 1, basis%subset(iset)%nexp 502 WRITE (UNIT=unit_nr, FMT="(T2,F24.14,30(1X,ES24.14))") & 503 basis%subset(iset)%exps(iexp), basis%subset(iset)%coeff(iexp, :) 504 END DO 505 END DO 506 END IF 507 508 END SUBROUTINE write_basis 509 510! ************************************************************************************************** 511!> \brief Initialize the derived basis sets and the vectors containing their 512!> assembly information ehich is used for regeneration of the sets. 513!> \param info_new ... 514!> \param info_ref ... 515!> \param basis ... 516!> \param basis_new ... 517!> \author Florian Schiffmann 518! ************************************************************************************************** 519 520 SUBROUTINE setup_used_parts_init_basis(info_new, info_ref, basis, basis_new) 521 TYPE(derived_basis_info) :: info_new, info_ref 522 TYPE(flex_basis_type) :: basis, basis_new 523 524 CHARACTER(len=*), PARAMETER :: routineN = 'setup_used_parts_init_basis', & 525 routineP = moduleN//':'//routineN 526 527 INTEGER :: i, jset, lind, nsets 528 529! copy the reference information on the new set 530 531 ALLOCATE (info_new%in_use_set(SIZE(info_ref%in_use_set))) 532 info_new%in_use_set(:) = info_ref%in_use_set 533 ALLOCATE (info_new%use_contr(SIZE(info_ref%in_use_set))) 534 DO i = 1, SIZE(info_ref%in_use_set) 535 ALLOCATE (info_new%use_contr(i)%in_use(SIZE(info_ref%use_contr(i)%in_use))) 536 info_new%use_contr(i)%in_use(:) = info_ref%use_contr(i)%in_use 537 END DO 538 DO i = 1, info_new%nsets 539 info_new%in_use_set(info_new%remove_set(i)) = .FALSE. 540 END DO 541 DO i = 1, info_new%ncontr 542 lind = convert_l_contr_to_entry(basis%subset(info_new%remove_contr(i, 1))%lmin, & 543 basis%subset(info_new%remove_contr(i, 1))%l, & 544 info_new%remove_contr(i, 3), info_new%remove_contr(i, 2)) 545 546 info_new%use_contr(info_new%remove_contr(i, 1))%in_use(lind) = .FALSE. 547 END DO 548 549 nsets = 0 550 DO i = 1, basis%nsets 551 IF (info_new%in_use_set(i)) nsets = nsets + 1 552 END DO 553 basis_new%nsets = nsets 554 ALLOCATE (basis_new%subset(nsets)) 555 jset = 0 556 DO i = 1, basis%nsets 557 IF (info_new%in_use_set(i)) THEN 558 jset = jset + 1 559 CALL create_new_subset(basis%subset(i), basis_new%subset(jset), info_new%use_contr(jset)%in_use) 560 END IF 561 END DO 562 563 END SUBROUTINE setup_used_parts_init_basis 564 565! ************************************************************************************************** 566!> \brief Fill the low level information of the derived basis set. 567!> \param subset ... 568!> \param subset_new ... 569!> \param in_use ... 570!> \author Florian Schiffmann 571! ************************************************************************************************** 572 573 SUBROUTINE create_new_subset(subset, subset_new, in_use) 574 TYPE(subset_type) :: subset, subset_new 575 LOGICAL, DIMENSION(:) :: in_use 576 577 CHARACTER(len=*), PARAMETER :: routineN = 'create_new_subset', & 578 routineP = moduleN//':'//routineN 579 580 INTEGER :: icon, iind, il 581 INTEGER, ALLOCATABLE, DIMENSION(:) :: tmp_l 582 583 ALLOCATE (tmp_l(SIZE(subset%l))) 584 tmp_l(:) = subset%l 585 subset_new%lmin = subset%lmin 586 subset_new%lmax = subset%lmin - 1 587 subset_new%nexp = subset%nexp 588 subset_new%n = subset%n 589 DO il = 1, SIZE(subset%l) 590 DO icon = 1, subset%l(il) 591 iind = convert_l_contr_to_entry(subset%lmin, subset%l, icon, subset%lmin + il - 1) 592 IF (.NOT. in_use(iind)) tmp_l(il) = tmp_l(il) - 1 593 END DO 594 IF (tmp_l(il) .GT. 0) subset_new%lmax = subset_new%lmax + 1 595 END DO 596 subset_new%nl = subset_new%lmax - subset_new%lmin + 1 597 subset_new%ncon_tot = SUM(tmp_l) 598 ALLOCATE (subset_new%l(subset_new%nl)) 599 ALLOCATE (subset_new%coeff(subset_new%nexp, subset_new%ncon_tot)) 600 ALLOCATE (subset_new%exps(subset_new%nexp)) 601 subset_new%exps(:) = subset%exps 602 DO il = 1, SIZE(subset%l) 603 IF (tmp_l(il) == 0) EXIT 604 subset_new%l(il) = tmp_l(il) 605 END DO 606 DEALLOCATE (tmp_l) 607 iind = 0 608 DO icon = 1, subset%ncon_tot 609 IF (in_use(icon)) THEN 610 iind = iind + 1 611 subset_new%coeff(:, iind) = subset%coeff(:, icon) 612 END IF 613 END DO 614 615 END SUBROUTINE create_new_subset 616 617! ************************************************************************************************** 618!> \brief for completeness generate the derived info for set 0(reference from file) 619!> \param info ... 620!> \param basis ... 621!> \author Florian Schiffmann 622! ************************************************************************************************** 623 624 SUBROUTINE init_deriv_info_ref(info, basis) 625 TYPE(derived_basis_info) :: info 626 TYPE(flex_basis_type) :: basis 627 628 CHARACTER(len=*), PARAMETER :: routineN = 'init_deriv_info_ref', & 629 routineP = moduleN//':'//routineN 630 631 INTEGER :: i 632 633 ALLOCATE (info%in_use_set(basis%nsets)) 634 info%in_use_set = .TRUE. 635 ALLOCATE (info%use_contr(basis%nsets)) 636 DO i = 1, basis%nsets 637 ALLOCATE (info%use_contr(i)%in_use(basis%subset(i)%ncon_tot)) 638 info%use_contr(i)%in_use = .TRUE. 639 END DO 640 641 END SUBROUTINE init_deriv_info_ref 642 643! ************************************************************************************************** 644!> \brief get the general information for the basis sets. E.g. what to optimize, 645!> Basis set name, constraints upon optimization and read the reference basis 646!> \param kind_section ... 647!> \param opt_bas ... 648!> \param para_env ... 649!> \author Florian Schiffmann 650! ************************************************************************************************** 651 652 SUBROUTINE generate_initial_basis(kind_section, opt_bas, para_env) 653 TYPE(section_vals_type), POINTER :: kind_section 654 TYPE(basis_optimization_type) :: opt_bas 655 TYPE(cp_para_env_type), POINTER :: para_env 656 657 CHARACTER(len=*), PARAMETER :: routineN = 'generate_initial_basis', & 658 routineP = moduleN//':'//routineN 659 660 INTEGER :: ikind, variable_counter 661 LOGICAL :: explicit 662 TYPE(section_vals_type), POINTER :: set_section 663 664 CALL section_vals_get(kind_section, n_repetition=opt_bas%nkind) 665 ALLOCATE (opt_bas%kind_basis(opt_bas%nkind)) 666 667 ! counter to get the number of free variables in optimization 668 variable_counter = 0 669 DO ikind = 1, opt_bas%nkind 670 CALL section_vals_val_get(kind_section, "_SECTION_PARAMETERS_", c_val=opt_bas%kind_basis(ikind)%element, & 671 i_rep_section=ikind) 672 CALL section_vals_val_get(kind_section, "BASIS_SET", c_val=opt_bas%kind_basis(ikind)%basis_name, & 673 i_rep_section=ikind) 674 set_section => section_vals_get_subs_vals(kind_section, "DERIVED_BASIS_SETS", & 675 i_rep_section=ikind) 676 CALL section_vals_get(set_section, n_repetition=opt_bas%kind_basis(ikind)%nbasis_deriv, explicit=explicit) 677 IF (.NOT. explicit) opt_bas%kind_basis(ikind)%nbasis_deriv = 0 678 ALLOCATE (opt_bas%kind_basis(ikind)%flex_basis(0:opt_bas%kind_basis(ikind)%nbasis_deriv)) 679 ALLOCATE (opt_bas%kind_basis(ikind)%deriv_info(0:opt_bas%kind_basis(ikind)%nbasis_deriv)) 680 681 CALL fill_basis_template(kind_section, opt_bas%kind_basis(ikind)%flex_basis(0), opt_bas%template_basis_file, & 682 opt_bas%kind_basis(ikind)%element, opt_bas%kind_basis(ikind)%basis_name, para_env, ikind) 683 684 CALL setup_exp_constraints(kind_section, opt_bas%kind_basis(ikind)%flex_basis(0)) 685 686 CALL parse_derived_basis(kind_section, opt_bas%kind_basis(ikind)%deriv_info, ikind) 687 688 variable_counter = variable_counter + opt_bas%kind_basis(ikind)%flex_basis(0)%nopt 689 END DO 690 691 ALLOCATE (opt_bas%x_opt(variable_counter)) 692 693 variable_counter = 0 694 DO ikind = 1, opt_bas%nkind 695 CALL assign_x_to_basis(opt_bas%x_opt, opt_bas%kind_basis(ikind)%flex_basis(0), variable_counter) 696 END DO 697 698 CPASSERT(variable_counter == SIZE(opt_bas%x_opt)) 699 700 END SUBROUTINE generate_initial_basis 701 702! ************************************************************************************************** 703!> \brief get low level information about how to construc new basis sets from reference 704!> \param kind_section ... 705!> \param deriv_info ... 706!> \param ikind ... 707!> \author Florian Schiffmann 708! ************************************************************************************************** 709 710 SUBROUTINE parse_derived_basis(kind_section, deriv_info, ikind) 711 TYPE(section_vals_type), POINTER :: kind_section 712 TYPE(derived_basis_info), DIMENSION(:) :: deriv_info 713 INTEGER :: ikind 714 715 CHARACTER(len=*), PARAMETER :: routineN = 'parse_derived_basis', & 716 routineP = moduleN//':'//routineN 717 718 INTEGER :: i_rep, iset, jset, n_rep, nsets 719 INTEGER, DIMENSION(:), POINTER :: i_vals 720 LOGICAL :: explicit 721 TYPE(section_vals_type), POINTER :: set1_section 722 723 nsets = SIZE(deriv_info) - 1 724 set1_section => section_vals_get_subs_vals(kind_section, "DERIVED_BASIS_SETS", & 725 i_rep_section=ikind) 726 DO jset = 1, nsets 727 ! stracnge but as derive info is allcated from 0 to n the count over here has to be shifted 728 iset = jset + 1 729 CALL section_vals_val_get(set1_section, "BASIS_SET_NAME", c_val=deriv_info(iset)%basis_name, & 730 i_rep_section=jset) 731 CALL section_vals_val_get(set1_section, "REFERENCE_SET", i_vals=i_vals, i_rep_section=jset) 732 deriv_info(iset)%reference_set = i_vals(1) 733 CALL section_vals_val_get(set1_section, "REMOVE_CONTRACTION", explicit=explicit, n_rep_val=n_rep, & 734 i_rep_section=jset) 735 deriv_info(iset)%ncontr = n_rep 736 IF (explicit) THEN 737 ALLOCATE (deriv_info(iset)%remove_contr(n_rep, 3)) 738 DO i_rep = 1, n_rep 739 CALL section_vals_val_get(set1_section, "REMOVE_CONTRACTION", i_rep_val=i_rep, i_vals=i_vals, & 740 i_rep_section=jset) 741 deriv_info(iset)%remove_contr(i_rep, :) = i_vals(:) 742 END DO 743 END IF 744 CALL section_vals_val_get(set1_section, "REMOVE_SET", explicit=explicit, n_rep_val=n_rep, & 745 i_rep_section=jset) 746 deriv_info(iset)%nsets = n_rep 747 IF (explicit) THEN 748 ALLOCATE (deriv_info(iset)%remove_set(n_rep)) 749 DO i_rep = 1, n_rep 750 CALL section_vals_val_get(set1_section, "REMOVE_SET", i_rep_val=i_rep, i_vals=i_vals, & 751 i_rep_section=jset) 752 deriv_info(iset)%remove_set(i_rep) = i_vals(1) 753 END DO 754 END IF 755 END DO 756 757 END SUBROUTINE parse_derived_basis 758 759! ************************************************************************************************** 760!> \brief get low level information about constraint on exponents 761!> \param kind1_section ... 762!> \param flex_basis ... 763!> \author Florian Schiffmann 764! ************************************************************************************************** 765 766 SUBROUTINE setup_exp_constraints(kind1_section, flex_basis) 767 TYPE(section_vals_type), POINTER :: kind1_section 768 TYPE(flex_basis_type) :: flex_basis 769 770 CHARACTER(len=*), PARAMETER :: routineN = 'setup_exp_constraints', & 771 routineP = moduleN//':'//routineN 772 773 INTEGER :: ipgf, irep, iset, nrep 774 INTEGER, DIMENSION(:), POINTER :: def_exp 775 LOGICAL :: is_bound, is_varlim 776 TYPE(section_vals_type), POINTER :: const_section 777 778 const_section => section_vals_get_subs_vals(kind1_section, "CONSTRAIN_EXPONENTS") 779 CALL section_vals_get(const_section, n_repetition=nrep) 780 DO irep = 1, nrep 781 CALL section_vals_val_get(const_section, "USE_EXP", i_vals=def_exp, i_rep_section=irep) 782 CALL section_vals_val_get(const_section, "BOUNDARIES", explicit=is_bound, i_rep_section=irep) 783 CALL section_vals_val_get(const_section, "MAX_VAR_FRACTION", explicit=is_varlim, i_rep_section=irep) 784 IF (is_bound .AND. is_varlim) & 785 CALL cp_abort(__LOCATION__, "Exponent has two constraints. "// & 786 "This is not possible at the moment. Please change input.") 787 IF (.NOT. is_bound .AND. .NOT. is_varlim) & 788 CALL cp_abort(__LOCATION__, "Exponent is declared to be constraint but none is given"// & 789 " Please change input.") 790 IF (def_exp(1) == -1) THEN 791 DO iset = 1, flex_basis%nsets 792 IF (def_exp(2) == -1) THEN 793 DO ipgf = 1, flex_basis%subset(iset)%nexp 794 CALL set_constraint(flex_basis, iset, ipgf, const_section, is_bound, is_varlim, irep) 795 END DO 796 ELSE 797 IF (def_exp(2) .LE. flex_basis%subset(iset)%nexp) & 798 CALL cp_abort(__LOCATION__, & 799 "Exponent declared in constraint is larger than number of exponents in the set"// & 800 " Please change input.") 801 CALL set_constraint(flex_basis, iset, def_exp(2), const_section, is_bound, is_varlim, irep) 802 END IF 803 END DO 804 ELSE 805 IF (.NOT. def_exp(1) .LE. flex_basis%nsets) & 806 CALL cp_abort(__LOCATION__, & 807 "Set number of constraint is larger than number of sets in the template basis set."// & 808 " Please change input.") 809 IF (def_exp(2) == -1) THEN 810 DO ipgf = 1, flex_basis%subset(iset)%nexp 811 CALL set_constraint(flex_basis, def_exp(1), ipgf, const_section, is_bound, is_varlim, irep) 812 END DO 813 ELSE 814 IF (.NOT. def_exp(2) .LE. flex_basis%subset(def_exp(1))%nexp) & 815 CALL cp_abort(__LOCATION__, & 816 "Exponent declared in constraint is larger than number of exponents in the set"// & 817 " Please change input.") 818 CALL set_constraint(flex_basis, def_exp(1), def_exp(2), const_section, is_bound, is_varlim, irep) 819 END IF 820 END IF 821 END DO 822 823 END SUBROUTINE setup_exp_constraints 824 825! ************************************************************************************************** 826!> \brief put the constraint information in type and process if requires 827!> BOUNDARIES constraint gets transformed into MAX_VAR_FRACTION constraint. 828!> \param flex_basis ... 829!> \param iset ... 830!> \param ipgf ... 831!> \param const_section ... 832!> \param is_bound ... 833!> \param is_varlim ... 834!> \param irep ... 835!> \author Florian Schiffmann 836! ************************************************************************************************** 837 838 SUBROUTINE set_constraint(flex_basis, iset, ipgf, const_section, is_bound, is_varlim, irep) 839 TYPE(flex_basis_type) :: flex_basis 840 INTEGER :: iset, ipgf 841 TYPE(section_vals_type), POINTER :: const_section 842 LOGICAL :: is_bound, is_varlim 843 INTEGER :: irep 844 845 CHARACTER(len=*), PARAMETER :: routineN = 'set_constraint', routineP = moduleN//':'//routineN 846 847 REAL(KIND=dp) :: r_val 848 REAL(KIND=dp), DIMENSION(:), POINTER :: r_vals 849 850 IF (flex_basis%subset(iset)%exp_has_const(ipgf)) & 851 CALL cp_abort(__LOCATION__, & 852 "Multiple constraints due to collision in CONSTRAIN_EXPONENTS."// & 853 " Please change input.") 854 flex_basis%subset(iset)%exp_has_const(ipgf) = .TRUE. 855 IF (is_bound) THEN 856 flex_basis%subset(iset)%exp_const(ipgf)%const_type = 0 857 CALL section_vals_val_get(const_section, "BOUNDARIES", r_vals=r_vals, i_rep_section=irep) 858 flex_basis%subset(iset)%exp_const(ipgf)%llim = MINVAL(r_vals) 859 flex_basis%subset(iset)%exp_const(ipgf)%ulim = MAXVAL(r_vals) 860 r_val = flex_basis%subset(iset)%exps(ipgf) 861 IF (flex_basis%subset(iset)%exps(ipgf) .GT. MAXVAL(r_vals) .OR. flex_basis%subset(iset)%exps(ipgf) .LT. MINVAL(r_vals)) & 862 CALL cp_abort(__LOCATION__, & 863 "Exponent "//cp_to_string(r_val)// & 864 " declared in constraint is out of bounds of constraint"//cp_to_string(MINVAL(r_vals))// & 865 " to"//cp_to_string(MAXVAL(r_vals))// & 866 " Please change input.") 867 flex_basis%subset(iset)%exp_const(ipgf)%init = SUM(r_vals)/2.0_dp 868 flex_basis%subset(iset)%exp_const(ipgf)%var_fac = MAXVAL(r_vals)/flex_basis%subset(iset)%exp_const(ipgf)%init - 1.0_dp 869 END IF 870 IF (is_varlim) THEN 871 flex_basis%subset(iset)%exp_const(ipgf)%const_type = 1 872 CALL section_vals_val_get(const_section, "MAX_VAR_FRACTION", r_vals=r_vals, i_rep_section=irep) 873 flex_basis%subset(iset)%exp_const(ipgf)%var_fac = r_vals(1) 874 flex_basis%subset(iset)%exp_const(ipgf)%init = flex_basis%subset(iset)%exps(ipgf) 875 END IF 876 877 END SUBROUTINE set_constraint 878 879! ************************************************************************************************** 880!> \brief Initialize the optimization vector with the values from the refernece sets 881!> \param x ... 882!> \param basis ... 883!> \param x_ind ... 884!> \author Florian Schiffmann 885! ************************************************************************************************** 886 887 SUBROUTINE assign_x_to_basis(x, basis, x_ind) 888 REAL(KIND=dp), DIMENSION(:) :: x 889 TYPE(flex_basis_type) :: basis 890 INTEGER :: x_ind 891 892 CHARACTER(len=*), PARAMETER :: routineN = 'assign_x_to_basis', & 893 routineP = moduleN//':'//routineN 894 895 INTEGER :: icont, ipgf, iset 896 897 DO iset = 1, basis%nsets 898 DO ipgf = 1, basis%subset(iset)%nexp 899 IF (basis%subset(iset)%opt_exps(ipgf)) THEN 900 x_ind = x_ind + 1 901 basis%subset(iset)%exp_x_ind(ipgf) = x_ind 902 x(x_ind) = basis%subset(iset)%exps(ipgf) 903 END IF 904 DO icont = 1, basis%subset(iset)%ncon_tot 905 IF (basis%subset(iset)%opt_coeff(ipgf, icont)) THEN 906 x_ind = x_ind + 1 907 basis%subset(iset)%coeff_x_ind(ipgf, icont) = x_ind 908 x(x_ind) = basis%subset(iset)%coeff(ipgf, icont) 909 END IF 910 END DO 911 END DO 912 END DO 913 914 END SUBROUTINE assign_x_to_basis 915 916! ************************************************************************************************** 917!> \brief Fill the reference set and get the free varialbles from input 918!> \param kind1_section ... 919!> \param flex_basis ... 920!> \param template_basis_file ... 921!> \param element ... 922!> \param basis_name ... 923!> \param para_env ... 924!> \param ikind ... 925!> \author Florian Schiffmann 926! ************************************************************************************************** 927 928 SUBROUTINE fill_basis_template(kind1_section, flex_basis, template_basis_file, element, basis_name, para_env, ikind) 929 TYPE(section_vals_type), POINTER :: kind1_section 930 TYPE(flex_basis_type) :: flex_basis 931 CHARACTER(LEN=default_path_length) :: template_basis_file 932 CHARACTER(LEN=default_string_length) :: element, basis_name 933 TYPE(cp_para_env_type), POINTER :: para_env 934 INTEGER :: ikind 935 936 CHARACTER(len=*), PARAMETER :: routineN = 'fill_basis_template', & 937 routineP = moduleN//':'//routineN 938 939 INTEGER :: icont, idof, ipgf, irep, iset, nrep 940 INTEGER, DIMENSION(:), POINTER :: switch 941 942 CALL parse_basis(flex_basis, template_basis_file, element, basis_name, para_env) 943 944 ! get the optimizable parameters. Many way to modify them but in the end only logical matrix 945 ! is either set or values get flipped according to the input 946 CALL section_vals_val_get(kind1_section, "INITIAL_DEGREES_OF_FREEDOM", i_val=idof, & 947 i_rep_section=ikind) 948 DO iset = 1, flex_basis%nsets 949 SELECT CASE (idof) 950 CASE (do_opt_none) 951 ! initialization in parse subset did the job 952 CASE (do_opt_all) 953 flex_basis%subset(iset)%opt_coeff = .TRUE. 954 flex_basis%subset(iset)%opt_exps = .TRUE. 955 CASE (do_opt_coeff) 956 flex_basis%subset(iset)%opt_coeff = .TRUE. 957 CASE (do_opt_exps) 958 flex_basis%subset(iset)%opt_exps = .TRUE. 959 CASE DEFAULT 960 CPABORT("No initialization available?????") 961 END SELECT 962 END DO 963 964 CALL section_vals_val_get(kind1_section, "SWITCH_CONTRACTION_STATE", n_rep_val=nrep, i_rep_section=ikind) 965 DO irep = 1, nrep 966 CALL section_vals_val_get(kind1_section, "SWITCH_CONTRACTION_STATE", i_rep_val=irep, & 967 i_rep_section=ikind, i_vals=switch) 968 icont = convert_l_contr_to_entry(flex_basis%subset(switch(1))%lmin, flex_basis%subset(switch(1))%l, switch(3), switch(2)) 969 DO ipgf = 1, flex_basis%subset(switch(1))%nexp 970 flex_basis%subset(switch(1))%opt_coeff(ipgf, icont) = .NOT. flex_basis%subset(switch(1))%opt_coeff(ipgf, icont) 971 END DO 972 END DO 973 974 CALL section_vals_val_get(kind1_section, "SWITCH_COEFF_STATE", n_rep_val=nrep, i_rep_section=ikind) 975 DO irep = 1, nrep 976 CALL section_vals_val_get(kind1_section, "SWITCH_COEFF_STATE", i_rep_val=irep, & 977 i_rep_section=ikind, i_vals=switch) 978 icont = convert_l_contr_to_entry(flex_basis%subset(switch(1))%lmin, flex_basis%subset(switch(1))%l, switch(3), switch(2)) 979 flex_basis%subset(switch(1))%opt_coeff(switch(4), icont) = & 980 .NOT. flex_basis%subset(switch(1))%opt_coeff(switch(4), icont) 981 END DO 982 983 CALL section_vals_val_get(kind1_section, "SWITCH_EXP_STATE", n_rep_val=nrep, i_rep_section=ikind) 984 DO irep = 1, nrep 985 CALL section_vals_val_get(kind1_section, "SWITCH_EXP_STATE", i_rep_val=irep, & 986 i_rep_section=ikind, i_vals=switch) 987 flex_basis%subset(switch(1))%opt_exps(switch(2)) = .NOT. flex_basis%subset(switch(1))%opt_exps(switch(2)) 988 END DO 989 990 CALL section_vals_val_get(kind1_section, "SWITCH_SET_STATE", n_rep_val=nrep, i_rep_section=ikind) 991 DO irep = 1, nrep 992 CALL section_vals_val_get(kind1_section, "SWITCH_SET_STATE", i_rep_val=irep, & 993 i_rep_section=ikind, i_vals=switch) 994 DO ipgf = 1, flex_basis%subset(switch(2))%nexp 995 SELECT CASE (switch(1)) 996 CASE (0) ! switch all states in the set 997 DO icont = 1, flex_basis%subset(switch(2))%ncon_tot 998 flex_basis%subset(switch(2))%opt_coeff(ipgf, icont) = & 999 .NOT. flex_basis%subset(switch(2))%opt_coeff(ipgf, icont) 1000 END DO 1001 flex_basis%subset(switch(2))%opt_exps(ipgf) = .NOT. flex_basis%subset(switch(2))%opt_exps(ipgf) 1002 CASE (1) ! switch only exp 1003 flex_basis%subset(switch(2))%opt_exps(ipgf) = .NOT. flex_basis%subset(switch(2))%opt_exps(ipgf) 1004 CASE (2) ! switch only coeff 1005 DO icont = 1, flex_basis%subset(switch(2))%ncon_tot 1006 flex_basis%subset(switch(2))%opt_coeff(ipgf, icont) = & 1007 .NOT. flex_basis%subset(switch(2))%opt_coeff(ipgf, icont) 1008 END DO 1009 CASE DEFAULT 1010 CPABORT("Invalid option in SWITCH_SET_STATE, 1st value has to be 0, 1 or 2") 1011 END SELECT 1012 END DO 1013 END DO 1014 1015 ! perform a final modification. If basis set is uncontracted coefficient will never have to be optimized 1016 DO irep = 1, flex_basis%nsets 1017 IF (flex_basis%subset(irep)%nexp == 1) THEN 1018 DO ipgf = 1, flex_basis%subset(irep)%nexp 1019 flex_basis%subset(irep)%opt_coeff(ipgf, 1) = .FALSE. 1020 END DO 1021 END IF 1022 END DO 1023 1024 ! finally count the total number of free parameters 1025 flex_basis%nopt = 0 1026 DO irep = 1, flex_basis%nsets 1027 DO ipgf = 1, flex_basis%subset(irep)%nexp 1028 DO icont = 1, flex_basis%subset(irep)%ncon_tot 1029 IF (flex_basis%subset(irep)%opt_coeff(ipgf, icont)) flex_basis%nopt = flex_basis%nopt + 1 1030 END DO 1031 IF (flex_basis%subset(irep)%opt_exps(ipgf)) flex_basis%nopt = flex_basis%nopt + 1 1032 END DO 1033 END DO 1034 1035 END SUBROUTINE fill_basis_template 1036 1037! ************************************************************************************************** 1038!> \brief Helper function to parse input. Converts l and index position of 1039!> a contraction to index in the contraction array of the set using lmin and nl 1040!> \param lmin ... 1041!> \param nl ... 1042!> \param icontr ... 1043!> \param l ... 1044!> \return ... 1045!> \author Florian Schiffmann 1046! ************************************************************************************************** 1047 1048 FUNCTION convert_l_contr_to_entry(lmin, nl, icontr, l) RESULT(ientry) 1049 INTEGER :: lmin 1050 INTEGER, DIMENSION(:) :: nl 1051 INTEGER :: icontr, l, ientry 1052 1053 INTEGER :: i, icon2l, iwork 1054 1055 iwork = l - lmin 1056 icon2l = 0 1057 DO i = 1, iwork 1058 icon2l = icon2l + nl(i) 1059 END DO 1060 ientry = icon2l + icontr 1061 1062 END FUNCTION convert_l_contr_to_entry 1063 1064! ************************************************************************************************** 1065!> \brief Read the reference basis sets from the template basis file 1066!> \param flex_basis ... 1067!> \param template_basis_file ... 1068!> \param element ... 1069!> \param basis_name ... 1070!> \param para_env ... 1071!> \author Florian Schiffmann 1072! ************************************************************************************************** 1073 1074 SUBROUTINE parse_basis(flex_basis, template_basis_file, element, basis_name, para_env) 1075 TYPE(flex_basis_type) :: flex_basis 1076 CHARACTER(LEN=default_path_length) :: template_basis_file 1077 CHARACTER(LEN=default_string_length) :: element, basis_name 1078 TYPE(cp_para_env_type), POINTER :: para_env 1079 1080 CHARACTER(len=*), PARAMETER :: routineN = 'parse_basis', routineP = moduleN//':'//routineN 1081 1082 CHARACTER(LEN=240) :: line 1083 CHARACTER(LEN=242) :: line2 1084 CHARACTER(LEN=LEN(basis_name)+2) :: basis_name2 1085 CHARACTER(LEN=LEN(element)+2) :: element2 1086 INTEGER :: iset, strlen1, strlen2 1087 LOGICAL :: basis_found, found, match 1088 TYPE(cp_parser_type), POINTER :: parser 1089 1090 basis_found = .FALSE. 1091 CALL uppercase(element) 1092 CALL uppercase(basis_name) 1093 NULLIFY (parser) 1094 CALL parser_create(parser, template_basis_file, para_env=para_env) 1095 1096 search_loop: DO 1097 CALL parser_search_string(parser, TRIM(basis_name), .TRUE., found, line) 1098 IF (found) THEN 1099 match = .FALSE. 1100 CALL uppercase(line) 1101 ! Check both the element symbol and the basis set name 1102 line2 = " "//line//" " 1103 element2 = " "//TRIM(element)//" " 1104 basis_name2 = " "//TRIM(basis_name)//" " 1105 strlen1 = LEN_TRIM(element2) + 1 1106 strlen2 = LEN_TRIM(basis_name2) + 1 1107 IF ((INDEX(line2, element2(:strlen1)) > 0) .AND. & 1108 (INDEX(line2, basis_name2(:strlen2)) > 0)) match = .TRUE. 1109 IF (match) THEN 1110 CALL parser_get_object(parser, flex_basis%nsets, newline=.TRUE.) 1111 ALLOCATE (flex_basis%subset(flex_basis%nsets)) 1112 DO iset = 1, flex_basis%nsets 1113 CALL parse_subset(parser, flex_basis%subset(iset)) 1114 END DO 1115 basis_found = .TRUE. 1116 EXIT 1117 END IF 1118 ELSE 1119 EXIT search_loop 1120 END IF 1121 END DO search_loop 1122 CALL parser_release(parser) 1123 1124 IF (.NOT. basis_found) CALL cp_abort(__LOCATION__, & 1125 "The requested basis set <"//TRIM(basis_name)// & 1126 "> for element <"//TRIM(element)//"> was not "// & 1127 "found in the template basis set file ") 1128 1129 END SUBROUTINE parse_basis 1130 1131! ************************************************************************************************** 1132!> \brief Read the subset information from the template basis file 1133!> \param parser ... 1134!> \param subset ... 1135!> \author Florian Schiffmann 1136! ************************************************************************************************** 1137 1138 SUBROUTINE parse_subset(parser, subset) 1139 TYPE(cp_parser_type), POINTER :: parser 1140 TYPE(subset_type) :: subset 1141 1142 CHARACTER(len=*), PARAMETER :: routineN = 'parse_subset', routineP = moduleN//':'//routineN 1143 1144 CHARACTER(len=20*default_string_length) :: line_att 1145 INTEGER :: icon1, icon2, il, ipgf, ishell, istart 1146 REAL(KIND=dp) :: gs_scale 1147 REAL(KIND=dp), POINTER :: r_val 1148 1149 line_att = "" 1150 CALL parser_get_object(parser, subset%n, newline=.TRUE.) 1151 CALL parser_get_object(parser, subset%lmin) 1152 CALL parser_get_object(parser, subset%lmax) 1153 CALL parser_get_object(parser, subset%nexp) 1154 subset%nl = subset%lmax - subset%lmin + 1 1155 ALLOCATE (r_val) 1156 ALLOCATE (subset%l(subset%nl)) 1157 ALLOCATE (subset%exps(subset%nexp)) 1158 ALLOCATE (subset%exp_has_const(subset%nexp)) 1159 subset%exp_has_const = .FALSE. 1160 ALLOCATE (subset%opt_exps(subset%nexp)) 1161 subset%opt_exps = .FALSE. 1162 ALLOCATE (subset%exp_const(subset%nexp)) 1163 ALLOCATE (subset%exp_x_ind(subset%nexp)) 1164 DO ishell = 1, subset%nl 1165 CALL parser_get_object(parser, subset%l(ishell)) 1166 END DO 1167 subset%ncon_tot = SUM(subset%l) 1168 ALLOCATE (subset%coeff(subset%nexp, subset%ncon_tot)) 1169 ALLOCATE (subset%opt_coeff(subset%nexp, subset%ncon_tot)) 1170 subset%opt_coeff = .FALSE. 1171 ALLOCATE (subset%coeff_x_ind(subset%nexp, subset%ncon_tot)) 1172 DO ipgf = 1, subset%nexp 1173 CALL parser_get_object(parser, r_val, newline=.TRUE.) 1174 subset%exps(ipgf) = r_val 1175 DO ishell = 1, subset%ncon_tot 1176 CALL parser_get_object(parser, r_val) 1177 subset%coeff(ipgf, ishell) = r_val 1178 END DO 1179 END DO 1180 1181 ! orthonormalize contraction coefficients using gram schmidt 1182 istart = 1 1183 DO il = 1, subset%nl 1184 DO icon1 = istart, istart + subset%l(il) - 2 1185 DO icon2 = icon1 + 1, istart + subset%l(il) - 1 1186 gs_scale = DOT_PRODUCT(subset%coeff(:, icon2), subset%coeff(:, icon1))/ & 1187 DOT_PRODUCT(subset%coeff(:, icon1), subset%coeff(:, icon1)) 1188 subset%coeff(:, icon2) = subset%coeff(:, icon2) - gs_scale*subset%coeff(:, icon1) 1189 END DO 1190 END DO 1191 istart = istart + subset%l(il) 1192 END DO 1193 1194 ! just to get an understandable basis normalize coefficients 1195 DO icon1 = 1, subset%ncon_tot 1196 subset%coeff(:, icon1) = subset%coeff(:, icon1)/ & 1197 SQRT(DOT_PRODUCT(subset%coeff(:, icon1), subset%coeff(:, icon1))) 1198 END DO 1199 DEALLOCATE (r_val) 1200 1201 END SUBROUTINE parse_subset 1202 1203! ************************************************************************************************** 1204!> \brief Initialize the variables for the powell optimizer 1205!> \param p_param ... 1206!> \param powell_section ... 1207!> \author Florian Schiffmann 1208! ************************************************************************************************** 1209 1210 SUBROUTINE init_powell_var(p_param, powell_section) 1211 TYPE(opt_state_type) :: p_param 1212 TYPE(section_vals_type), POINTER :: powell_section 1213 1214 CHARACTER(len=*), PARAMETER :: routineN = 'init_powell_var', & 1215 routineP = moduleN//':'//routineN 1216 1217 p_param%state = 0 1218 p_param%nvar = 0 1219 p_param%iprint = 0 1220 p_param%unit = default_output_unit 1221 CALL section_vals_val_get(powell_section, "ACCURACY", r_val=p_param%rhoend) 1222 CALL section_vals_val_get(powell_section, "STEP_SIZE", r_val=p_param%rhobeg) 1223 CALL section_vals_val_get(powell_section, "MAX_FUN", i_val=p_param%maxfun) 1224 1225 END SUBROUTINE init_powell_var 1226 1227END MODULE optimize_basis_utils 1228