TABLE OF CONTENTS
- ABINIT/m_results_out
- m_results_out/copy_results_out
- m_results_out/destroy_results_out
- m_results_out/gather_results_out
- m_results_out/init_results_out
- m_results_out/results_out_type
ABINIT/m_results_out [ Modules ]
NAME
m_results_out
FUNCTION
This module provides the definition of the results_out_type used to store results from GS calculations.
COPYRIGHT
Copyright (C) 2008-2018 ABINIT group (MT) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
PARENTS
CHILDREN
TODO
One should replace the 'pointer' by 'allocatable'. This was tried, in October 2014, but Petrus_nag complained (test v67mbpt t31...t34), and also max2 (paral#08 np=10).
SOURCE
26 #if defined HAVE_CONFIG_H 27 #include "config.h" 28 #endif 29 30 #include "abi_common.h" 31 32 MODULE m_results_out 33 34 use defs_basis 35 use defs_abitypes 36 use m_errors 37 use m_profiling_abi 38 use m_xmpi 39 40 implicit none 41 42 private 43 44 ! public procedures. 45 public :: init_results_out 46 public :: destroy_results_out 47 public :: copy_results_out 48 public :: gather_results_out
m_results_out/copy_results_out [ Functions ]
[ Top ] [ m_results_out ] [ Functions ]
NAME
copy_results_out
FUNCTION
Copy a results_out datastructure into another
INPUTS
results_out_in=<type(results_out_type)>=input results_out datastructure
OUTPUT
results_out_out=<type(results_out_type)>=output results_out datastructure
PARENTS
m_results_out
CHILDREN
copy_results_out,init_results_out,xmpi_allgather,xmpi_allgatherv xmpi_gatherv
SOURCE
465 subroutine copy_results_out(results_out_in,results_out_out) 466 467 468 !This section has been created automatically by the script Abilint (TD). 469 !Do not modify the following lines by hand. 470 #undef ABI_FUNC 471 #define ABI_FUNC 'copy_results_out' 472 !End of the abilint section 473 474 implicit none 475 476 !Arguments ------------------------------------ 477 !arrays 478 type(results_out_type),intent(in) :: results_out_in 479 type(results_out_type),intent(out) :: results_out_out 480 !Local variables------------------------------- 481 !scalars 482 integer :: natom_,natom_out,nimage_,nimage_out,nkpt_,nkpt_out,npsp_,npsp_out,nocc_,nocc_out,ntypat_,ntypat_out 483 484 !************************************************************************ 485 486 !@results_out_type 487 488 nimage_=size(results_out_in%etotal) 489 natom_ =size(results_out_in%fcart,2) 490 nkpt_ =size(results_out_in%npwtot,1) 491 nocc_ =size(results_out_in%occ,1) 492 npsp_ =size(results_out_in%mixalch,1) 493 ntypat_=size(results_out_in%mixalch,2) 494 nimage_out=0;if (associated(results_out_out%etotal))nimage_out=size(results_out_out%etotal) 495 natom_out =0;if (associated(results_out_out%fcart)) natom_out =size(results_out_out%fcart,2) 496 nkpt_out =0;if (associated(results_out_out%npwtot))nkpt_out =size(results_out_out%npwtot,1) 497 nocc_out =0;if (associated(results_out_out%occ)) nocc_out =size(results_out_out%occ,1) 498 npsp_out =0;if (associated(results_out_out%mixalch))npsp_out =size(results_out_out%mixalch,1) 499 ntypat_out=0;if (associated(results_out_out%mixalch))ntypat_out=size(results_out_out%mixalch,2) 500 501 if (nimage_>nimage_out) then 502 if (associated(results_out_out%acell)) then 503 ABI_DEALLOCATE(results_out_out%acell) 504 end if 505 if (associated(results_out_out%etotal)) then 506 ABI_DEALLOCATE(results_out_out%etotal) 507 end if 508 if (associated(results_out_out%rprim)) then 509 ABI_DEALLOCATE(results_out_out%rprim) 510 end if 511 if (associated(results_out_out%strten)) then 512 ABI_DEALLOCATE(results_out_out%strten) 513 end if 514 if (associated(results_out_out%vel_cell)) then 515 ABI_DEALLOCATE(results_out_out%vel_cell) 516 end if 517 ABI_ALLOCATE(results_out_out%acell,(3,nimage_)) 518 ABI_ALLOCATE(results_out_out%etotal,(nimage_)) 519 ABI_ALLOCATE(results_out_out%rprim,(3,3,nimage_)) 520 ABI_ALLOCATE(results_out_out%strten,(6,nimage_)) 521 ABI_ALLOCATE(results_out_out%vel_cell,(3,3,nimage_)) 522 end if 523 if (nimage_>nimage_out.or.natom_>natom_out) then 524 if (associated(results_out_out%fcart)) then 525 ABI_DEALLOCATE(results_out_out%fcart) 526 end if 527 if (associated(results_out_out%fred)) then 528 ABI_DEALLOCATE(results_out_out%fred) 529 end if 530 if (associated(results_out_out%vel)) then 531 ABI_DEALLOCATE(results_out_out%vel) 532 end if 533 if (associated(results_out_out%xred)) then 534 ABI_DEALLOCATE(results_out_out%xred) 535 end if 536 ABI_ALLOCATE(results_out_out%fcart,(3,natom_,nimage_)) 537 ABI_ALLOCATE(results_out_out%fred,(3,natom_,nimage_)) 538 ABI_ALLOCATE(results_out_out%vel,(3,natom_,nimage_)) 539 ABI_ALLOCATE(results_out_out%xred,(3,natom_,nimage_)) 540 end if 541 if (nimage_>nimage_out.or.nkpt_>nkpt_out) then 542 if (associated(results_out_out%npwtot)) then 543 ABI_DEALLOCATE(results_out_out%npwtot) 544 end if 545 ABI_ALLOCATE(results_out_out%npwtot,(nkpt_,nimage_)) 546 end if 547 if (nimage_>nimage_out.or.nocc_>nocc_out) then 548 if (associated(results_out_out%occ)) then 549 ABI_DEALLOCATE(results_out_out%occ) 550 end if 551 ABI_ALLOCATE(results_out_out%occ,(nocc_,nimage_)) 552 end if 553 if (ntypat_>ntypat_out) then 554 if (associated(results_out_out%amu)) then 555 ABI_DEALLOCATE(results_out_out%amu) 556 end if 557 ABI_ALLOCATE(results_out_out%amu,(ntypat_,nimage_)) 558 end if 559 560 if (npsp_>npsp_out.or.ntypat_>ntypat_out) then 561 if (associated(results_out_out%mixalch)) then 562 ABI_DEALLOCATE(results_out_out%mixalch) 563 end if 564 ABI_ALLOCATE(results_out_out%mixalch,(npsp_,ntypat_,nimage_)) 565 end if 566 567 results_out_out%nimage=results_out_in%nimage 568 results_out_out%natom =results_out_in%natom 569 results_out_out%nkpt =results_out_in%nkpt 570 results_out_out%nocc =results_out_in%nocc 571 results_out_out%acell(1:3,1:nimage_) =results_out_in%acell(1:3,1:nimage_) 572 results_out_out%amu(1:ntypat_,1:nimage_) =results_out_in%amu(1:ntypat_,1:nimage_) 573 results_out_out%etotal(1:nimage_) =results_out_in%etotal(1:nimage_) 574 results_out_out%fcart(1:3,1:natom_,1:nimage_)=results_out_in%fcart(1:3,1:natom_,1:nimage_) 575 results_out_out%fred(1:3,1:natom_,1:nimage_) =results_out_in%fred(1:3,1:natom_,1:nimage_) 576 results_out_out%mixalch(1:npsp_,1:ntypat_,1:nimage_)=results_out_in%mixalch(1:npsp_,1:ntypat_,1:nimage_) 577 results_out_out%npwtot(1:nkpt_,1:nimage_) =results_out_in%npwtot(1:nkpt_,1:nimage_) 578 results_out_out%occ(1:nocc_,1:nimage_) =results_out_in%occ(1:nocc_,1:nimage_) 579 results_out_out%rprim(1:3,1:3,1:nimage_) =results_out_in%rprim(1:3,1:3,1:nimage_) 580 results_out_out%strten(1:6,1:nimage_) =results_out_in%strten(1:6,1:nimage_) 581 results_out_out%xred(1:3,1:natom_,1:nimage_) =results_out_in%xred(1:3,1:natom_,1:nimage_) 582 results_out_out%vel(1:3,1:natom_,1:nimage_) =results_out_in%vel(1:3,1:natom_,1:nimage_) 583 results_out_out%vel_cell(1:3,1:3,1:nimage_) =results_out_in%vel_cell(1:3,1:3,1:nimage_) 584 585 end subroutine copy_results_out
m_results_out/destroy_results_out [ Functions ]
[ Top ] [ m_results_out ] [ Functions ]
NAME
destroy_results_out
FUNCTION
Clean and destroy an array of results_out datastructures
SIDE EFFECTS
results_out(:)=<type(results_out_type)>=results_out datastructure array
PARENTS
abinit,driver
CHILDREN
copy_results_out,init_results_out,xmpi_allgather,xmpi_allgatherv xmpi_gatherv
SOURCE
364 subroutine destroy_results_out(results_out) 365 366 367 !This section has been created automatically by the script Abilint (TD). 368 !Do not modify the following lines by hand. 369 #undef ABI_FUNC 370 #define ABI_FUNC 'destroy_results_out' 371 !End of the abilint section 372 373 implicit none 374 375 !Arguments ------------------------------------ 376 !arrays 377 type(results_out_type),intent(inout) :: results_out(:) 378 !Local variables------------------------------- 379 !scalars 380 integer :: idt1,idt2,ii,results_out_size 381 382 !************************************************************************ 383 384 !@results_out_type 385 386 results_out_size=size(results_out) 387 if (results_out_size>0) then 388 389 idt1=lbound(results_out,1);idt2=ubound(results_out,1) 390 do ii=idt1,idt2 391 results_out(ii)%nimage=0 392 results_out(ii)%natom=0 393 results_out(ii)%nkpt=0 394 results_out(ii)%nocc=0 395 if (associated(results_out(ii)%acell)) then 396 ABI_DEALLOCATE(results_out(ii)%acell) 397 end if 398 if (associated(results_out(ii)%amu)) then 399 ABI_DEALLOCATE(results_out(ii)%amu) 400 end if 401 if (associated(results_out(ii)%etotal)) then 402 ABI_DEALLOCATE(results_out(ii)%etotal) 403 end if 404 if (associated(results_out(ii)%fcart)) then 405 ABI_DEALLOCATE(results_out(ii)%fcart) 406 end if 407 if (associated(results_out(ii)%fred)) then 408 ABI_DEALLOCATE(results_out(ii)%fred) 409 end if 410 if (associated(results_out(ii)%mixalch)) then 411 ABI_DEALLOCATE(results_out(ii)%mixalch) 412 end if 413 if (associated(results_out(ii)%npwtot)) then 414 ABI_DEALLOCATE(results_out(ii)%npwtot) 415 end if 416 if (associated(results_out(ii)%occ)) then 417 ABI_DEALLOCATE(results_out(ii)%occ) 418 end if 419 if (associated(results_out(ii)%rprim)) then 420 ABI_DEALLOCATE(results_out(ii)%rprim) 421 end if 422 if (associated(results_out(ii)%strten)) then 423 ABI_DEALLOCATE(results_out(ii)%strten) 424 end if 425 if (associated(results_out(ii)%vel)) then 426 ABI_DEALLOCATE(results_out(ii)%vel) 427 end if 428 if (associated(results_out(ii)%vel_cell)) then 429 ABI_DEALLOCATE(results_out(ii)%vel_cell) 430 end if 431 if (associated(results_out(ii)%xred)) then 432 ABI_DEALLOCATE(results_out(ii)%xred) 433 end if 434 end do 435 436 end if 437 438 end subroutine destroy_results_out
m_results_out/gather_results_out [ Functions ]
[ Top ] [ m_results_out ] [ Functions ]
NAME
gather_results_out
FUNCTION
Gather results_out datastructure array using communicator over images (replicas) of the cell. Each contribution of single processor is gathered into a big array on master processor
INPUTS
allgather= --optional, default=false-- if TRUE do ALL_GATHER instead of GATHER dtsets(:)= <type datafiles_type> contains all input variables, master= --optional, default=0-- index of master proc receiving gathered data (if allgather=false) mpi_enregs=informations about MPI parallelization only_one_per_img= --optional, default=true-- if TRUE, the gather operation is only done by one proc per image (master of the comm_cell) results_out(:)=<type(results_out_type)>=results_out datastructure array on each proc use_results_all=true if results_out_all datastructure is allocated for current proc
SIDE EFFECTS
=== f use_results_all=true === results_out_all(:)=<type(results_out_type)>=global (gathered) results_out datastructure array
PARENTS
abinit,driver
CHILDREN
copy_results_out,init_results_out,xmpi_allgather,xmpi_allgatherv xmpi_gatherv
SOURCE
621 subroutine gather_results_out(dtsets,mpi_enregs,results_out,results_out_all,use_results_all,& 622 & master,allgather,only_one_per_img) ! optional arguments 623 624 625 !This section has been created automatically by the script Abilint (TD). 626 !Do not modify the following lines by hand. 627 #undef ABI_FUNC 628 #define ABI_FUNC 'gather_results_out' 629 !End of the abilint section 630 631 implicit none 632 633 !Arguments ------------------------------------ 634 !scalars 635 integer,optional,intent(in) :: master 636 logical,optional,intent(in) :: allgather,only_one_per_img 637 logical,intent(in) :: use_results_all 638 !arrays 639 type(dataset_type),intent(in) :: dtsets(:) 640 type(results_out_type),intent(in) :: results_out(:) 641 type(results_out_type),intent(inout) :: results_out_all(:) 642 type(MPI_type), intent(inout) :: mpi_enregs(:) 643 !Local variables------------------------------- 644 !scalars 645 integer :: dtsets_size 646 integer :: ibufi,ibufr 647 integer :: idt1,idt2,ierr,ii,iproc,jj 648 integer :: isize,isize_img 649 integer :: master_all,master_img,master_one_img 650 integer :: mpi_enregs_size,mxnatom,mxnband,mxnkpt,mxnpsp,mxnsppol,mxntypat 651 integer :: natom_,nkpt_,nocc_,npsp_,ntypat_,nimage,nimagetot 652 integer :: results_out_size,results_out_all_size 653 integer :: rsize,rsize_img 654 logical :: do_allgather,one_per_img 655 character(len=500) :: msg 656 ! type(MPI_type):: mpi_img 657 !arrays 658 integer,allocatable :: ibuffer(:),ibuffer_all(:),ibufshft(:) 659 integer,allocatable :: iimg(:),isize_img_all(:),nimage_all(:) 660 integer,allocatable :: rbufshft(:),rsize_img_all(:) 661 real(dp),allocatable :: rbuffer(:),rbuffer_all(:) 662 663 !************************************************************************ 664 665 !@results_out_type 666 667 one_per_img=.true.;if (present(only_one_per_img)) one_per_img=only_one_per_img 668 do_allgather=.false.;if (present(allgather)) do_allgather=allgather 669 master_all=0;if (present(master)) master_all=master 670 671 ! call init_mpi_enreg(mpi_img,init_mpi=.false.) 672 master_img=0;master_one_img=0 673 ! i_am_master=(mpi_img%me==master_all) 674 ! use_results_all= & 675 !& ((( do_allgather).and.( one_per_img).and.(mpi_img%me_cell==master_one_img)) .or. & 676 !& (( do_allgather).and.(.not.one_per_img)) .or. & 677 !& ((.not.do_allgather).and.( one_per_img).and.(mpi_img%me==master_all)) .or. & 678 !& ((.not.do_allgather).and.(.not.one_per_img).and.(mpi_img%me_img==master_img))) 679 680 dtsets_size=size(dtsets);results_out_size=size(results_out) 681 mpi_enregs_size=size(mpi_enregs) 682 if (dtsets_size/=results_out_size) then 683 msg=' Wrong sizes for dtsets and results_out datastructures !' 684 MSG_BUG(msg) 685 end if 686 if (mpi_enregs_size/=results_out_size) then 687 msg=' Wrong sizes for dtsets and results_out datastructures !' 688 MSG_BUG(msg) 689 end if 690 691 if (use_results_all) then 692 results_out_all_size=size(results_out_all) 693 if (results_out_size/=results_out_all_size) then 694 msg=' Wrong size for results_out_all datastructure !' 695 MSG_BUG(msg) 696 end if 697 end if 698 699 if (results_out_size>0) then 700 701 idt1=lbound(results_out,1);idt2=ubound(results_out,1) 702 703 ! Create global results_out_all datastructure 704 if (use_results_all) then 705 mxnatom=1;mxnband=1;mxnkpt=1;mxnpsp=1;mxntypat=1 706 do ii=idt1,idt2 707 isize=size(results_out(ii)%fcart,2) ;if (isize>mxnatom) mxnatom=isize 708 isize=size(results_out(ii)%occ,1) ;if (isize>mxnband) mxnband=isize 709 isize=size(results_out(ii)%mixalch,1);if(isize>mxnpsp) mxnpsp=isize 710 isize=size(results_out(ii)%npwtot,1);if (isize>mxnkpt) mxnkpt=isize 711 isize=size(results_out(ii)%mixalch,2);if(isize>mxntypat) mxntypat=isize 712 end do 713 mxnband=mxnband/mxnkpt;mxnsppol=1 714 call init_results_out(dtsets,2,0,mpi_enregs,mxnatom,mxnband,mxnkpt,mxnpsp,mxnsppol,mxntypat,results_out_all) 715 end if 716 717 ! Loop over results_out components (datasets) 718 do ii=idt1,idt2 719 720 ! Simple copy in case of 1 image 721 if (dtsets(ii)%npimage<=1) then 722 if (use_results_all) then 723 call copy_results_out(results_out(ii),results_out_all(ii)) 724 end if 725 else 726 727 ! Retrieve MPI informations for this dataset 728 729 if ((.not.one_per_img).or.(mpi_enregs(ii)%me_cell==master_one_img)) then 730 731 ! Gather number of images treated by each proc 732 ABI_ALLOCATE(nimage_all,(mpi_enregs(ii)%nproc_img)) 733 nimage_all=0 734 nimage=results_out(ii)%nimage 735 call xmpi_allgather(nimage,nimage_all,mpi_enregs(ii)%comm_img,ierr) 736 nimagetot=sum(nimage_all) 737 738 ! Copy scalars from distributed results_out to gathered one 739 if (use_results_all) then 740 results_out_all(ii)%nimage=nimagetot 741 results_out_all(ii)%natom =results_out(ii)%natom 742 results_out_all(ii)%nkpt =results_out(ii)%nkpt 743 results_out_all(ii)%nocc =results_out(ii)%nocc 744 results_out_all(ii)%npsp =results_out(ii)%npsp 745 results_out_all(ii)%ntypat=results_out(ii)%ntypat 746 end if 747 748 ! Compute number of integers/reals needed by current 749 ! results_out structure for current proc 750 isize=results_out(ii)%nkpt 751 rsize=28+12*results_out(ii)%natom+results_out(ii)%nocc+results_out(ii)%npsp*results_out(ii)%ntypat+results_out(ii)%ntypat 752 isize_img=results_out(ii)%nimage*isize 753 rsize_img=results_out(ii)%nimage*rsize 754 ABI_ALLOCATE(isize_img_all,(mpi_enregs(ii)%nproc_img)) 755 ABI_ALLOCATE(rsize_img_all,(mpi_enregs(ii)%nproc_img)) 756 isize_img_all(:)=isize*nimage_all(:) 757 rsize_img_all(:)=rsize*nimage_all(:) 758 ABI_DEALLOCATE(nimage_all) 759 760 ! Compute shifts in buffer arrays for each proc 761 ABI_ALLOCATE(ibufshft,(mpi_enregs(ii)%nproc_img)) 762 ibufshft(1)=0 763 ABI_ALLOCATE(rbufshft,(mpi_enregs(ii)%nproc_img)) 764 rbufshft(1)=0 765 do jj=2,mpi_enregs(ii)%nproc_img 766 ibufshft(jj)=ibufshft(jj-1)+isize_img_all(jj-1) 767 rbufshft(jj)=rbufshft(jj-1)+rsize_img_all(jj-1) 768 end do 769 770 ! Load buffers 771 ABI_ALLOCATE(ibuffer,(isize_img)) 772 ABI_ALLOCATE(rbuffer,(rsize_img)) 773 ibufi=0;ibufr=0 774 natom_=results_out(ii)%natom 775 nkpt_ =results_out(ii)%nkpt 776 nocc_ =results_out(ii)%nocc 777 npsp_ =results_out(ii)%npsp 778 ntypat_ =results_out(ii)%ntypat 779 do jj=1,results_out(ii)%nimage 780 ibuffer(ibufi+1:ibufi+nkpt_)=results_out(ii)%npwtot(1:nkpt_,jj) 781 ibufi=ibufi+nkpt_ 782 rbuffer(ibufr+1:ibufr+3)=results_out(ii)%acell(1:3,jj) 783 ibufr=ibufr+3 784 rbuffer(ibufr+1:ibufr+ntypat_)=results_out(ii)%amu(1:ntypat_,jj) 785 ibufr=ibufr+ntypat_ 786 rbuffer(ibufr+1)=results_out(ii)%etotal(jj) 787 ibufr=ibufr+1 788 rbuffer(ibufr+1:ibufr+3*natom_)=reshape(results_out(ii)%fcart(1:3,1:natom_,jj),(/3*natom_/)) 789 ibufr=ibufr+3*natom_ 790 rbuffer(ibufr+1:ibufr+3*natom_)=reshape(results_out(ii)%fred(1:3,1:natom_,jj),(/3*natom_/)) 791 ibufr=ibufr+3*natom_ 792 rbuffer(ibufr+1:ibufr+npsp_*ntypat_)=& 793 & reshape(results_out(ii)%mixalch(1:npsp_,1:ntypat_,jj),(/npsp_*ntypat_/) ) 794 ibufr=ibufr+npsp_*ntypat_ 795 rbuffer(ibufr+1:ibufr+nocc_)=results_out(ii)%occ(1:nocc_,jj) 796 ibufr=ibufr+nocc_ 797 rbuffer(ibufr+1:ibufr+9)=reshape(results_out(ii)%rprim(1:3,1:3,jj),(/9/)) 798 ibufr=ibufr+9 799 rbuffer(ibufr+1:ibufr+9)=reshape(results_out(ii)%vel_cell(1:3,1:3,jj),(/9/)) 800 ibufr=ibufr+9 801 rbuffer(ibufr+1:ibufr+6)=results_out(ii)%strten(1:6,jj) 802 ibufr=ibufr+6 803 rbuffer(ibufr+1:ibufr+3*natom_)=reshape(results_out(ii)%vel(1:3,1:natom_,jj),(/3*natom_/)) 804 ibufr=ibufr+3*natom_ 805 rbuffer(ibufr+1:ibufr+3*natom_)=reshape(results_out(ii)%xred(1:3,1:natom_,jj),(/3*natom_/)) 806 ibufr=ibufr+3*natom_ 807 end do 808 if (ibufi/=isize_img.or.ibufr/=rsize_img) then 809 msg=' wrong buffer sizes !' 810 MSG_BUG(msg) 811 end if 812 813 ! Gather all data 814 if (use_results_all) then 815 ABI_ALLOCATE(ibuffer_all,(isize*nimagetot)) 816 ABI_ALLOCATE(rbuffer_all,(rsize*nimagetot)) 817 end if 818 if (.not.use_results_all) then 819 ABI_ALLOCATE(ibuffer_all,(0)) 820 ABI_ALLOCATE(rbuffer_all,(0)) 821 end if 822 if (do_allgather) then 823 call xmpi_allgatherv(ibuffer,isize_img,ibuffer_all,isize_img_all,ibufshft,& 824 & mpi_enregs(ii)%comm_img,ierr) 825 call xmpi_allgatherv(rbuffer,rsize_img,rbuffer_all,rsize_img_all,rbufshft,& 826 & mpi_enregs(ii)%comm_img,ierr) 827 else 828 call xmpi_gatherv(ibuffer,isize_img,ibuffer_all,isize_img_all,ibufshft,& 829 & master_img,mpi_enregs(ii)%comm_img,ierr) 830 call xmpi_gatherv(rbuffer,rsize_img,rbuffer_all,rsize_img_all,rbufshft,& 831 & master_img,mpi_enregs(ii)%comm_img,ierr) 832 end if 833 ABI_DEALLOCATE(isize_img_all) 834 ABI_DEALLOCATE(rsize_img_all) 835 ABI_DEALLOCATE(ibuffer) 836 ABI_DEALLOCATE(rbuffer) 837 838 ! Transfer buffers into gathered results_out_all (master proc only) 839 if (use_results_all) then 840 ABI_ALLOCATE(iimg,(mpi_enregs(ii)%nproc_img)) 841 iimg=0 842 natom_=results_out_all(ii)%natom 843 nkpt_=results_out_all(ii)%nkpt 844 nocc_=results_out_all(ii)%nocc 845 npsp_ =results_out_all(ii)%npsp 846 ntypat_ =results_out_all(ii)%ntypat 847 do jj=1,nimagetot 848 ! The following line supposes that images are sorted by increasing index 849 iproc=mpi_enregs(ii)%distrb_img(jj)+1;iimg(iproc)=iimg(iproc)+1 850 ibufi=ibufshft(iproc)+(iimg(iproc)-1)*isize 851 ibufr=rbufshft(iproc)+(iimg(iproc)-1)*rsize 852 results_out_all(ii)%npwtot(1:nkpt_,jj)=ibuffer_all(ibufi+1:ibufi+nkpt_) 853 ibufi=ibufi+nkpt_ 854 results_out_all(ii)%acell(1:3,jj)=rbuffer_all(ibufr+1:ibufr+3) 855 ibufr=ibufr+3 856 results_out_all(ii)%amu(1:ntypat_,jj)=rbuffer_all(ibufr+1:ibufr+ntypat_) 857 ibufr=ibufr+ntypat_ 858 results_out_all(ii)%etotal(jj)=rbuffer_all(ibufr+1) 859 ibufr=ibufr+1 860 results_out_all(ii)%fcart(1:3,1:natom_,jj)= & 861 & reshape(rbuffer_all(ibufr+1:ibufr+3*natom_),(/3,natom_/)) 862 ibufr=ibufr+3*natom_ 863 results_out_all(ii)%fred(1:3,1:natom_,jj)= & 864 & reshape(rbuffer_all(ibufr+1:ibufr+3*natom_),(/3,natom_/)) 865 ibufr=ibufr+3*natom_ 866 results_out_all(ii)%mixalch(1:npsp_,1:ntypat_,jj)= & 867 & reshape(rbuffer_all(ibufr+1:ibufr+npsp_*ntypat_),(/npsp_,ntypat_/)) 868 ibufr=ibufr+npsp_*ntypat_ 869 results_out_all(ii)%occ(1:nocc_,jj)=rbuffer_all(ibufr+1:ibufr+nocc_) 870 ibufr=ibufr+nocc_ 871 results_out_all(ii)%rprim(1:3,1:3,jj)=reshape(rbuffer_all(ibufr+1:ibufr+9),(/3,3/)) 872 ibufr=ibufr+9 873 results_out_all(ii)%vel_cell(1:3,1:3,jj)=reshape(rbuffer_all(ibufr+1:ibufr+9),(/3,3/)) 874 ibufr=ibufr+9 875 results_out_all(ii)%strten(1:6,jj)=rbuffer_all(ibufr+1:ibufr+6) 876 ibufr=ibufr+6 877 results_out_all(ii)%vel(1:3,1:natom_,jj)= & 878 & reshape(rbuffer_all(ibufr+1:ibufr+3*natom_),(/3,natom_/)) 879 ibufr=ibufr+3*natom_ 880 results_out_all(ii)%xred(1:3,1:natom_,jj)= & 881 & reshape(rbuffer_all(ibufr+1:ibufr+3*natom_),(/3,natom_/)) 882 ibufr=ibufr+3*natom_ 883 end do 884 ABI_DEALLOCATE(iimg) 885 end if 886 887 ! Free memory 888 ABI_DEALLOCATE(ibufshft) 889 ABI_DEALLOCATE(rbufshft) 890 ABI_DEALLOCATE(ibuffer_all) 891 ABI_DEALLOCATE(rbuffer_all) 892 893 end if 894 end if 895 end do 896 end if 897 898 end subroutine gather_results_out
m_results_out/init_results_out [ Functions ]
[ Top ] [ m_results_out ] [ Functions ]
NAME
init_results_out
FUNCTION
Init all scalars and pointers in an array of results_out datastructures
INPUTS
dtsets(:)= <type datafiles_type> contains all input variables, option_alloc=0: only allocate datastructure 1: allocate and initialize the whole datastructure 2: allocate datastructure and initialize only first member option_size=0: allocate results_out with a global number images (use mxnimage=max(dtset%nimage)) 1: allocate results_out with a number of images per processor (use mxnimage=max(mpi_enreg%my_nimage)) mpi_enregs=information about MPI parallelization mxnimage=-optional- maximal value of nimage over datasets if this argument is present, it is used for allocations if it is not present, allocations are automatic natom= number of atoms nband= number of bands nkpt= number of k-points nsppol= number of independant spin components
SIDE EFFECTS
results_out(:)=<type(results_out_type)>=results_out datastructure array
PARENTS
abinit,m_results_out
CHILDREN
copy_results_out,init_results_out,xmpi_allgather,xmpi_allgatherv xmpi_gatherv
SOURCE
188 subroutine init_results_out(dtsets,option_alloc,option_size,mpi_enregs,& 189 & mxnatom,mxnband,mxnkpt,mxnpsp,mxnsppol,mxntypat,results_out) 190 191 192 !This section has been created automatically by the script Abilint (TD). 193 !Do not modify the following lines by hand. 194 #undef ABI_FUNC 195 #define ABI_FUNC 'init_results_out' 196 !End of the abilint section 197 198 implicit none 199 200 !Arguments ------------------------------------ 201 !scalars 202 integer,intent(in) :: option_alloc,option_size 203 integer,intent(in) :: mxnatom,mxnband,mxnkpt,mxnpsp,mxnsppol,mxntypat 204 !arrays 205 type(dataset_type),intent(in) :: dtsets(:) 206 type(results_out_type),intent(inout) :: results_out(:) 207 type(MPI_type), intent(in) :: mpi_enregs(:) 208 !Local variables------------------------------- 209 !scalars 210 integer :: dtsets_size,idt1,idt2,idt3,ii,jj,kk 211 integer :: mpi_enregs_size,mxnimage_,natom_,nkpt_,nocc_ 212 integer :: results_out_size 213 ! type(MPI_type) :: mpi_img 214 !arrays 215 integer,allocatable :: img(:,:),nimage(:) 216 real(dp),allocatable :: tmp(:,:) 217 218 !************************************************************************ 219 220 !@results_out_type 221 222 dtsets_size=size(dtsets) 223 results_out_size=size(results_out) 224 mpi_enregs_size=size(mpi_enregs) 225 if (dtsets_size/=mpi_enregs_size .or. dtsets_size/=results_out_size) then 226 MSG_ERROR("init_results_out: wrong sizes (2)!") 227 endif 228 229 if (results_out_size>0) then 230 231 idt1=lbound(results_out,1);idt2=ubound(results_out,1) 232 idt3=idt2;if (option_alloc==2) idt3=idt1 233 ABI_ALLOCATE(nimage,(idt1:idt2)) 234 nimage=0 235 mxnimage_=1 236 if (option_size==0) then 237 do ii=idt1,idt2 238 nimage(ii)=dtsets(ii)%nimage 239 if (nimage(ii)>mxnimage_) mxnimage_=nimage(ii) 240 end do 241 if (option_alloc>0) then 242 ABI_ALLOCATE(img,(mxnimage_,idt1:idt3)) 243 img=0 244 do ii=idt1,idt3 245 do jj=1,nimage(ii) 246 img(jj,ii)=jj 247 end do 248 end do 249 end if 250 else 251 do ii=idt1,idt2 252 nimage(ii)=mpi_enregs(ii)%my_nimage 253 if (nimage(ii)>mxnimage_) mxnimage_=nimage(ii) 254 end do 255 if (option_alloc>0) then 256 ABI_ALLOCATE(img,(mxnimage_,idt1:idt3)) 257 img=0 258 do ii=idt1,idt3 259 do jj=1,nimage(ii) 260 img(jj,ii)=mpi_enregs(ii)%my_imgtab(jj) 261 end do 262 end do 263 end if 264 end if 265 266 do ii=idt1,idt2 267 268 ABI_ALLOCATE(results_out(ii)%acell,(3,mxnimage_)) 269 ABI_ALLOCATE(results_out(ii)%amu,(mxntypat,mxnimage_)) 270 ABI_ALLOCATE(results_out(ii)%etotal,(mxnimage_)) 271 ABI_ALLOCATE(results_out(ii)%fcart,(3,mxnatom,mxnimage_)) 272 ABI_ALLOCATE(results_out(ii)%fred,(3,mxnatom,mxnimage_)) 273 ABI_ALLOCATE(results_out(ii)%mixalch,(mxnpsp,mxntypat,mxnimage_)) 274 ABI_ALLOCATE(results_out(ii)%npwtot,(mxnkpt,mxnimage_)) 275 ABI_ALLOCATE(results_out(ii)%occ,(mxnband*mxnkpt*mxnsppol,mxnimage_)) 276 ABI_ALLOCATE(results_out(ii)%rprim,(3,3,mxnimage_)) 277 ABI_ALLOCATE(results_out(ii)%strten,(6,mxnimage_)) 278 ABI_ALLOCATE(results_out(ii)%vel,(3,mxnatom,mxnimage_)) 279 ABI_ALLOCATE(results_out(ii)%vel_cell,(3,3,mxnimage_)) 280 ABI_ALLOCATE(results_out(ii)%xred,(3,mxnatom,mxnimage_)) 281 282 if ((option_alloc==1).or.(option_alloc==2.and.ii==idt3)) then 283 results_out(ii)%nimage=nimage(ii) 284 results_out(ii)%natom =mxnatom 285 results_out(ii)%nkpt =mxnkpt 286 results_out(ii)%npsp =mxnpsp 287 results_out(ii)%ntypat =mxntypat 288 results_out(ii)%nocc =mxnband*mxnkpt*mxnsppol 289 natom_=dtsets(ii)%natom 290 nkpt_=dtsets(ii)%nkpt;if(ii==0) nkpt_=mxnkpt 291 nocc_=mxnband*dtsets(ii)%nkpt*dtsets(ii)%nsppol 292 results_out(ii)%nimage=nimage(ii) 293 results_out(ii)%natom=natom_ 294 results_out(ii)%nkpt=nkpt_ 295 results_out(ii)%nocc=nocc_ 296 results_out(ii)%acell=zero 297 results_out(ii)%amu=zero 298 results_out(ii)%etotal(:)=zero 299 results_out(ii)%fred(:,:,:)=zero 300 results_out(ii)%fcart(:,:,:)=zero 301 results_out(ii)%mixalch(:,:,:)=zero 302 results_out(ii)%occ=zero 303 results_out(ii)%rprim=zero 304 results_out(ii)%strten(:,:)=zero 305 results_out(ii)%vel=zero 306 results_out(ii)%vel_cell=zero 307 results_out(ii)%xred=zero 308 results_out(ii)%npwtot(:,:)=0 309 if (nimage(ii)>0) then 310 do jj=1,nimage(ii) 311 kk=img(jj,ii) 312 results_out(ii)%acell(:,jj) =dtsets(ii)%acell_orig(:,kk) 313 results_out(ii)%amu(:,jj) =dtsets(ii)%amu_orig(:,kk) 314 results_out(ii)%rprim(:,:,jj) =dtsets(ii)%rprim_orig(:,:,kk) 315 results_out(ii)%vel_cell(:,:,jj)=dtsets(ii)%vel_cell_orig(:,:,kk) 316 results_out(ii)%mixalch(:,:,jj) =dtsets(ii)%mixalch_orig(:,:,kk) 317 if (natom_>0) then 318 ABI_ALLOCATE(tmp,(3,natom_)) 319 tmp(1:3,1:natom_)=dtsets(ii)%vel_orig(1:3,1:natom_,kk) 320 results_out(ii)%vel(1:3,1:natom_,jj)=tmp(1:3,1:natom_) 321 tmp(1:3,1:natom_)=dtsets(ii)%xred_orig(1:3,1:natom_,kk) 322 results_out(ii)%xred(1:3,1:natom_,jj)=tmp(1:3,1:natom_) 323 ABI_DEALLOCATE(tmp) 324 end if 325 if (nocc_>0) then 326 results_out(ii)%occ(1:nocc_,jj)=dtsets(ii)%occ_orig(1:nocc_) 327 end if 328 end do 329 end if 330 end if 331 332 end do 333 ABI_DEALLOCATE(nimage) 334 !if (option_size/=0.and.option_alloc==1) then 335 if (allocated(img)) then 336 ABI_DEALLOCATE(img) 337 end if 338 end if 339 340 end subroutine init_results_out
m_results_out/results_out_type [ Types ]
[ Top ] [ m_results_out ] [ Types ]
NAME
results_out_type
FUNCTION
This structured datatype contains a subset of the results of a GS calculation, needed to perform the so-called "internal tests", and to perform the timing analysis
SOURCE
62 type, public :: results_out_type 63 64 ! WARNING : if you modify this datatype, please check whether there might be creation/destruction/copy routines, 65 ! declared in another part of ABINIT, that might need to take into account your modification. 66 67 ! Integer scalar 68 69 integer :: natom 70 ! The number of atoms for this dataset 71 integer :: nimage 72 ! The number of images of the cell for this dataset (treated by current proc) 73 integer :: nkpt 74 ! The number of k-pints for this dataset 75 integer :: nocc 76 ! The number of occupations for this dataset 77 integer :: npsp 78 ! The number of pseudopotentials 79 integer :: ntypat 80 ! The number of types of atoms 81 82 ! Integer arrays 83 84 integer, pointer :: npwtot(:,:) 85 ! npw(mxnkpt,nimage) Full number of plane waves for each 86 ! k point, computed with the "true" rprimd 87 ! Not taking into account the decrease due to istwfk 88 ! Not taking into account the spread of pws on different procs 89 90 ! Real (real(dp)) arrays 91 92 real(dp), pointer :: acell(:,:) 93 ! acell(3,nimage) 94 ! Length of primitive vectors 95 96 real(dp), pointer :: amu(:,:) 97 ! amu(ntypat,nimage) 98 ! Mass of the atomic type 99 100 real(dp), pointer :: etotal(:) 101 ! etotal(nimage) 102 ! Total energy (Hartree) 103 104 real(dp), pointer :: fcart(:,:,:) 105 ! fcart(3,natom,nimage) Cartesian forces (Hartree/Bohr) 106 ! Forces in cartesian coordinates (Hartree) 107 108 real(dp), pointer :: fred(:,:,:) 109 ! fred(3,natom,nimage) 110 ! Forces in reduced coordinates (Hartree) 111 ! Actually, gradient of the total energy with respect 112 ! to change of reduced coordinates 113 114 real(dp), pointer :: mixalch(:,:,:) 115 ! mixalch(npsp,ntypat,nimage) [note that in psps datastructure, the dimensioning is npspalch,ntypalch] 116 ! Mixing coefficients going from the input pseudopotentials (those for alchemical mixing) to the alchemical atoms 117 118 real(dp), pointer :: occ(:,:) 119 ! occ(mxmband_upper*mxnkpt*mxnsppol,nimage) 120 ! Electronic occupations 121 122 real(dp), pointer :: rprim(:,:,:) 123 ! rprim(3,3,nimage) 124 ! Dimensionless real space primitive translations 125 126 real(dp), pointer :: strten(:,:) 127 ! strten(6,nimage) 128 ! Stress tensor 129 130 real(dp), pointer :: vel(:,:,:) 131 ! vel(3,natom,nimage) 132 ! Atomic velocities 133 134 real(dp), pointer :: vel_cell(:,:,:) 135 ! vel_cell(3,3,nimage) 136 ! Cell velocities 137 ! Time derivatives of dimensional primitive translations 138 139 real(dp), pointer :: xred(:,:,:) 140 ! xred(3,natom,nimage) 141 ! Atomic positions in reduced coordinates 142 143 end type results_out_type