1!dalton_copyright_start 2! 3! 4!dalton_copyright_end 5 6module lucita_ci_task_interface 7 8 use lucita_integral_density_interface 9 10#ifdef VAR_MPI 11#ifdef USE_MPI_MOD_F90 12 use mpi 13 implicit none 14#else 15 implicit none 16#include "mpif.h" 17#endif 18#else 19 implicit none 20#endif 21 22 public CI_task_list_interface 23 public create_CI_task_list 24 25 private 26 27#ifdef VAR_MPI 28 integer(kind=MPI_INTEGER_KIND) :: my_MPI_REAL8 = MPI_REAL8 29 integer(kind=MPI_INTEGER_KIND) :: ierr_mpi, len_mpi, root_comm=0 30#endif 31 32contains 33 34!********************************************************************** 35 36 subroutine CI_task_list_interface(ci_task_list, & 37 cref, & 38 hc, & 39 resolution_mat, & 40 int1_or_rho1, & 41 int2_or_rho2, & 42 block_list, & 43 par_dist_block_list, & 44 proclist, & 45 grouplist, & 46 fh_array, & 47 rcctos, & 48 nbatch, & 49 nblock, & 50 print_lvl) 51!------------------------------------------------------------------------------- 52! 53! purpose: interface routine to all individual CI tasks defined in a proper order 54! in the list "ci_task_list". 55! 56!------------------------------------------------------------------------------- 57 use lucita_mcscf_ci_cfg 58 character (len=12), intent(in) :: ci_task_list(*) 59!------------- optional input depending on MCSCF/CI run ------------------------ 60 real(8), intent(inout) :: cref(*) 61 real(8), intent(inout) :: hc(*) 62 real(8), intent(inout) :: resolution_mat(*) 63 real(8), intent(inout) :: int1_or_rho1(*) 64 real(8), intent(inout) :: int2_or_rho2(*) 65!------------- end of optional input depending on MCSCF/CI run ----------------- 66 integer, intent(in) :: nbatch 67 integer, intent(in) :: nblock 68 integer, intent(in) :: block_list(nblock) 69 integer, intent(in) :: par_dist_block_list(nblock) 70 integer, intent(in) :: proclist(*) 71 integer, intent(in) :: grouplist(*) 72 integer, intent(in) :: fh_array(*) 73 integer, intent(in) :: rcctos(nblock) 74 integer, intent(in) :: print_lvl 75!------------------------------------------------------------------------------- 76 integer :: ci_task_ticket 77!------------------------------------------------------------------------------- 78 79! initialize ci_task_ticket counter 80 ci_task_ticket = 0 81 82 select case(ci_task_list(ci_task_ticket+1)) 83 84 case ('return CIdia', 'perform Davi', 'return sigma', 'return rotVC', 'return densM', 'ijkl resort ', & 85 'fci dump ') 86 87! fill integral indices pointers and possibly read in integrals 88! ------------------------------------------------------------- 89 call lucita_pointer_integral_driver(ci_task_list(ci_task_ticket+1), & 90 int1_or_rho1, & 91 int2_or_rho2, & 92 mcscf_ci_update_ijkl, & 93 mcscf_orbital_trial_vector, & 94 mcscf_ci_trial_vector, & 95 integrals_from_mcscf_env, & 96 print_lvl) 97 98 case default ! 'report CIspc', 'report CIana' 99 100! print *, ' skipped integral pointer construction and vector block allocation' 101 102 end select 103 104! 105 do ! infinte loop over CI task tickets 106 107 ci_task_ticket = ci_task_ticket + 1 108 109 select case(ci_task_list(ci_task_ticket)) 110 111 112 case ('report CIspc') 113 exit ! return because there is nothing else to do 114 case ('return CIdia') 115 call calculate_CI_diagonal(cref,print_lvl,block_list,par_dist_block_list) 116 case ('perform Davi') 117 call davidson_ci_driver(block_list,par_dist_block_list,proclist, & 118 grouplist,fh_array,rcctos,nblock,print_lvl, & 119 cref,hc,resolution_mat,docisrdft_mc2lu) 120 case ('return sigma') 121 call return_sigma_vector(print_lvl,nbatch,block_list,par_dist_block_list, & 122 rcctos,grouplist,proclist, & 123 cref,hc,resolution_mat,int1_or_rho1,int2_or_rho2) 124 case ('return rotVC') 125 call return_bvec_transformed_2_new_mo_basis(cref,hc,resolution_mat,int1_or_rho1, & 126 par_dist_block_list,block_list, & 127 rcctos,grouplist,proclist) 128 case ('report CIana') 129 call report_CI_vector_analysis(cref,print_lvl) 130 case ('return densM') 131 call return_Xp_density_matrix(print_lvl,nbatch,block_list,par_dist_block_list, & 132 rcctos,grouplist,proclist, & 133 cref,hc,resolution_mat,int1_or_rho1,int2_or_rho2) 134 case default ! standard loop exit criterion 135 exit 136 end select 137 138 end do 139 140 end subroutine CI_task_list_interface 141!********************************************************************** 142 143 subroutine create_CI_task_list(lucita_ci_run_id, & 144 report_ci_analysis, & 145 return_density_matrices, & 146 ci_task_list, & 147 max_ci_tasks) 148!------------------------------------------------------------------------------- 149! 150! purpose: translate external (MCSCF) CI run ID's in LUCITA internal 151! CI tasks. add additional CI tasks as requested by input 152! parameters. 153! 154!------------------------------------------------------------------------------- 155 character (len=12), intent(in) :: lucita_ci_run_id 156 integer , intent(in) :: report_ci_analysis 157 integer , intent(in) :: return_density_matrices 158 integer , intent(in) :: max_ci_tasks 159 character (len=12), intent(inout) :: ci_task_list(*) 160!------------------------------------------------------------------------------- 161 integer :: ci_task_ticket 162!------------------------------------------------------------------------------- 163 164! initialize 165 do ci_task_ticket = 1, max_ci_tasks 166 ci_task_list(ci_task_ticket) = 'undefined ' 167 end do 168 169 ci_task_ticket = 1 170 171 select case(lucita_ci_run_id) 172 173 case('return CIdim') ! calculate # of determinants per symmetry irrep 174 ci_task_list(ci_task_ticket) = 'report CIspc' 175 176 case('return CIdia') ! calculate the diagonal part of the CI Hamiltonian matrix 177 ci_task_list(ci_task_ticket) = 'return CIdia' 178 179 case('sigma vec ') ! calculate a sigma vector 180 ci_task_list(ci_task_ticket) = 'return sigma' 181 182 case('Xp-density m') ! calculate the 1-/2-particle density matrix 183 ci_task_list(ci_task_ticket) = 'return densM' 184 185 case('analyze Cvec') ! analyze CI vector 186 ci_task_list(ci_task_ticket) = 'report CIana' 187 188 case('rotate Cvec') ! rotate CI vector 189 ci_task_list(ci_task_ticket) = 'return rotVC' 190 191 case('srdft ci ') ! srdft ci 192 ci_task_list(ci_task_ticket) = 'return CIdia' 193 ci_task_ticket = ci_task_ticket + 1 194 ci_task_list(ci_task_ticket) = 'perform Davi' 195 196 if(report_ci_analysis > 0)then 197 ci_task_ticket = ci_task_ticket + 1 198 ci_task_list(ci_task_ticket) = 'report CIana' 199 end if 200 201 if(return_density_matrices > 0)then 202 ci_task_ticket = ci_task_ticket + 1 203 ci_task_list(ci_task_ticket) = 'return densM' 204 end if 205 206 case('standard ci ', 'initial ci ') ! perform Davidson CI run 207 208 ci_task_list(ci_task_ticket) = 'return CIdia' 209 ci_task_ticket = ci_task_ticket + 1 210 ci_task_list(ci_task_ticket) = 'perform Davi' 211 212 if(report_ci_analysis > 0)then 213 ci_task_ticket = ci_task_ticket + 1 214 ci_task_list(ci_task_ticket) = 'report CIana' 215 end if 216 217 if(return_density_matrices > 0)then 218 ci_task_ticket = ci_task_ticket + 1 219 ci_task_list(ci_task_ticket) = 'return densM' 220 end if 221 222 case('ijkl resort ') ! resort integrals to lucita format (prior to a large-scale CI) 223 ci_task_list(ci_task_ticket) = 'ijkl resort ' 224 225 case('fci dump ') ! resort integrals to lucita format (prior to a large-scale CI) 226 ci_task_list(ci_task_ticket) = 'fci dump ' 227 228 case default 229 call quit('error in create_CI_task_list: undefined CI run id.') 230 end select 231 232 end subroutine create_CI_task_list 233!********************************************************************** 234 235 subroutine calculate_CI_diagonal(cref, & 236 print_lvl, & 237 block_list, & 238 par_dist_block_list) 239 240 use lucita_energy_types 241 use file_type_module 242 use ttss_block_module 243 244#ifdef VAR_MPI 245 use par_mcci_io 246#endif 247 248#include "parluci.h" 249#include "mxpdim.inc" 250#include "crun.inc" 251#include "oper.inc" 252#include "clunit.inc" 253!------------------------------------------------------------------------------- 254 real(8), intent(inout) :: cref(*) 255 integer, intent(in) :: print_lvl 256 integer, intent(in) :: block_list(num_blocks) 257 integer, intent(in) :: par_dist_block_list(num_blocks) 258!------------------------------------------------------------------------------- 259#ifdef VAR_MPI 260#include "maxorb.h" 261#include "infpar.h" 262 integer(MPI_INTEGER_KIND):: my_STATUS(MPI_STATUS_SIZE) 263 integer(MPI_OFFSET_KIND) :: offset 264 integer :: block_length 265 integer :: isblk 266 integer :: iproc 267 integer :: my_MPI_COMM_WORLD = MPI_COMM_WORLD 268#endif 269!------------------------------------------------------------------------------- 270 271! construct the diagonal part of the Hamiltonian 272 if(idiag == 2) ludia = lusc1 273 if(icistr >= 2) call rewino(ludia) 274 275 call gasdiat(cref,ludia,ecore_orig-ecore,icistr,i12, & 276 ttss_info%ttss_block_type,num_blocks,ttss_info%ttss_vec_split, & 277 par_dist_block_list) 278 279 if(nocsf == 1 .and. icistr == 1)then 280 call rewino(ludia) 281 call todsc_luci(cref,l_combi,l_combi,ludia) 282 end if 283 284 if(luci_nmproc > 1)then 285#ifdef VAR_MPI 286 rewind ludia 287 288 call isetvc(file_info%iluxlist(1,file_info%current_file_nr_diag),1, & 289 file_info%max_list_length) 290 291! collect H diagonal 292 call mcci_cp_vcd_mpi_2_seq_io_interface(cref,ludia, & 293 file_info%fh_lu(file_info% & 294 current_file_nr_diag), & 295 file_info%file_offsets( & 296 file_info%current_file_nr_diag), & 297 file_info%iluxlist(1, & 298 file_info%current_file_nr_diag), & 299 par_dist_block_list, & 300 block_list, & 301 my_MPI_COMM_WORLD, & 302 num_blocks,1,1,2) 303#endif 304 end if 305 306#ifdef LUCI_DEBUG 307! debug printing (parallel run) 308 if(print_lvl .ge. 10)then 309#ifdef VAR_MPI 310 if(luci_nmproc > 1)then 311 offset = file_info%file_offsets(file_info%current_file_nr_diag) 312 do isblk = 1, num_blocks 313 call get_block_proc(par_dist_block_list,isblk,iproc) 314 if(luci_myproc == iproc)then 315 call get_block_length(block_list,isblk,block_length) 316 call mpi_file_read_at(file_info%fh_lu(file_info%current_file_nr_diag),& 317 offset,cref,block_length,my_MPI_REAL8,my_STATUS,ierr_mpi) 318 write(luwrt,*) ' # of diagonal elements ==> ',block_length 319 call wrtmatmn(cref,1,block_length,1,block_length,luwrt) 320 offset = offset + block_length 321 end if 322 end do 323 end if 324#endif 325 end if 326#endif 327 328 end subroutine calculate_CI_diagonal 329!********************************************************************** 330 331 subroutine report_CI_vector_analysis(cref,print_lvl) 332 333 use file_type_module 334 use ttss_block_module 335 336#include "parluci.h" 337#include "mxpdim.inc" 338#include "crun.inc" 339#include "cstate.inc" 340#include "clunit.inc" 341!------------------------------------------------------------------------------- 342 real(8), intent(inout) :: cref(*) 343 integer, intent(in) :: print_lvl 344!------------------------------------------------------------------------------- 345 integer :: eigen_state_id 346 integer :: imzero, iampack 347 integer :: luc_vector_file 348 integer :: lusc_vector_file 349!------------------------------------------------------------------------------- 350 351 if(luci_myproc == luci_master)then 352 353 luc_vector_file = luc 354! MCSCF - CI task run: file handle set in MCSCF/CI interface routine 355 if(file_info%current_file_fh_seqf(1) > 0) luc_vector_file = file_info%current_file_fh_seqf(1) 356 357 call rewino(luc_vector_file) 358 359 do eigen_state_id = 1, nroot 360 361 if(icistr == 1)then 362 call frmdsc_luci(cref,l_combi,l_combi,luc_vector_file,imzero,iampack) 363 else 364 lusc_vector_file = luc_vector_file 365 if(nroot > 1)then 366 call rewino(lusc1) 367 call copvcd(luc_vector_file,lusc1,cref,0,-1) 368 lusc_vector_file = lusc1 369 end if 370 end if 371! 372 write(luwrt,'(/a)') ' ******************************************************' 373 write(luwrt,'(a,i3)') ' Analysis of leading coefficients for eigen state = ',eigen_state_id 374 write(luwrt,'(a)') ' ******************************************************' 375 376 377 call rewino(lusc_vector_file) 378 call gasana(l0block,num_blocks,ttss_info%ttss_vec_split,ttss_info%ttss_block_type,lusc_vector_file,icistr) 379 end do 380 end if 381 382 end subroutine report_CI_vector_analysis 383!********************************************************************** 384 385 subroutine return_Xp_density_matrix(print_lvl, & 386 nbatch, & 387 block_list, & 388 par_dist_block_list, & 389 rcctos, & 390 grouplist, & 391 proclist, & 392 cref, & 393 hc, & 394 resolution_mat, & 395 int1_or_rho1, & 396 int2_or_rho2) 397 use file_type_module 398 use lucita_cfg 399 use lucita_mcscf_ci_cfg 400 use lucita_mcscf_srdftci_cfg 401#ifdef VAR_MPI 402 use par_mcci_io 403#endif 404 405 406#include "priunit.h" 407#include "parluci.h" 408#include "mxpdim.inc" 409#include "crun.inc" 410#include "cstate.inc" 411#include "cprnt.inc" 412#include "clunit.inc" 413#include "oper.inc" 414#include "orbinp.inc" 415#include "lucinp.inc" 416#include "cicisp.inc" 417#include "cands.inc" 418#include "glbbas.inc" 419!------------------------------------------------------------------------------- 420 integer, intent(in) :: print_lvl 421 integer, intent(in) :: nbatch 422 integer, intent(in) :: block_list(num_blocks) 423 integer, intent(in) :: par_dist_block_list(num_blocks) 424 integer, intent(in) :: rcctos(num_blocks) 425 integer, intent(in) :: grouplist(luci_nmproc) 426 integer, intent(in) :: proclist(luci_nmproc) 427!------------------------------------------------------------------------------- 428 real(8), intent(inout) :: cref(*) 429 real(8), intent(inout) :: hc(*) 430 real(8), intent(inout) :: resolution_mat(*) 431 real(8), intent(inout) :: int1_or_rho1(*) 432 real(8), intent(inout) :: int2_or_rho2(*) 433!------------------------------------------------------------------------------- 434 real(8) :: work 435#include "wrkspc.inc" 436 real(8) :: test_energy 437 real(8) :: exps2 438 real(8) :: cv_dummy, hcv_dummy 439 integer, parameter :: isigden = 2 ! density matrix switch for sigden_ci 440 integer :: eigen_state_id 441 integer :: imzero, iampack 442 integer :: lusc_vector_file 443 integer :: luhc_vector_file 444 integer :: luc_internal 445 integer :: luhc_internal 446 integer :: twopart_densdim 447 integer :: rhotype 448 integer :: rhotype_spin1 449 integer :: idum 450 integer(8) :: kdum 451 integer(8) :: k_scratchsbatch = 1 ! to avoid compiler warnings 452 integer(8) :: k_dens2_scratch 453 integer(8) :: k_scratch1 454 integer(8) :: k_scratch2 455 integer(8) :: k_scratch3 456 integer(8) :: k_scratch4 457#ifdef VAR_MPI 458 integer, allocatable :: blocks_per_batch(:) 459 integer, allocatable :: batch_length(:) 460 integer, allocatable :: block_offset_batch(:) 461 integer, allocatable :: block_info_batch(:,:) 462 integer :: iatp, ibtp 463 integer :: nbatch_par, nblock_par 464 integer :: iloop, nbatch_par_max 465 integer(kind=MPI_INTEGER_KIND) :: mynew_comm_mpi, my_luci_master 466 integer(kind=MPI_OFFSET_KIND) :: my_lu4_off_tmp 467 integer :: my_MPI_COMM_WORLD = MPI_COMM_WORLD 468#endif 469!------------------------------------------------------------------------------- 470 471! set level of particle-density matrix calculation 472 i12 = idensi 473 474! initialize s**2 expectation value 475 exps2 = 0.0d0 476 477! set rhotype 478 rhotype = 1 479! remember: lucita-internal routines for natorbs expects an unpacked density matrix 480 if(lucita_cfg_natural_orb_occ_nr) rhotype = 0 481 if(lucita_cfg_transition_densm) rhotype = 3 482 483! MCSCF - CI task run: file handle set in MCSCF/CI interface routine 484 luc_internal = luc 485 luhc_internal = luhc 486 if(file_info%current_file_fh_seqf(1) > 0) luc_internal = file_info%current_file_fh_seqf(1) 487 if(file_info%current_file_fh_seqf(2) > 0) luhc_internal = file_info%current_file_fh_seqf(2) 488 489#ifdef VAR_MPI 490 my_luci_master = luci_master 491 if(luci_nmproc > 1)then 492 493 if(icsm /= issm) & 494 call quit('*** return_Xp_density_matrix: property/response density matrix not parallelized yet. ***') 495 496! sequential --> MPI I/O 497! ---------------------- 498 allocate(blocks_per_batch(mxntts)) 499 allocate(batch_length(mxntts)) 500 allocate(block_offset_batch(mxntts)) 501 allocate(block_info_batch(8,mxntts)) 502 blocks_per_batch = 0 503 batch_length = 0 504 block_offset_batch = 0 505 block_info_batch = 0 506 iatp = 1 507 ibtp = 2 508 call Z_BLKFO_partitioning_parallel(icspc,icsm,iatp,ibtp, & 509 blocks_per_batch,batch_length, & 510 block_offset_batch, & 511 block_info_batch,nbatch_par, & 512 nblock_par,par_dist_block_list) 513 514 if(luci_myproc == luci_master .and. rhotype == 1)then 515 idum = 0 516 call memman(idum,idum,'MARK ',idum,'Xpden0') 517 nbatch_par_max = 0 518 do iloop = 1, nbatch_par 519 nbatch_par_max = max(nbatch_par_max,batch_length(iloop)) 520 end do 521! write(lupri,*) 'allocate scratch space for sigma vec: ',nbatch_par_max 522 call memman(k_scratchsbatch,nbatch_par_max,'ADDL ',2,'SBATCH') 523 call dzero(work(k_scratchsbatch),nbatch_par_max) 524 end if 525 526 527 if(.not.file_info%file_type_mc)then ! not within MCSCF 528 529! hardwired to ILU1 + ILU2 for CI 530 file_info%current_file_nr_active1 = 2 531 file_info%current_file_nr_active2 = 2 532 call rewino(luc_internal) 533 534 call izero(file_info%iluxlist(1,file_info%current_file_nr_active1), & 535 file_info%max_list_length) 536 call izero(file_info%iluxlist(1,file_info%current_file_nr_active2), & 537 file_info%max_list_length) 538 539! step 1: the rhs vector 540 call mcci_cp_vcd_mpi_2_seq_io_interface(cref, & 541 luc_internal, & 542 file_info%fh_lu(file_info% & 543 current_file_nr_active1), & 544 file_info%file_offsets( & 545 file_info%current_file_nr_active1), & 546 file_info%iluxlist(1, & 547 file_info%current_file_nr_active1), & 548 par_dist_block_list, & 549 block_list, & 550 my_MPI_COMM_WORLD, & 551 NUM_BLOCKS,NROOT,1,1) 552 end if 553 554 if(lucita_cfg_transition_densm)then 555 556 if(.not.file_info%file_type_mc)then ! not within MCSCF 557 558! step 2: the lhs vector 559 file_info%current_file_nr_active2 = 3 560 call rewino(luhc_internal) 561 562 call mcci_cp_vcd_mpi_2_seq_io_interface(hc, & 563 luhc_internal, & 564 file_info%fh_lu(file_info% & 565 current_file_nr_active2), & 566 file_info%file_offsets( & 567 file_info%current_file_nr_active2), & 568 file_info%iluxlist(1, & 569 file_info%current_file_nr_active2), & 570 par_dist_block_list, & 571 block_list, & 572 my_MPI_COMM_WORLD, & 573 NUM_BLOCKS,NROOT,1,1) 574 end if 575 576 else 577 578 call icopy(file_info%max_list_length,file_info%iluxlist(1,file_info%current_file_nr_active1),1, & 579 file_info%iluxlist(1,file_info%current_file_nr_active2),1) 580 end if 581 582 luhc_vector_file = file_info%fh_lu(file_info%current_file_nr_active2) 583! save and temporarily re-set the file offset ("ilu4" is used inside sigden_ci) 584 my_lu4_off_tmp = my_lu4_off 585 my_lu4_off = file_info%file_offsets(file_info%current_file_nr_active2) 586 587 lusc_vector_file = file_info%fh_lu(file_info%current_file_nr_bvec) 588 end if 589#endif 590 591 call rewino(luc_internal) 592 call rewino(luhc_internal) 593 594!#define LUCI_DEBUG 595#ifdef LUCI_DEBUG 596 iprden = 100 597 write(luwrt,'(/a)') ' IPRDEN raised explicitly in return_Xp_density_matrix ' 598#endif 599 600 do eigen_state_id = 1, nroot 601 602 write(luwrt,'(/a)') ' **************************************************************' 603 write(luwrt,'(a,i3)') ' calculating 1-/2-particle density matrix for eigen state = ', eigen_state_id 604 write(luwrt,'(a)') ' **************************************************************' 605 606 if(luci_nmproc == 1)then 607 if(icistr == 1)then 608 call frmdsc_luci(cref,l_combi,l_combi,luc_internal,imzero,iampack) 609 if(lucita_cfg_transition_densm)then 610! transition density 611 call frmdsc_luci(hc,l_combi,l_combi,luhc_internal,imzero,iampack) 612 else 613 call dcopy(l_combi,cref,1,hc,1) 614 end if 615 else 616 if(lucita_cfg_transition_densm)then 617! transition density 618 lusc_vector_file = luc_internal 619 luhc_vector_file = luhc_internal 620 else 621 call rewino(lusc1) 622 call copvcd(luc_internal,lusc1,cref,0,-1) 623 624 call copvcd(lusc1,lusc2,cref,1,-1) 625 call rewino(lusc1) 626 call rewino(lusc2) 627!#define LUCI_DEBUG 628#ifdef LUCI_DEBUG 629 CALL REWINE(lusc1,-1) 630 WRITE(LUWRT,*) ' final solution vector for root ==> ',eigen_state_id 631 CALL WRTVCD(cref,LUsC1,0,-1) 632 CALL REWINE(LUsC1,-1) 633#undef LUCI_DEBUG 634#endif 635 lusc_vector_file = lusc1 636 luhc_vector_file = lusc2 637 end if 638 end if 639 else 640#ifdef VAR_MPI 641! MPI I/O --> MPI I/O node-master collection file 642! ----------------------------------------------- 643 mynew_comm_mpi = mynew_comm 644 call mpi_barrier(mynew_comm_mpi,ierr_mpi) ! probably not needed if collection of densities is in action 645 call izero(file_info%ilublist,file_info%max_list_length_bvec) 646 647 call mcci_cp_vcd_batch(file_info%fh_lu(file_info%current_file_nr_active1), & 648 file_info%fh_lu(file_info%current_file_nr_bvec), & 649 cref,nbatch_par,blocks_per_batch, & 650 batch_length,block_offset_batch,block_info_batch, & 651 block_list, & 652 file_info%file_offsets(file_info%current_file_nr_active1),& 653 file_info%file_offsets(file_info%current_file_nr_bvec), & 654 file_info%iluxlist(1,file_info%current_file_nr_active1), & 655 file_info%ilublist, & 656 my_vec2_ioff,my_act_blk2,eigen_state_id-1) 657! offset in MPI I/O file for lhs-vector 658 jvec_sf = eigen_state_id - 1 659#endif 660 end if 661 662! 663! stefan - jan 2012: the following construct is a consequence of the physical memory 664! identity of cref and hc on the master node in case of a parallel mcscf 665! and a regular density matrix run 666! 667 if(luci_nmproc > 1 .and. luci_myproc == luci_master .and. rhotype == 1)then 668 call sigden_ci(cref,work(k_scratchsbatch),resolution_mat,lusc_vector_file,luhc_vector_file,& 669 cv_dummy,hcv_dummy,isigden,rhotype,exps2 & 670#ifdef VAR_MPI 671 ,file_info%ilublist,file_info%iluxlist(1,file_info%current_file_nr_active2), & 672 block_list,par_dist_block_list, & 673 rcctos,grouplist,proclist & 674#endif 675 ) 676 else 677 call sigden_ci(cref,hc,resolution_mat,lusc_vector_file,luhc_vector_file, & 678 cv_dummy,hcv_dummy,isigden,rhotype,exps2 & 679#ifdef VAR_MPI 680 ,file_info%ilublist,file_info%iluxlist(1,file_info%current_file_nr_active2), & 681 block_list,par_dist_block_list, & 682 rcctos,grouplist,proclist & 683#endif 684 ) 685 end if 686 687 688 if(luci_nmproc > 1)then 689#ifdef VAR_MPI 690! collect density matrices 691! ------------------------ 692 if(luci_myproc == luci_master)then 693 len_mpi = nacob**2 694 call mpi_reduce(mpi_in_place,work(krho1),len_mpi,my_mpi_real8, & 695 mpi_sum,my_luci_master,mpi_comm_world,ierr_mpi) 696 if(ispnden > 0)then 697 call mpi_reduce(mpi_in_place,work(ksrho1),len_mpi,my_mpi_real8, & 698 mpi_sum,my_luci_master,mpi_comm_world,ierr_mpi) 699 call mpi_reduce(mpi_in_place,work(ksrho1a),len_mpi,my_mpi_real8, & 700 mpi_sum,my_luci_master,mpi_comm_world,ierr_mpi) 701 call mpi_reduce(mpi_in_place,work(ksrho1b),len_mpi,my_mpi_real8, & 702 mpi_sum,my_luci_master,mpi_comm_world,ierr_mpi) 703 end if 704 if(i12 > 1)then 705 len_mpi = nacob**2*(nacob**2+1)/2 706 call mpi_reduce(mpi_in_place,work(krho2),len_mpi, & 707 my_mpi_real8,mpi_sum,my_luci_master,mpi_comm_world,ierr_mpi) 708 len_mpi = 1 709 call mpi_reduce(mpi_in_place,exps2,len_mpi, & 710 my_mpi_real8,mpi_sum,my_luci_master,mpi_comm_world,ierr_mpi) 711 end if 712 else 713 len_mpi = nacob**2 714 call mpi_reduce(work(krho1),mpi_in_place,len_mpi,my_mpi_real8, & 715 mpi_sum,my_luci_master,mpi_comm_world,ierr_mpi) 716 if(ispnden > 0)then 717 call mpi_reduce(work(ksrho1),mpi_in_place,len_mpi,my_mpi_real8, & 718 mpi_sum,my_luci_master,mpi_comm_world,ierr_mpi) 719 call mpi_reduce(work(ksrho1a),mpi_in_place,len_mpi,my_mpi_real8, & 720 mpi_sum,my_luci_master,mpi_comm_world,ierr_mpi) 721 call mpi_reduce(work(ksrho1b),mpi_in_place,len_mpi,my_mpi_real8, & 722 mpi_sum,my_luci_master,mpi_comm_world,ierr_mpi) 723 end if 724 if(i12 > 1)then 725 len_mpi = nacob**2*(nacob**2+1)/2 726 call mpi_reduce(work(krho2),mpi_in_place,len_mpi, & 727 my_mpi_real8,mpi_sum,my_luci_master,mpi_comm_world,ierr_mpi) 728 len_mpi = 1 729 call mpi_reduce(exps2,mpi_in_place,len_mpi, & 730 my_mpi_real8,mpi_sum,my_luci_master,mpi_comm_world,ierr_mpi) 731 end if 732 end if 733#endif 734 end if ! luci_nmproc > 1 735 736 if(i12 > 1)then 737 write(luwrt,'(/a )') ' ------------------------------------------------' 738 write(luwrt,'(a,1f10.3)') ' expectation value of operator <S**2> =', exps2 739 write(luwrt,'(a/ )') ' ------------------------------------------------' 740 end if 741 742! export 1-/2-particle density matrix to MCSCF format 743! --------------------------------------------------- 744! 745!#ifdef largeFCInotneeded 746 idum = 0 747 call memman(idum,idum,'MARK ',idum,'Xpden1') 748 if(i12 > 1) call memman(k_dens2_scratch,nacob**4,'ADDL ',2,'PVfull') 749 if(i12 == 1) call memman(k_dens2_scratch, 0,'ADDL ',2,'PVfull') 750 751! set rhotype for spin-densities after reordering 752 rhotype_spin1 = 1 753 if(lucita_cfg_natural_orb_occ_nr) rhotype_spin1 = rhotype 754 755! activate for MCSCF/improved CI 756 if(integrals_from_mcscf_env)then 757 758 !> reoder first spin-densities in order to not overwrite the outgoing 1p density matrix in int1_or_rho1 759 if(ispnden == 1)then 760 call lucita_spinputdens_1p(work(kSRHO1a),work(krho2),int1_or_rho1,int2_or_rho2, & 761 work(k_dens2_scratch),i12,isigden,rhotype_spin1,eigen_state_id,1) 762 call lucita_spinputdens_1p(work(kSRHO1b),work(krho2),int1_or_rho1,int2_or_rho2, & 763 work(k_dens2_scratch),i12,isigden,rhotype_spin1,eigen_state_id,2) 764 end if 765 call lucita_putdens_generic(work(krho1),work(krho2),int1_or_rho1,int2_or_rho2, & 766 work(k_dens2_scratch),i12,isigden,rhotype,eigen_state_id) 767 768 else 769 twopart_densdim = 0 770 if(i12 > 1) twopart_densdim = (nacob*(nacob+1)/2)**2 771 call memman(k_scratch1 ,nacob**2 ,'ADDL ',2,'scrat1') 772 call memman(k_scratch2 ,twopart_densdim ,'ADDL ',2,'scrat2') 773 774 call lucita_putdens_generic(work(krho1),work(krho2),work(k_scratch1),work(k_scratch2), & 775 work(k_dens2_scratch),i12,isigden,rhotype,eigen_state_id) 776 if(ispnden >= 1)then 777 call lucita_spinputdens_1p(work(ksrho1a),work(krho2),work(k_scratch1),work(k_scratch2), & 778 work(k_dens2_scratch), 1,isigden,rhotype_spin1,eigen_state_id,1) 779 call lucita_spinputdens_1p(work(ksrho1b),work(krho2),work(k_scratch1),work(k_scratch2), & 780 work(k_dens2_scratch), 1,isigden,rhotype_spin1,eigen_state_id,2) 781 end if 782 end if 783! 784 call memman(kdum ,idum,'FLUSM ',2,'Xpden1') 785!#endif 786 787! natural orbital occupation numbers 788! ---------------------------------- 789 if(lucita_cfg_natural_orb_occ_nr)then 790 idum = 0 791 call memman(idum,idum,'MARK ',idum,'Xpden2') 792 call memman(k_scratch1,nacob**2 ,'ADDL ',2,'scrat1') 793 call memman(k_scratch2,nacob**2 ,'ADDL ',2,'scrat2') 794 call memman(k_scratch3,nacob ,'ADDL ',2,'scrat3') 795 call memman(k_scratch4,nacob*(nacob+1)/2,'ADDL ',2,'scrat4') 796! if (ispnden == 1 .and. lucita_cfg_is_spin_multiplett /= 1) then 797 if (ispnden >= 1)then 798 write(luwrt,'(/a)') ' Natural spin-orbital occupation numbers for alpha spin-orbitals' 799 call lnatorb(work(ksrho1a),nsmob,ntoobs,nacobs,ninobs, & 800 ireost,work(k_scratch1),work(k_scratch2),work(k_scratch3),nacob, & 801 work(k_scratch4),1) 802 write(luwrt,'(/a)') ' Natural spin-orbital occupation numbers for beta spin-orbitals' 803 call lnatorb(work(ksrho1b),nsmob,ntoobs,nacobs,ninobs, & 804 ireost,work(k_scratch1),work(k_scratch2),work(k_scratch3),nacob, & 805 work(k_scratch4),1) 806 end if 807 808 write(luwrt,'(/a)') ' Natural orbital occupation numbers' 809 call lnatorb(work(krho1),nsmob,ntoobs,nacobs,ninobs, & 810 ireost,work(k_scratch1),work(k_scratch2),work(k_scratch3),nacob, & 811 work(k_scratch4),1) 812 813 call memman(kdum ,idum,'FLUSM ',2,'Xpden2') 814 end if 815 816 if(i12 > 1 .and. .not. srdft_ci_with_lucita) call en_from_dens(test_energy,i12) 817 818 end do ! loop over eigen states 819#ifdef LUCI_DEBUG 820 write(luwrt,'(/a)') ' IPRDEN lowered explicitly in return_Xp_density_matrix ' 821 iprden = 1 822#endif 823 824#ifdef VAR_MPI 825 if(luci_nmproc > 1)then 826 if(luci_myproc == luci_master .and. rhotype == 1)then 827 call memman(kdum ,idum,'FLUSM ',2,'Xpden0') 828 end if 829! reset file off-set and free memory 830 if(.not.lucita_cfg_transition_densm) my_lu4_off = my_lu4_off_tmp 831 deallocate(block_info_batch) 832 deallocate(block_offset_batch) 833 deallocate(batch_length) 834 deallocate(blocks_per_batch) 835 end if 836#endif 837 838 end subroutine return_Xp_density_matrix 839!********************************************************************** 840 841 subroutine return_sigma_vector(print_lvl, & 842 nbatch, & 843 block_list, & 844 par_dist_block_list, & 845 rcctos, & 846 grouplist, & 847 proclist, & 848 cref, & 849 hc, & 850 resolution_mat, & 851 int1_or_rho1, & 852 int2_or_rho2) 853 use file_type_module 854 use lucita_cfg 855 use lucita_mcscf_ci_cfg 856#ifdef VAR_MPI 857 use par_mcci_io 858#endif 859 860 861#include "parluci.h" 862#include "mxpdim.inc" 863#include "crun.inc" 864#include "cstate.inc" 865#include "clunit.inc" 866#include "oper.inc" 867#include "orbinp.inc" 868#include "lucinp.inc" 869#include "cicisp.inc" 870#include "cands.inc" 871#include "glbbas.inc" 872!------------------------------------------------------------------------------- 873 integer, intent(in) :: print_lvl 874 integer, intent(in) :: nbatch 875 integer, intent(in) :: block_list(num_blocks) 876 integer, intent(in) :: par_dist_block_list(num_blocks) 877 integer, intent(in) :: rcctos(num_blocks) 878 integer, intent(in) :: grouplist(luci_nmproc) 879 integer, intent(in) :: proclist(luci_nmproc) 880!------------------------------------------------------------------------------- 881 real(8), intent(inout) :: cref(*) 882 real(8), intent(inout) :: hc(*) 883 real(8), intent(inout) :: resolution_mat(*) 884 real(8), intent(inout) :: int1_or_rho1(*) 885 real(8), intent(inout) :: int2_or_rho2(*) 886!------------------------------------------------------------------------------- 887 real(8) :: cv_dummy, hcv_dummy,exps2 888 integer :: eigen_state_id 889 integer :: imzero, iampack 890 integer :: lusc_vector_file 891 integer :: luhc_vector_file 892 integer :: luc_internal 893 integer :: luhc_internal 894 integer, parameter :: isigden = 1 ! sigma vector switch for sigden_ci 895#ifdef VAR_MPI 896 integer, allocatable :: blocks_per_batch(:) 897 integer, allocatable :: batch_length(:) 898 integer, allocatable :: block_offset_batch(:) 899 integer, allocatable :: block_info_batch(:,:) 900 integer :: iatp, ibtp 901 integer :: nbatch_par, nblock_par 902 integer(kind=MPI_INTEGER_KIND) :: mynew_comm_mpi 903 integer(kind=mpi_offset_kind) :: my_lu4_off_tmp 904 integer :: my_MPI_COMM_WORLD = MPI_COMM_WORLD 905#endif 906!------------------------------------------------------------------------------- 907 908! MCSCF - CI task run: file handle set in MCSCF/CI interface routine 909 luc_internal = luc 910 luhc_internal = luhc 911 if(file_info%current_file_fh_seqf(1) > 0) luc_internal = file_info%current_file_fh_seqf(1) 912 913#ifdef VAR_MPI 914 if(luci_nmproc > 1)then 915 916 if(icsm /= issm) & 917 call quit('*** return_sigma_vector: property/response E2b calculation not parallelized yet. ***') 918 919! sequential --> MPI I/O 920! ---------------------- 921 922 allocate(blocks_per_batch(mxntts)) 923 allocate(batch_length(mxntts)) 924 allocate(block_offset_batch(mxntts)) 925 allocate(block_info_batch(8,mxntts)) 926 blocks_per_batch = 0 927 batch_length = 0 928 block_offset_batch = 0 929 block_info_batch = 0 930 iatp = 1 931 ibtp = 2 932 call Z_BLKFO_partitioning_parallel(icspc,icsm,iatp,ibtp, & 933 blocks_per_batch,batch_length, & 934 block_offset_batch, & 935 block_info_batch,nbatch_par, & 936 nblock_par,par_dist_block_list) 937 938! set HC vector file target 939 file_info%current_file_nr_active2 = 3 940 941 if(.not.file_info%file_type_mc)then ! not within MCSCF 942 943! hardwired to ILU1 for CI 944 file_info%current_file_nr_active1 = 2 945 946 call izero(file_info%iluxlist(1,file_info%current_file_nr_active1), & 947 file_info%max_list_length) 948 call izero(file_info%iluxlist(1,file_info%current_file_nr_active2), & 949 file_info%max_list_length) 950 951 call rewino(luc_internal) 952! step 1: the rhs vector 953 call mcci_cp_vcd_mpi_2_seq_io_interface(cref, & 954 luc_internal, & 955 file_info%fh_lu(file_info% & 956 current_file_nr_active1), & 957 file_info%file_offsets( & 958 file_info%current_file_nr_active1), & 959 file_info%iluxlist(1, & 960 file_info%current_file_nr_active1), & 961 par_dist_block_list, & 962 block_list, & 963 my_MPI_COMM_WORLD, & 964 NUM_BLOCKS,NROOT,1,1) 965 end if ! not within MCSCF 966 end if 967#endif 968 969 call rewino(luhc_internal) 970 call rewino(luc_internal) 971 972 do eigen_state_id = 1, nroot 973 974 write(luwrt,'(/a)') ' **********************************************' 975 write(luwrt,'(a,i3)') ' calculating E2b lhs vector for b state = ', eigen_state_id 976 write(luwrt,'(a)') ' **********************************************' 977 978 if(luci_nmproc == 1)then 979 if(icistr == 1)then 980 call frmdsc_luci(cref,l_combi,l_combi,luc_internal,imzero,iampack) 981 else 982 call rewino(lusc1) 983 call copvcd(luc_internal,lusc1,cref,0,-1) 984!#define LUCI_DEBUG 985#ifdef LUCI_DEBUG 986 call wrtvcd(cref,lusc1,-1,-1) 987#endif 988 call rewino(lusc1) 989 end if 990 lusc_vector_file = lusc1 991 luhc_vector_file = luhc_internal 992 else 993#ifdef VAR_MPI 994! MPI I/O --> MPI I/O node-master collection file 995! ----------------------------------------------- 996 mynew_comm_mpi = mynew_comm 997 call mpi_barrier(mynew_comm_mpi,ierr_mpi) 998 call izero(file_info%ilublist,file_info%max_list_length_bvec) 999 1000 call mcci_cp_vcd_batch(file_info%fh_lu(file_info%current_file_nr_active1), & 1001 file_info%fh_lu(file_info%current_file_nr_bvec), & 1002 hc,nbatch_par,blocks_per_batch, & 1003 batch_length,block_offset_batch,block_info_batch, & 1004 block_list, & 1005 file_info%file_offsets(file_info%current_file_nr_active1),& 1006 file_info%file_offsets(file_info%current_file_nr_bvec), & 1007 file_info%iluxlist(1,file_info%current_file_nr_active1), & 1008 file_info%ilublist, & 1009 my_vec2_ioff,my_act_blk2,eigen_state_id-1) 1010 1011 lusc_vector_file = file_info%fh_lu(file_info%current_file_nr_bvec) 1012 luhc_vector_file = file_info%fh_lu(file_info%current_file_nr_active2) 1013 1014! save and temporarily re-set the file offset ("ilu4" is used inside sigden_ci) 1015 my_lu4_off_tmp = my_lu4_off 1016 my_lu4_off = file_info%file_offsets(file_info%current_file_nr_active2) 1017 1018! offset in MPI I/O file for lhs-vector 1019 jvec_sf = eigen_state_id - 1 1020#endif 1021 end if 1022! 1023 call sigden_ci(cref,hc,resolution_mat,lusc_vector_file,luhc_vector_file, & 1024 cv_dummy,hcv_dummy,isigden,-1,exps2 & 1025#ifdef VAR_MPI 1026 ,file_info%ilublist,file_info%iluxlist(1,file_info%current_file_nr_active2),& 1027 block_list,par_dist_block_list, & 1028 rcctos,grouplist,proclist & 1029#endif 1030 ) 1031 1032 1033 end do ! loop over eigen states 1034 1035 if(luci_nmproc > 1)then 1036#ifdef VAR_MPI 1037 if(.not.file_info%file_type_mc)then ! not within MCSCF 1038! collect e2b vector(s) 1039! -------------------- 1040 call rewino(luhc_internal) 1041 1042! the lhs vector 1043 call mcci_cp_vcd_mpi_2_seq_io_interface(hc,luhc_internal,luhc_vector_file, & 1044 file_info%file_offsets( & 1045 file_info%current_file_nr_active2), & 1046 file_info%iluxlist(1, & 1047 file_info%current_file_nr_active2), & 1048 par_dist_block_list, & 1049 block_list, & 1050 my_MPI_COMM_WORLD, & 1051 num_blocks,nroot,1,2) 1052 1053 end if ! not within MCSCF 1054 1055 my_lu4_off = my_lu4_off_tmp 1056 1057 deallocate(block_info_batch) 1058 deallocate(block_offset_batch) 1059 deallocate(batch_length) 1060 deallocate(blocks_per_batch) 1061#endif 1062 end if ! luci_nmproc > 1 1063 1064 end subroutine return_sigma_vector 1065!********************************************************************** 1066 1067 subroutine return_bvec_transformed_2_new_mo_basis(vec1,vec2,c2,mo2mo_mat, & 1068 par_dist_block_list, & 1069 block_list,rcctos, & 1070 grouplist,proclist) 1071! 1072! purpose: master routine for transforming a CI vector to new orbital basis 1073! 1074! Jeppe Olsen, January 98 1075! re-factored and parallelized by Stefan Knecht, November 2011 1076! 1077!------------------------------------------------------------------------------- 1078 use file_type_module 1079!------------------------------------------------------------------------------- 1080#include "mxpdim.inc" 1081#include "lucinp.inc" 1082#include "clunit.inc" 1083#include "cgas.inc" 1084#include "csm.inc" 1085#include "crun.inc" 1086#include "cstate.inc" 1087#include "glbbas.inc" 1088#include "orbinp.inc" 1089#include "cicisp.inc" 1090#include "cands.inc" 1091#include "parluci.h" 1092 INTEGER ADASX,ASXAD,ADSXA,SXSXDX,SXDXSX 1093 COMMON/CSMPRD/ADASX(MXPOBS,MXPOBS),ASXAD(MXPOBS,2*MXPOBS),& 1094 ADSXA(MXPOBS,2*MXPOBS), & 1095 SXSXDX(2*MXPOBS,2*MXPOBS),SXDXSX(2*MXPOBS,4*MXPOBS) 1096 real(8), intent(inout) :: vec1(*) 1097 real(8), intent(inout) :: vec2(*) 1098 real(8), intent(inout) :: c2(*) 1099 real(8), intent(inout) :: mo2mo_mat(*) 1100 integer, intent(in) :: par_dist_block_list(*) 1101 integer, intent(in) :: block_list(*) 1102 integer, intent(in) :: proclist(*) 1103 integer, intent(in) :: grouplist(*) 1104 integer, intent(in) :: rcctos(*) 1105!------------------------------------------------------------------------------- 1106 real(8) :: work 1107#include "wrkspc.inc" 1108 integer :: my_in_fh, my_out_fh, my_sc1_fh, my_sc2_fh 1109 integer :: my_BVC_fh 1110 integer :: len_ilu1 1111 integer :: len_ilu2 1112 integer :: len_iluc 1113 integer :: iatp, ibtp 1114 integer :: nbatch, nblock 1115 integer :: lu_ref, lu_refout 1116 integer :: luc_internal, luhc_internal 1117 integer(kind=8) :: my_in_off 1118 integer(kind=8) :: my_out_off 1119 integer(kind=8) :: my_scr_off 1120 integer(kind=8) :: my_BVC_off 1121 1122 integer, allocatable :: blocks_per_batch(:) 1123 integer, allocatable :: batch_length(:) 1124 integer, allocatable :: block_offset_batch(:) 1125 integer, allocatable :: block_info_batch(:,:) 1126 integer, allocatable :: blocktype(:) 1127 integer(8) :: my_lu4_off_tmp 1128!------------------------------------------------------------------------------- 1129 1130!#define LUCI_DEBUG 1131#ifdef LUCI_DEBUG 1132 WRITE(luwrt,*) ' Welcome to return_bvec_transformed_2_new_mo_basis' 1133 WRITE(luwrt,*) ' =================================================' 1134 call flshfo(luwrt) 1135#endif 1136 1137! MCSCF - CI task run: file handle set in MCSCF/CI interface routine 1138 luc_internal = luc 1139 luhc_internal = luhc 1140 if(file_info%current_file_fh_seqf(1) > 0) luc_internal = file_info%current_file_fh_seqf(1) 1141 1142 len_ilu1 = 0 1143 len_ilu2 = 0 1144 len_iluc = 0 1145 my_in_off = 0 1146 my_out_off = 0 1147 my_scr_off = 0 1148 my_BVC_off = 0 1149 my_BVC_fh = 0 1150#ifdef VAR_MPI 1151 file_info%current_file_nr_active2 = 3 1152 1153 my_in_fh = file_info%fh_lu(file_info%current_file_nr_active1) 1154 my_out_fh = file_info%fh_lu(file_info%current_file_nr_active1) 1155 my_sc1_fh = file_info%fh_lu(file_info%current_file_nr_active2) 1156 my_sc2_fh = file_info%fh_lu(file_info%current_file_nr_active2) 1157 my_BVC_fh = file_info%fh_lu(file_info%current_file_nr_bvec) 1158 1159 my_in_off = file_info%file_offsets( & 1160 file_info%current_file_nr_active1) 1161 my_out_off = file_info%file_offsets( & 1162 file_info%current_file_nr_active1) 1163 my_scr_off = file_info%file_offsets( & 1164 file_info%current_file_nr_active2) 1165 my_BVC_off = file_info%file_offsets( & 1166 file_info%current_file_nr_bvec) 1167 1168 len_ilu1 = file_info%max_list_length 1169 len_ilu2 = file_info%max_list_length 1170 len_iluc = file_info%max_list_length_bvec 1171 1172! save and temporarily re-set the file offset ("ilu4" is used inside sigden_ci) 1173 my_lu4_off_tmp = my_lu4_off 1174 my_lu4_off = file_info%file_offsets(file_info%current_file_nr_active2) 1175 1176! update co-workers with single-orbital transformation matrix 1177 len_mpi = nacob**2 1178 call mpi_bcast(mo2mo_mat,len_mpi,my_mpi_real8,root_comm,mpi_comm_world,ierr_mpi) 1179 1180#else 1181 my_in_fh = lusc1 1182 my_out_fh = luhc_internal 1183 my_sc1_fh = lusc2 1184 my_sc2_fh = lusc3 1185 1186 call copvcd(luc_internal,my_in_fh,vec1,1,-1) 1187 call rewine(my_in_fh,-1) 1188 CALL rewine(my_out_fh,-1) 1189#endif 1190 lu_ref = luc_internal 1191 lu_refout = luhc_internal 1192 1193! set up block and batch structure of vector 1194 1195 allocate(blocks_per_batch(mxntts)) 1196 allocate(batch_length(mxntts)) 1197 allocate(block_offset_batch(mxntts)) 1198 allocate(block_info_batch(8,mxntts)) 1199 allocate(blocktype(nsmst)) 1200 1201 blocks_per_batch = 0 1202 batch_length = 0 1203 block_offset_batch = 0 1204 block_info_batch = 0 1205 blocktype = 0 1206 1207 IATP = 1 1208 IBTP = 2 1209 CALL Z_BLKFO_partitioning_parallel(ICSPC,ICSM,iatp,ibtp, & 1210 blocks_per_batch, & 1211 batch_length, & 1212 block_offset_batch, & 1213 block_info_batch, & 1214 NBATCH,NBLOCK, & 1215 par_dist_block_list) 1216 1217 CALL ZBLTP_IDC(ISMOST(1,ICSM),NSMST,IDC,blocktype) 1218 1219! transform CI vector 1220 call tracid(mo2mo_mat,work(kint1),lu_ref,lu_refout, & 1221 my_in_fh,my_out_fh, & 1222 my_sc1_fh,my_sc2_fh, & 1223 my_BVC_fh, & 1224 vec1,vec2,c2, & 1225 NBATCH,NBLOCK,blocks_per_batch, & 1226 batch_length,block_offset_batch, & 1227 block_info_batch,blocktype, & 1228 par_dist_block_list,block_list, & 1229 rcctos,grouplist,proclist, & 1230 file_info%iluxlist(1,file_info% & 1231 current_file_nr_active1), & 1232 file_info%iluxlist(1,file_info% & 1233 current_file_nr_active2), & 1234 file_info%ilublist, & 1235 len_ilu1,len_ilu2,len_iluc, & 1236 my_in_off,my_out_off,my_scr_off, & 1237 my_BVC_off) 1238 1239#ifdef VAR_MPI 1240! reset 1241 my_lu4_off = my_lu4_off_tmp 1242#endif 1243 1244 deallocate(blocks_per_batch) 1245 deallocate(batch_length) 1246 deallocate(block_offset_batch) 1247 deallocate(block_info_batch) 1248 deallocate(blocktype) 1249 1250 end subroutine return_bvec_transformed_2_new_mo_basis 1251!********************************************************************** 1252 1253end module 1254