TABLE OF CONTENTS


ABINIT/m_xfpack [ Modules ]

[ Top ] [ 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 ]

[ Top ] [ 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 ]

[ Top ] [ 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 ]

[ Top ] [ 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 ]

[ Top ] [ 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 ]

[ Top ] [ 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