TABLE OF CONTENTS
- ABINIT/m_xfpack
- ABINIT/xfh_recover_new
- ABINIT/xfh_update
- ABINIT/xfpack_f2vout
- ABINIT/xfpack_vin2x
- ABINIT/xfpack_x2vin
ABINIT/m_xfpack [ Modules ]
NAME
m_xfpack
FUNCTION
COPYRIGHT
Copyright (C) 1998-2022 ABINIT group (XG, MJV, DCA, GMR, JCC, SE) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
15 #if defined HAVE_CONFIG_H 16 #include "config.h" 17 #endif 18 19 #include "abi_common.h" 20 21 module m_xfpack 22 23 use defs_basis 24 use m_errors 25 use m_abicore 26 use m_abimover 27 28 use m_symtk, only : matr3inv 29 use m_geometry, only : mkradim, mkrdim, metric, strainsym 30 use m_results_gs , only : results_gs_type 31 use m_bfgs, only : hessupdt 32 33 implicit none 34 35 private
ABINIT/xfh_recover_new [ Functions ]
NAME
xfh_recover_new
FUNCTION
Update the contents of the history xfhist taking values from xred, acell, rprim, gred_corrected and strten
INPUTS
OUTPUT
SOURCE
637 subroutine xfh_recover_new(ab_xfh,ab_mover,acell,cycl_main,gred,& 638 & hessin,ndim,rprim,rprimd0,strten,ucvol,ucvol0,vin,vin_prev,vout,& 639 & vout_prev,xred) 640 641 !Arguments ------------------------------------ 642 !scalars 643 644 integer,intent(in) :: ndim 645 integer,intent(out) :: cycl_main 646 real(dp),intent(inout) :: ucvol,ucvol0 647 type(ab_xfh_type),intent(inout) :: ab_xfh 648 type(abimover),intent(in) :: ab_mover 649 650 651 !arrays 652 real(dp),intent(inout) :: acell(3) 653 real(dp),intent(inout) :: hessin(:,:) 654 real(dp),intent(inout) :: xred(3,ab_mover%natom) 655 real(dp),intent(inout) :: rprim(3,3) 656 real(dp),intent(inout) :: rprimd0(3,3) 657 real(dp),intent(inout) :: gred(3,ab_mover%natom) 658 real(dp),intent(inout) :: strten(6) 659 real(dp),intent(inout) :: vin(:) 660 real(dp),intent(inout) :: vin_prev(:) 661 real(dp),intent(inout) :: vout(:) 662 real(dp),intent(inout) :: vout_prev(:) 663 664 !Local variables------------------------------- 665 !scalars 666 integer :: ixfh ! kk,jj 667 668 !********************************************************************* 669 670 if(ab_xfh%nxfh/=0)then 671 ! Loop over previous time steps 672 do ixfh=1,ab_xfh%nxfh 673 674 ! For that time step, get new (x,f) from xfhist 675 xred(:,:) =ab_xfh%xfhist(:,1:ab_mover%natom ,1,ixfh) 676 rprim(1:3,1:3)=ab_xfh%xfhist(:,ab_mover%natom+2:ab_mover%natom+4,1,ixfh) 677 acell(:) =ab_xfh%xfhist(:,ab_mover%natom+1,1,ixfh) 678 gred(:,:) =ab_xfh%xfhist(:,1:ab_mover%natom,2,ixfh) 679 ! This use of results_gs is unusual 680 strten(1:3) =ab_xfh%xfhist(:,ab_mover%natom+2,2,ixfh) 681 strten(4:6) =ab_xfh%xfhist(:,ab_mover%natom+3,2,ixfh) 682 683 ! !DEBUG 684 ! write (ab_out,*) '---READED FROM XFHIST---' 685 686 ! write (ab_out,*) 'XRED' 687 ! do kk=1,ab_mover%natom 688 ! write (ab_out,*) xred(:,kk) 689 ! end do 690 ! write (ab_out,*) 'FRED' 691 ! do kk=1,ab_mover%natom 692 ! write (ab_out,*) gred(:,kk) 693 ! end do 694 ! write(ab_out,*) 'RPRIM' 695 ! do kk=1,3 696 ! write(ab_out,*) rprim(:,kk) 697 ! end do 698 ! write(ab_out,*) 'ACELL' 699 ! write(ab_out,*) acell(:) 700 ! !DEBUG 701 702 ! Transfer it in vin, vout 703 call xfpack_x2vin(acell,ab_mover%natom,& 704 & ndim,ab_mover%nsym,ab_mover%optcell,rprim,rprimd0,& 705 & ab_mover%symrel,ucvol,ucvol0,vin,xred) 706 call xfpack_f2vout(gred,ab_mover%natom,& 707 & ndim,ab_mover%optcell,ab_mover%strtarget,strten,& 708 & ucvol,vout) 709 ! Get old time step, if any, and update inverse hessian 710 if(ixfh/=1)then 711 xred(:,:) =ab_xfh%xfhist(:,1:ab_mover%natom,1,ixfh-1) 712 rprim(1:3,1:3)=& 713 & ab_xfh%xfhist(:,ab_mover%natom+2:ab_mover%natom+4,1,ixfh-1) 714 acell(:)=ab_xfh%xfhist(:,ab_mover%natom+1,1,ixfh-1) 715 gred(:,:)=ab_xfh%xfhist(:,1:ab_mover%natom,2,ixfh-1) 716 ! This use of results_gs is unusual 717 strten(1:3)=ab_xfh%xfhist(:,ab_mover%natom+2,2,ixfh-1) 718 strten(4:6)=ab_xfh%xfhist(:,ab_mover%natom+3,2,ixfh-1) 719 ! Tranfer it in vin_prev, vout_prev 720 call xfpack_x2vin(acell,ab_mover%natom,& 721 & ndim,ab_mover%nsym,ab_mover%optcell,rprim,rprimd0,& 722 & ab_mover%symrel,ucvol,ucvol0,vin_prev,xred) 723 call xfpack_f2vout(gred,ab_mover%natom,& 724 & ndim,ab_mover%optcell,ab_mover%strtarget,strten,& 725 & ucvol,vout_prev) 726 727 ! write(ab_out,*) 'Hessian matrix before update',ndim,'x',ndim 728 ! write(ab_out,*) 'ixfh=',ixfh 729 ! do kk=1,ndim 730 ! do jj=1,ndim,3 731 ! if (jj+2<=ndim)then 732 ! write(ab_out,*) jj,hessin(jj:jj+2,kk) 733 ! else 734 ! write(ab_out,*) jj,hessin(jj:ndim,kk) 735 ! end if 736 ! end do 737 ! end do 738 739 call hessupdt(hessin,ab_mover%iatfix,ab_mover%natom,ndim,& 740 & vin,vin_prev,vout,vout_prev) 741 742 ! !DEBUG 743 ! write(ab_out,*) 'Hessian matrix after update',ndim,'x',ndim 744 ! do kk=1,ndim 745 ! do jj=1,ndim,3 746 ! if (jj+2<=ndim)then 747 ! write(ab_out,*) jj,hessin(jj:jj+2,kk) 748 ! else 749 ! write(ab_out,*) jj,hessin(jj:ndim,kk) 750 ! end if 751 ! end do 752 ! end do 753 ! !DEBUG 754 755 end if !if(ab_xfh%nxfh/=0) 756 end do ! End loop over previous time steps 757 758 ! The hessian has been generated, 759 ! as well as the latest vin and vout 760 ! so will cycle the main loop 761 cycl_main=1 762 end if 763 764 end subroutine xfh_recover_new
ABINIT/xfh_update [ Functions ]
NAME
xfh_update
FUNCTION
Update the contents of the history xfhist taking values from xred, acell, rprim, gred_corrected and strten
INPUTS
OUTPUT
SOURCE
781 subroutine xfh_update(ab_xfh,acell,gred_corrected,natom,rprim,strten,xred) 782 783 !Arguments ------------------------------------ 784 !scalars 785 type(ab_xfh_type),intent(inout) :: ab_xfh 786 integer,intent(in) :: natom 787 788 !arrays 789 real(dp),intent(in) :: acell(3) 790 real(dp),intent(in) :: xred(3,natom) 791 real(dp),intent(in) :: rprim(3,3) 792 real(dp),intent(in) :: gred_corrected(3,natom) 793 real(dp),intent(in) :: strten(6) 794 795 !Local variables------------------------------- 796 !scalars 797 !integer :: kk 798 799 !********************************************************************* 800 801 !DEBUG 802 !write (ab_out,*) '---WROTE TO XFHIST---' 803 804 !write (ab_out,*) 'XRED' 805 !do kk=1,natom 806 !write (ab_out,*) xred(:,kk) 807 !end do 808 !write (ab_out,*) 'FRED' 809 !do kk=1,natom 810 !write (ab_out,*) gred_corrected(:,kk) 811 !end do 812 !write(ab_out,*) 'RPRIM' 813 !do kk=1,3 814 !write(ab_out,*) rprim(:,kk) 815 !end do 816 !write(ab_out,*) 'ACELL' 817 !write(ab_out,*) acell(:) 818 !DEBUG 819 820 ab_xfh%nxfh=ab_xfh%nxfh+1 821 822 ab_xfh%xfhist(:,1:natom,1,ab_xfh%nxfh)=xred(:,:) 823 ab_xfh%xfhist(:,natom+1,1,ab_xfh%nxfh)=acell(:) 824 ab_xfh%xfhist(:,natom+2:natom+4,1,ab_xfh%nxfh)=rprim(:,:) 825 ab_xfh%xfhist(:,1:natom,2,ab_xfh%nxfh)=gred_corrected(:,:) 826 ab_xfh%xfhist(:,natom+2,2,ab_xfh%nxfh)=strten(1:3) 827 ab_xfh%xfhist(:,natom+3,2,ab_xfh%nxfh)=strten(4:6) 828 829 end subroutine xfh_update
ABINIT/xfpack_f2vout [ Functions ]
NAME
xfpack_f2vout
FUNCTION
Old option=3, transfer gred and strten to vout
INPUTS
natom=number of atoms in cell ndim=dimension of vout arrays optcell=option for the optimisation of the unit cell. Described in abinit_help. Depending on its value, different part of strten are contained in vout. strtarget(6)=target stresses ; they will be subtracted from strten when vout is computed. ucvol=unit cell volume (bohr^3), needed for some values of optcell.
OUTPUT
(see side effects)
SIDE EFFECTS
Input/Output variables gred(3,natom)=grads of Etot wrt reduced coordinates (hartree) strten(6)=components of the stress tensor (hartree/bohr^3) vout(ndim)=vector that contains gred and some quantity derived from strten, depending on the value of optcell, and taking care ot strtarget
SOURCE
517 subroutine xfpack_f2vout(gred,natom,ndim,optcell,strtarget,strten,ucvol,vout) 518 519 !Arguments ------------------------------------ 520 !scalars 521 integer,intent(in) :: natom,ndim,optcell 522 real(dp),intent(in) :: ucvol 523 !arrays 524 real(dp),intent(in) :: strtarget(6) 525 real(dp),intent(in) :: gred(3,natom),strten(6) 526 real(dp),intent(out) :: vout(ndim) 527 528 !Local variables------------------------------- 529 !scalars 530 real(dp) :: strdiag 531 character(len=500) :: message 532 !arrays 533 real(dp) :: dstr(6) 534 535 ! ************************************************************************* 536 537 !!DEBUG 538 !write(ab_out,*) '' 539 !write(ab_out,*) 'xfpack_f2vout' 540 !write(ab_out,*) 'natom=',natom 541 !write(ab_out,*) 'ndim=',ndim 542 !write(ab_out,*) 'optcell=',optcell 543 !write(ab_out,*) 'ucvol=',ucvol 544 !!DEBUG 545 546 547 !########################################################## 548 !### 1. Test for compatible ndim 549 550 if(optcell==0 .and. ndim/=3*natom)then 551 write(message,'(a,a,a,i4,a,i4,a)' )& 552 & ' When optcell=0, ndim MUST be equal to 3*natom,',ch10,& 553 & ' while ndim=',ndim,' and 3*natom=',3*natom,'.' 554 ABI_BUG(message) 555 end if 556 557 if( optcell==1 .and. ndim/=3*natom+1)then 558 write(message,'(a,a,a,i4,a,i4,a)' )& 559 & ' When optcell=1, ndim MUST be equal to 3*natom+1,',ch10,& 560 & ' while ndim=',ndim,' and 3*natom+1=',3*natom+1,'.' 561 ABI_BUG(message) 562 end if 563 564 if( (optcell==2 .or. optcell==3) & 565 & .and. ndim/=3*natom+6)then 566 write(message,'(a,a,a,i4,a,i4,a)' )& 567 & ' When optcell=2 or 3, ndim MUST be equal to 3*natom+6,',ch10,& 568 & ' while ndim=',ndim,' and 3*natom+6=',3*natom+6,'.' 569 ABI_BUG(message) 570 end if 571 572 if( optcell>=4 .and. ndim/=3*natom+3)then 573 write(message,'(a,a,a,i4,a,i4,a)' )& 574 & ' When optcell=4,5,6,7,8 or 9, ndim MUST be equal to 3*natom+3,',ch10,& 575 & ' while ndim=',ndim,' and 3*natom+3=',3*natom+3,'.' 576 ABI_BUG(message) 577 end if 578 ! 579 !Get vout from gred and strten 580 ! 581 vout(1:3*natom)= reshape(gred(:,:), (/3*natom/) ) 582 dstr(:)=strten(:)-strtarget(:) 583 584 if(optcell==1)then 585 586 vout(3*natom+1)=( dstr(1)+dstr(2)+dstr(3))*ucvol 587 588 else if(optcell>=2)then 589 ! Eventually take away the trace 590 strdiag=0.0_dp 591 if(optcell==3) strdiag=(dstr(1)+dstr(2)+dstr(3))/3.0_dp 592 if(optcell==2 .or. optcell==3)then 593 vout(3*natom+1:3*natom+3)=(dstr(1:3)-strdiag)*ucvol 594 ! For non-diagonal derivatives, must take into account 595 ! that eps(i,j) AND eps(j,i) are varied at the same time. Thus, derivative 596 ! is twice larger 597 vout(3*natom+4:3*natom+6)=dstr(4:6)*ucvol*2.0_dp 598 else if(optcell==7 .or. optcell==8 .or. optcell==9)then 599 ! Similar to case optcell==2 or optcell==3, but in 2 dimensions. 600 vout(3*natom+1:3*natom+3)=dstr(1:3)*ucvol 601 vout(3*natom+optcell-6) =dstr(optcell-3)*ucvol*2.0_dp 602 else if (optcell==4)then 603 vout(3*natom+1) = dstr(1)*ucvol 604 vout(3*natom+2) = dstr(5)*ucvol 605 vout(3*natom+3) = dstr(6)*ucvol 606 else if (optcell==5)then 607 vout(3*natom+1) = dstr(2)*ucvol 608 vout(3*natom+2) = dstr(4)*ucvol 609 vout(3*natom+3) = dstr(6)*ucvol 610 else if (optcell==6)then 611 vout(3*natom+1) = dstr(3)*ucvol 612 vout(3*natom+2) = dstr(4)*ucvol 613 vout(3*natom+3) = dstr(5)*ucvol 614 end if 615 616 end if 617 618 end subroutine xfpack_f2vout
ABINIT/xfpack_vin2x [ Functions ]
NAME
xfpack_vin2x
FUNCTION
Old option=2, transfer vin to xred, acell and rprim
INPUTS
acell0(3)=reference length scales of primitive translations (bohr), needed for some values of optcell. natom=number of atoms in cell ndim=dimension of vin array nsym=order of group. rprimd0(3,3)=reference real space primitive translations, needed for some values of optcell. optcell=option for the optimisation of the unit cell. Described in abinit_help. Depending on its value, different part of acell and rprim are contained in vin. symrel(3,3,nsym)=symmetry operators in terms of action on primitive translations ucvol=unit cell volume (bohr^3), needed for some values of optcell. ucvol0=reference unit cell volume (bohr^3), needed for some values of optcell.
OUTPUT
(see side effects)
SIDE EFFECTS
Input/Output variables acell(3)=length scales of primitive translations (bohr) rprim(3,3)=dimensionless real space primitive translations vin(ndim)=vector that contains xred and some quantity derived from acell and rprim, depending on the value of optcell. xred(3,natom)=reduced dimensionless atomic coordinates
SOURCE
82 subroutine xfpack_vin2x(acell,acell0,natom,ndim,nsym,optcell,& 83 & rprim,rprimd0,symrel,ucvol,ucvol0,vin,xred) 84 85 !Arguments ------------------------------------ 86 !scalars 87 integer,intent(in) :: natom,ndim,nsym,optcell 88 real(dp),intent(in) :: ucvol0 89 real(dp),intent(out) :: ucvol 90 !arrays 91 integer,intent(in) :: symrel(3,3,nsym) 92 real(dp),intent(in) :: acell0(3),rprimd0(3,3) 93 real(dp),intent(inout) :: acell(3),rprim(3,3) 94 real(dp),intent(in) :: vin(ndim) 95 real(dp),intent(out) :: xred(3,natom) 96 97 !Local variables------------------------------- 98 !scalars 99 integer :: ii,jj,kk 100 real(dp) :: scale 101 character(len=500) :: message 102 logical :: equal=.TRUE. 103 !arrays 104 real(dp) :: gmet(3,3),gprimd(3,3),rmet(3,3),rprimd(3,3) 105 real(dp) :: rprimd_symm(3,3),scaling(3,3) 106 107 ! ************************************************************************* 108 109 !!DEBUG 110 !write(ab_out,*) '' 111 !write(ab_out,*) 'xfpack_vin2x' 112 !write(ab_out,*) 'natom=',natom 113 !write(ab_out,*) 'ndim=',ndim 114 !write(ab_out,*) 'nsym=',nsym 115 !write(ab_out,*) 'optcell=',optcell 116 !write(ab_out,*) 'ucvol=',ucvol 117 !write(ab_out,*) 'xred=' 118 !do kk=1,natom 119 !write(ab_out,*) xred(:,kk) 120 !end do 121 !write(ab_out,*) 'VECTOR INPUT (vin) xfpack_vin2x INPUT' 122 !do ii=1,ndim,3 123 !if (ii+2<=ndim)then 124 !write(ab_out,*) ii,vin(ii:ii+2) 125 !else 126 !write(ab_out,*) ii,vin(ii:ndim) 127 !end if 128 !end do 129 !!DEBUG 130 131 132 !########################################################## 133 !### 1. Test for compatible ndim 134 135 if(optcell==0 .and. ndim/=3*natom)then 136 write(message,'(a,a,a,i4,a,i4,a)' )& 137 & ' When optcell=0, ndim MUST be equal to 3*natom,',ch10,& 138 & ' while ndim=',ndim,' and 3*natom=',3*natom,'.' 139 ABI_BUG(messagE) 140 end if 141 142 if( (optcell==1) .and. ndim/=3*natom+1)then 143 write(message,'(a,a,a,i4,a,i4,a)' )& 144 & ' When optcell=1 ndim MUST be equal to 3*natom+1,',ch10,& 145 & ' while ndim=',ndim,' and 3*natom+1=',3*natom+1,'.' 146 ABI_BUG(message) 147 end if 148 149 if( (optcell==2 .or. optcell==3) & 150 & .and. ndim/=3*natom+6)then 151 write(message,'(a,a,a,i4,a,i4,a)' )& 152 & ' When optcell=2 or 3, ndim MUST be equal to 3*natom+6,',ch10,& 153 & ' while ndim=',ndim,' and 3*natom+6=',3*natom+6,'.' 154 ABI_BUG(message) 155 end if 156 157 if( optcell>=4 .and. ndim/=3*natom+3)then 158 write(message,'(a,a,a,i4,a,i4,a)' )& 159 & ' When optcell=4,5,6,7,8 or 9, ndim MUST be equal to 3*natom+3,',ch10,& 160 & ' while ndim=',ndim,' and 3*natom+3=',3*natom+3,'.' 161 ABI_BUG(message) 162 end if 163 164 !########################################################## 165 !### 3. option=2, transfer vin to xred, acell and rprim 166 167 !Get xred, and eventually acell and rprim from vin 168 xred(:,:)=reshape( vin(1:3*natom), (/3,natom/) ) 169 170 if(optcell==1)then 171 172 ! acell(:)=acell0(:)*vin(3*natom+1)/(ucvol0**third) 173 acell(:)=acell0(:)*vin(3*natom+1) 174 175 else if (optcell>=2)then 176 177 scaling(:,:)=0.0_dp 178 scaling(1,1)=1.0_dp ; scaling(2,2)=1.0_dp ; scaling(3,3)=1.0_dp 179 180 if(optcell==2 .or. optcell==3)then 181 scaling(1,1)=vin(3*natom+1) 182 scaling(2,2)=vin(3*natom+2) 183 scaling(3,3)=vin(3*natom+3) 184 scaling(2,3)=vin(3*natom+4) ; scaling(3,2)=vin(3*natom+4) 185 scaling(1,3)=vin(3*natom+5) ; scaling(3,1)=vin(3*natom+5) 186 scaling(1,2)=vin(3*natom+6) ; scaling(2,1)=vin(3*natom+6) 187 else if(optcell==4)then 188 scaling(1,1)=vin(3*natom+1) 189 if (abs(scaling(1,1) - 1.0_dp) <1.E-14) scaling(1,1)=1.0_dp 190 scaling(3,1)=vin(3*natom+2) 191 if (abs(scaling(3,1)) <1.E-14) scaling(3,1)=0.0_dp 192 scaling(2,1)=vin(3*natom+3) 193 if (abs(scaling(2,1)) <1.E-14) scaling(2,1)=0.0_dp 194 else if(optcell==5)then 195 scaling(2,2)=vin(3*natom+1) 196 if (abs(scaling(2,2) - 1.0_dp) <1.E-14) scaling(2,2)=1.0_dp 197 scaling(3,2)=vin(3*natom+2) 198 if (abs(scaling(3,2)) <1.E-14) scaling(3,2)=0.0_dp 199 scaling(1,2)=vin(3*natom+3) 200 if (abs(scaling(1,2)) <1.E-14) scaling(1,2)=0.0_dp 201 else if(optcell==6)then 202 scaling(3,3)=vin(3*natom+1) 203 if (abs(scaling(3,3) - 1.0_dp) <1.E-14) scaling(3,3)=1.0_dp 204 scaling(2,3)=vin(3*natom+2) 205 if (abs(scaling(2,3)) <1.E-14) scaling(2,3)=0.0_dp 206 scaling(1,3)=vin(3*natom+3) 207 if (abs(scaling(1,3)) <1.E-14) scaling(1,3)=0.0_dp 208 else if(optcell==7)then 209 scaling(2,2)=vin(3*natom+2) ; scaling(3,3)=vin(3*natom+3) 210 scaling(2,3)=vin(3*natom+1) ; scaling(3,2)=vin(3*natom+1) 211 else if(optcell==8)then 212 scaling(1,1)=vin(3*natom+1) ; scaling(3,3)=vin(3*natom+3) 213 scaling(1,3)=vin(3*natom+2) ; scaling(3,1)=vin(3*natom+2) 214 else if(optcell==9)then 215 scaling(1,1)=vin(3*natom+1) ; scaling(2,2)=vin(3*natom+2) 216 scaling(1,2)=vin(3*natom+3) ; scaling(2,1)=vin(3*natom+3) 217 end if 218 if(optcell<=3 .or. optcell>=7)then 219 do ii=1,3 220 do jj=1,3 221 rprimd(ii,jj)=0.0_dp 222 do kk=1,3 223 rprimd(ii,jj)=rprimd(ii,jj)+scaling(ii,kk)*rprimd0(kk,jj) 224 end do 225 end do 226 end do 227 ! for optcell=4,5,6, implementing search for all 3 components of the vector to be relaxed according to Eq.10 of J. Chem. Phys. 228 ! 136, 074103 (2012), i.e. search direction given by rprimd0 * stress 229 else if(optcell==4)then 230 rprimd(:,2) = rprimd0(:,2) 231 rprimd(:,3) = rprimd0(:,3) 232 rprimd(:,1) = 0.0_dp 233 do ii=1,3 234 do kk=1,3 235 rprimd(ii,1) = rprimd(ii,1) + scaling(kk,1)*rprimd0(ii,kk) 236 end do 237 end do 238 else if(optcell==5)then 239 rprimd(:,1) = rprimd0(:,1) 240 rprimd(:,3) = rprimd0(:,3) 241 rprimd(:,2) = 0.0_dp 242 do ii=1,3 243 do kk=1,3 244 rprimd(ii,2) = rprimd(ii,2) + scaling(kk,2)*rprimd0(ii,kk) 245 end do 246 end do 247 else if(optcell==6)then 248 rprimd(:,1) = rprimd0(:,1) 249 rprimd(:,2) = rprimd0(:,2) 250 rprimd(:,3) = 0.0_dp 251 do ii=1,3 252 do kk=1,3 253 rprimd(ii,3) = rprimd(ii,3) + scaling(kk,3)*rprimd0(ii,kk) 254 end do 255 end do 256 end if 257 258 ! Rescale if the volume must be preserved 259 if(optcell==3)then 260 call metric(gmet,gprimd,-1,rmet,rprimd,ucvol) 261 scale=(ucvol0/ucvol)**third 262 rprimd(:,:)=scale*rprimd(:,:) 263 end if 264 call strainsym(nsym,rprimd0,rprimd,rprimd_symm,symrel) 265 do jj=1,3 266 do ii=1,3 267 ! write(ab_out,*) 'DIFF',ii,jj,abs(rprimd0(ii,jj)-rprimd_symm(ii,jj)) 268 if (abs(rprimd0(ii,jj)-rprimd_symm(ii,jj))>1.E-14)& 269 & equal=.FALSE. 270 end do 271 end do 272 273 if (equal)then 274 acell(:)=acell0(:) 275 rprimd(:,:)=rprimd0(:,:) 276 else 277 ! Use a representation based on normalised rprim vectors 278 call mkradim(acell,rprim,rprimd_symm) 279 end if 280 281 end if 282 283 end subroutine xfpack_vin2x
ABINIT/xfpack_x2vin [ Functions ]
NAME
xfpack_x2vin
FUNCTION
Old option=1, transfer xred, acell, and rprim to vin
INPUTS
acell0(3)=reference length scales of primitive translations (bohr), needed for some values of optcell. natom=number of atoms in cell ndim=dimension of vin arrays nsym=order of group. rprimd0(3,3)=reference real space primitive translations, needed for some values of optcell. optcell=option for the optimisation of the unit cell. Described in abinit_help. Depending on its value, different part of acell and rprim are contained in vin. symrel(3,3,nsym)=symmetry operators in terms of action on primitive translations ucvol=unit cell volume (bohr^3), needed for some values of optcell. ucvol0=reference unit cell volume (bohr^3), needed for some values of optcell.
OUTPUT
(see side effects)
SIDE EFFECTS
Input/Output variables acell(3)=length scales of primitive translations (bohr) rprim(3,3)=dimensionless real space primitive translations vin(ndim)=vector that contains xred and some quantity derived from acell and rprim, depending on the value of optcell. xred(3,natom)=reduced dimensionless atomic coordinates
SOURCE
320 subroutine xfpack_x2vin(acell,natom,ndim,nsym,optcell,& 321 & rprim,rprimd0,symrel,ucvol,ucvol0,vin,xred) 322 323 !Arguments ------------------------------------ 324 !scalars 325 integer,intent(in) :: natom,ndim,nsym,optcell 326 real(dp),intent(in) :: ucvol0 327 real(dp),intent(inout) :: ucvol !vz_i 328 !arrays 329 integer,intent(in) :: symrel(3,3,nsym) 330 real(dp),intent(in) :: rprimd0(3,3) 331 real(dp),intent(in) :: acell(3),rprim(3,3) 332 real(dp),intent(in) :: xred(3,natom) 333 real(dp),intent(out) :: vin(ndim) 334 335 !Local variables------------------------------- 336 !scalars 337 integer :: ii,jj,kk 338 real(dp) :: scale 339 character(len=500) :: message 340 !arrays 341 real(dp) :: gmet(3,3),gprimd(3,3),gprimd0(3,3),rmet(3,3),rprimd(3,3) 342 real(dp) :: rprimd_symm(3,3),scaling(3,3) 343 344 ! ************************************************************************* 345 346 !!DEBUG 347 !write(ab_out,*) '' 348 !write(ab_out,*) 'xfpack_x2vin' 349 !write(ab_out,*) 'natom=',natom 350 !write(ab_out,*) 'ndim=',ndim 351 !write(ab_out,*) 'nsym=',nsym 352 !write(ab_out,*) 'optcell=',optcell 353 !write(ab_out,*) 'ucvol=',ucvol 354 !write(ab_out,*) 'xred=' 355 !do kk=1,natom 356 !write(ab_out,*) xred(:,kk) 357 !end do 358 !write(ab_out,*) 'VECTOR INPUT (vin) xfpack_x2vin INPUT' 359 !do ii=1,ndim,3 360 !if (ii+2<=ndim)then 361 !write(ab_out,*) ii,vin(ii:ii+2) 362 !else 363 !write(ab_out,*) ii,vin(ii:ndim) 364 !end if 365 !end do 366 !!DEBUG 367 368 369 !########################################################## 370 !### 1. Test for compatible ndim 371 372 if(optcell==0 .and. ndim/=3*natom)then 373 write(message,'(a,a,a,i4,a,i4,a)' )& 374 & ' When optcell=0, ndim MUST be equal to 3*natom,',ch10,& 375 & ' while ndim=',ndim,' and 3*natom=',3*natom,'.' 376 ABI_BUG(message) 377 end if 378 379 if( optcell==1 .and. ndim/=3*natom+1)then 380 write(message,'(a,a,a,i4,a,i4,a)' )& 381 & ' When optcell=1, ndim MUST be equal to 3*natom+1,',ch10,& 382 & ' while ndim=',ndim,' and 3*natom+1=',3*natom+1,'.' 383 ABI_BUG(message) 384 end if 385 386 if( (optcell==2 .or. optcell==3) & 387 & .and. ndim/=3*natom+6)then 388 write(message,'(a,a,a,i4,a,i4,a)' )& 389 & ' When optcell=2,3,4,5,6, ndim MUST be equal to 3*natom+6,',ch10,& 390 & ' while ndim=',ndim,' and 3*natom+6=',3*natom+6,'.' 391 ABI_BUG(message) 392 end if 393 394 if( optcell>=4 .and. ndim/=3*natom+3)then 395 write(message,'(a,a,a,i4,a,i4,a)' )& 396 & ' When optcell=4,5,6,7,8 or 9, ndim MUST be equal to 3*natom+3,',ch10,& 397 & ' while ndim=',ndim,' and 3*natom+3=',3*natom+3,'.' 398 ABI_BUG(message) 399 end if 400 401 !########################################################## 402 !### 2. option=1, transfer xred, acell, and rprim to vin 403 404 !Get vin from xred, acell, and rprim 405 vin(1:3*natom)= reshape(xred(:,:), (/3*natom/) ) 406 407 if(optcell/=0)then 408 call mkrdim(acell,rprim,rprimd) 409 call strainsym(nsym,rprimd0,rprimd,rprimd_symm,symrel) 410 call metric(gmet,gprimd,-1,rmet,rprimd_symm,ucvol) 411 412 if(optcell==1)then 413 414 ! vin(3*natom+1)=ucvol**third 415 vin(3*natom+1)=(ucvol/ucvol0)**third 416 417 else if(optcell>=2)then 418 419 ! Generates gprimd0 420 call matr3inv(rprimd0,gprimd0) 421 if (optcell==2 .or. optcell==3 .or. optcell>=7)then 422 do ii=1,3 423 do jj=1,3 424 scaling(ii,jj)=0.0_dp 425 do kk=1,3 426 scaling(ii,jj)=scaling(ii,jj)+rprimd_symm(ii,kk)*gprimd0(jj,kk) 427 end do 428 end do 429 end do 430 ! Rescale if the volume must be preserved 431 if(optcell==3)then 432 scale=(ucvol0/ucvol)**third 433 scaling(:,:)=scale*scaling(:,:) 434 end if 435 if(optcell==2 .or. optcell==3)then 436 vin(3*natom+1)=scaling(1,1) ; vin(3*natom+4)=(scaling(2,3)+scaling(3,2))*0.5_dp 437 vin(3*natom+2)=scaling(2,2) ; vin(3*natom+5)=(scaling(1,3)+scaling(3,1))*0.5_dp 438 vin(3*natom+3)=scaling(3,3) ; vin(3*natom+6)=(scaling(1,2)+scaling(2,1))*0.5_dp 439 else if(optcell>=7)then 440 vin(3*natom+1)=scaling(1,1) 441 vin(3*natom+2)=scaling(2,2) 442 vin(3*natom+3)=scaling(3,3) 443 if(optcell==7)vin(3*natom+1)=(scaling(2,3)+scaling(3,2))*0.5_dp 444 if(optcell==8)vin(3*natom+2)=(scaling(1,3)+scaling(3,1))*0.5_dp 445 if(optcell==9)vin(3*natom+3)=(scaling(1,2)+scaling(2,1))*0.5_dp 446 end if 447 end if 448 449 if (optcell==4)then 450 scaling(:,:) = 0.0_dp 451 do ii=1,3 452 do kk=1,3 453 scaling(ii,1) = scaling(ii,1) + gprimd0(kk,ii)*rprimd_symm(kk,1) 454 end do 455 end do 456 vin(3*natom+1) = scaling(1,1) 457 vin(3*natom+2) = scaling(3,1) 458 vin(3*natom+3) = scaling(2,1) 459 else if (optcell==5)then 460 scaling(:,:) = 0.0_dp 461 do ii=1,3 462 do kk=1,3 463 scaling(ii,2) = scaling(ii,2) + gprimd0(kk,ii)*rprimd_symm(kk,2) 464 end do 465 end do 466 vin(3*natom+1) = scaling(2,2) 467 vin(3*natom+2) = scaling(3,2) 468 vin(3*natom+3) = scaling(1,2) 469 else if (optcell==6)then 470 scaling(:,:) = 0.0_dp 471 do ii=1,3 472 do kk=1,3 473 scaling(ii,3) = scaling(ii,3) + gprimd0(kk,ii)*rprimd_symm(kk,3) 474 end do 475 end do 476 vin(3*natom+1) = scaling(3,3) 477 vin(3*natom+2) = scaling(2,3) 478 vin(3*natom+3) = scaling(1,3) 479 end if 480 481 end if 482 483 end if 484 485 end subroutine xfpack_x2vin