TABLE OF CONTENTS

ABINIT/m_crystal [ Modules ]

[ Top ] [ Modules ]

NAME

 m_crystal

FUNCTION

 Module containing the definition of the crystal_t data type and methods used to handle it.

COPYRIGHT

  Copyright (C) 2008-2024 ABINIT group (MG, YP, MJV, GA)
  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

16 #if defined HAVE_CONFIG_H
17 #include "config.h"
18 #endif
19 
20 #include "abi_common.h"
21 
22 module m_crystal
23 
24  use defs_basis
25  use m_errors
26  use m_abicore
27  use m_atomdata
28  use m_xmpi
29  use m_nctk
30  use, intrinsic :: iso_c_binding
31 #ifdef HAVE_NETCDF
32  use netcdf
33 #endif
34 
35  use m_io_tools,       only : file_exists
36  use m_numeric_tools,  only : set2unit
37  use m_hide_lapack,    only : matrginv
38  use m_fstrings,       only : int2char10, sjoin, yesno, itoa, strcat
39  use m_symtk,          only : mati3inv, sg_multable, symatm, print_symmetries
40  use m_spgdata,        only : spgdata
41  use m_geometry,       only : metric, xred2xcart, xcart2xred, remove_inversion, getspinrot, symredcart, normv
42  use m_io_tools,       only : open_file
43 
44  implicit none
45 
46  private

m_crystal/adata_type [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

  adata_type

FUNCTION

  Return atomic data from the itypat index

SOURCE

1275 type(atomdata_t) function adata_type(crystal, itypat) result(atom)
1276 
1277 !Arguments ------------------------------------
1278 !scalars
1279  integer,intent(in) :: itypat
1280  class(crystal_t),intent(in) :: crystal
1281 
1282 ! *************************************************************************
1283 
1284  call atomdata_from_znucl(atom, crystal%znucl(itypat))
1285 
1286 end function adata_type

m_crystal/crystal_bcast [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

  crystal_bcast

FUNCTION

  Master broadcasts data and others allocate their arrays.

SOURCE

791 subroutine crystal_bcast(Cryst, comm)
792 
793 !Arguments ------------------------------------
794  class(crystal_t),intent(inout) :: Cryst
795  integer, intent(in) :: comm
796 
797 !Local variables -------------------------
798  integer, parameter :: master=0
799  integer :: ierr
800 
801 ! *********************************************************************
802 
803  if (xmpi_comm_size(comm) == 1) return
804 
805  DBG_ENTER("COLL")
806 
807  ! Integers
808  call xmpi_bcast(Cryst%natom, master, comm, ierr)
809  call xmpi_bcast(Cryst%nsym, master, comm, ierr)
810  call xmpi_bcast(Cryst%ntypat, master, comm, ierr)
811  call xmpi_bcast(Cryst%nirredat, master, comm, ierr)
812  call xmpi_bcast(Cryst%npsp, master, comm, ierr)
813  call xmpi_bcast(Cryst%space_group, master, comm, ierr)
814  call xmpi_bcast(Cryst%timrev, master, comm, ierr)
815  call xmpi_bcast(Cryst%use_antiferro, master, comm, ierr)
816 
817  if (xmpi_comm_rank(comm) /= master) then
818    call Cryst%free()
819    call Cryst%malloc()
820  end if
821 
822  ! Floats
823  call xmpi_bcast(Cryst%ucvol, master, comm, ierr)
824 
825  ! Arrays
826  call xmpi_bcast(Cryst%angdeg, master, comm, ierr)
827  call xmpi_bcast(Cryst%gmet, master, comm, ierr)
828  call xmpi_bcast(Cryst%gprimd, master, comm, ierr)
829  call xmpi_bcast(Cryst%rmet, master, comm, ierr)
830  call xmpi_bcast(Cryst%rprimd, master, comm, ierr)
831  call xmpi_bcast(Cryst%indsym, master, comm, ierr)
832  call xmpi_bcast(Cryst%symafm, master, comm, ierr)
833  call xmpi_bcast(Cryst%symrec, master, comm, ierr)
834  call xmpi_bcast(Cryst%symrel, master, comm, ierr)
835  call xmpi_bcast(Cryst%symrel_cart, master, comm, ierr)
836  call xmpi_bcast(Cryst%atindx, master, comm, ierr)
837  call xmpi_bcast(Cryst%atindx1, master, comm, ierr)
838  call xmpi_bcast(Cryst%typat, master, comm, ierr)
839  call xmpi_bcast(Cryst%nattyp, master, comm, ierr)
840  call xmpi_bcast(Cryst%tnons, master, comm, ierr)
841  call xmpi_bcast(Cryst%xcart, master, comm, ierr)
842  call xmpi_bcast(Cryst%xred, master, comm, ierr)
843  call xmpi_bcast(Cryst%spinrot, master, comm, ierr)
844  call xmpi_bcast(Cryst%amu, master, comm, ierr)
845  call xmpi_bcast(Cryst%zion, master, comm, ierr)
846  call xmpi_bcast(Cryst%znucl, master, comm, ierr)
847 
848 
849  call xmpi_bcast(Cryst%title, master, comm, ierr)
850 
851  ! It is not always allocated on master node,
852  ! and it can be computed afterward on each node.
853  !call xmpi_bcast(Cryst%irredatindx, master, comm, ierr)
854 
855  DBG_EXIT("COLL")
856 
857 end subroutine crystal_bcast

m_crystal/crystal_compare [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

  crystal_compare

FUNCTION

   Compare two crystalline structures,
   write warning messages to stdout if they differ, return exit status

INPUTS

  [header]=Optional header message.

OUTPUT

SOURCE

877 integer function crystal_compare(self, other, header) result(ierr)
878 
879 !Arguments ------------------------------------
880 !scalars
881  class(crystal_t),intent(in) :: self, other
882  character(len=*),optional,intent(in) :: header
883 
884 !Local variables-------------------------------
885  !integer :: isym, iat, itypat
886 ! *********************************************************************
887 
888  if (present(header)) call wrtout(std_out, header)
889  ierr = 0
890 
891  ! Test basic dimensions and metadata.
892  ABI_CHECK_IEQ_IERR(self%natom, other%natom, "Different natom" , ierr)
893  ABI_CHECK_IEQ_IERR(self%ntypat, other%ntypat, "Different ntypat" , ierr)
894  ABI_CHECK_IEQ_IERR(self%npsp, other%npsp, "Different npsp" , ierr)
895  ABI_CHECK_IEQ_IERR(self%nsym, other%nsym, "Different nsym" , ierr)
896  ABI_CHECK_IEQ_IERR(self%timrev, other%timrev, "Different timrev" , ierr)
897 
898  if (ierr /= 0) goto 10
899  ! After this point, we know that basic dimensions agree with each other.
900  ! Yes, I use GOTO and I'm proud of that!
901 
902  ! Check direct lattice
903  if (any(abs(self%rprimd - other%rprimd) > tol6)) then
904    ABI_WARNING("Found critical diffs in rprimd lattice vectors.")
905    ierr = ierr + 1
906  end if
907 
908  ! Check Symmetries
909  if (any(self%symrel /= other%symrel)) then
910    ABI_WARNING("Found critical diffs in symrel symmetries.")
911    ierr = ierr + 1
912  end if
913  if (any(abs(self%tnons - other%tnons) > tol3)) then
914    ABI_WARNING("Found critical diffs in fractional translations tnons.")
915    ierr = ierr + 1
916  end if
917  if (self%use_antiferro .neqv. other%use_antiferro) then
918    ABI_WARNING("Different values of use_antiferro")
919    ierr = ierr + 1
920  end if
921 
922  ! Atoms
923  if (any(self%typat /= other%typat)) then
924    ABI_WARNING("Found critical diffs in typat.")
925    ierr = ierr + 1
926  end if
927  if (any(abs(self%zion - other%zion) > tol3)) then
928    ABI_WARNING("Found critical diffs in zion.")
929    ierr = ierr + 1
930  end if
931  if (any(abs(self%znucl - other%znucl) > tol3)) then
932    ABI_WARNING("Found critical diffs in znucl.")
933    ierr = ierr + 1
934  end if
935  if (any(abs(self%amu - other%amu) > tol3)) then
936    ABI_WARNING("Found critical diffs in amu.")
937    ierr = ierr + 1
938  end if
939  if (any(abs(self%xred - other%xred) > tol6)) then
940    ABI_WARNING("Found critical diffs in xred.")
941    ierr = ierr + 1
942  end if
943 
944  if (ierr /= 0) goto 10
945  return
946 
947  ! Print structure to aid debugging. Caller will handle exit status.
948 10 call wrtout(std_out, " Comparing crystal1 and crystal2 for possible differences before returning ierr /= 0!")
949    call self%print(header="crystal1")
950    call wrtout(std_out, "")
951    call other%print(header="crystal2")
952    call wrtout(std_out, "")
953 
954 end function crystal_compare

m_crystal/crystal_compute_geometry [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

  crystal_compute_geometry

FUNCTION

  Compute the different metrics and the angle between primitive vectors.
  Also compute cartesian coordinates of atoms.

SOURCE

578 subroutine crystal_compute_geometry(Cryst)
579 
580 !Arguments ------------------------------------
581  class(crystal_t),intent(inout) :: Cryst
582 
583 ! *********************************************************************
584 
585  call metric(Cryst%gmet,Cryst%gprimd,-1,Cryst%rmet,Cryst%rprimd,Cryst%ucvol)
586 
587  Cryst%angdeg(1)=ACOS(Cryst%rmet(2,3)/SQRT(Cryst%rmet(2,2)*Cryst%rmet(3,3)))/two_pi*360.0d0
588  Cryst%angdeg(2)=ACOS(Cryst%rmet(1,3)/SQRT(Cryst%rmet(1,1)*Cryst%rmet(3,3)))/two_pi*360.0d0
589  Cryst%angdeg(3)=ACOS(Cryst%rmet(1,2)/SQRT(Cryst%rmet(1,1)*Cryst%rmet(2,2)))/two_pi*360.0d0
590 
591  call xred2xcart(Cryst%natom,Cryst%rprimd,Cryst%xcart,Cryst%xred)
592 
593 end subroutine crystal_compute_geometry

m_crystal/crystal_compute_sym [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

  crystal_compute_sym

FUNCTION

  Get symmetries in cartesian coordinates, construct rotation tables
  for atoms with and without spinor, and construct list of reducible atoms.

SOURCE

492 subroutine crystal_compute_sym(Cryst)
493 
494 !Arguments ------------------------------------
495  class(crystal_t),intent(inout) :: Cryst
496 
497 !Local variables-------------------------------
498 !scalars
499  integer :: iat,indx,isym
500  real(dp) :: tolsym8
501  logical, allocatable :: irredat_tmp(:)
502 !arrays
503  integer :: symrec(3,3)
504 
505 ! *********************************************************************
506 
507  ! Get symmetries in reciprocal space
508  do isym=1,Cryst%nsym
509    call mati3inv(Cryst%symrel(:,:,isym),symrec)
510    Cryst%symrec(:,:,isym)=symrec
511  end do
512 
513  ! Get symmetries in cartesian coordinates
514  do isym =1,Cryst%nsym
515    call symredcart(Cryst%rprimd, Cryst%gprimd, Cryst%symrel_cart(:,:,isym), Cryst%symrel(:,:,isym))
516    ! purify operations in cartesian coordinates.
517    where (abs(Cryst%symrel_cart(:,:,isym)) < tol14)
518      Cryst%symrel_cart(:,:,isym) = zero
519    end where
520  end do
521 
522  ! === Obtain a list of rotated atoms ===
523  ! $ R^{-1} (xred(:,iat)-\tau) = xred(:,iat_sym) + R_0 $
524  ! * indsym(4,  isym,iat) gives iat_sym in the original unit cell.
525  ! * indsym(1:3,isym,iat) gives the lattice vector $R_0$.
526  !
527  tolsym8=tol8
528  call symatm(Cryst%indsym, Cryst%natom, Cryst%nsym, Cryst%symrec, Cryst%tnons, tolsym8, Cryst%typat, Cryst%xred)
529 
530  ! Rotations in spinor space
531  do isym=1,Cryst%nsym
532    call getspinrot(Cryst%rprimd, Cryst%spinrot(:,isym), Cryst%symrel(:,:,isym))
533  end do
534 
535 ! Find list of irreducible atoms by using the indsym
536  ABI_MALLOC(irredat_tmp, (Cryst%natom))
537  irredat_tmp = .TRUE.
538 
539  Cryst%nirredat = 0
540  do iat = 1,Cryst%natom
541    if(irredat_tmp(iat))then
542       Cryst%nirredat = Cryst%nirredat + 1
543       do isym = 1,Cryst%nsym
544          if (Cryst%indsym(4,isym,iat) /= iat)then
545            irredat_tmp(Cryst%indsym(4,isym,iat)) = .FALSE.
546          endif
547       enddo
548    endif
549  enddo
550 
551  ! Write indexes of irreducible atoms
552  ABI_MALLOC(Cryst%irredatindx,(Cryst%nirredat))
553  indx = 0
554  do iat = 1,Cryst%natom
555     if (irredat_tmp(iat)) then
556       indx = indx + 1
557       cryst%irredatindx(indx) = iat
558     endif
559  enddo
560 
561  ABI_SFREE(irredat_tmp)
562 
563 end subroutine crystal_compute_sym

m_crystal/crystal_copy [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

  crystal_copy

FUNCTION

  Copy the object.

 OUTPUTS
  new = A new crystal instance

SOURCE

725 subroutine crystal_copy(Cryst, new)
726 
727 !Arguments ------------------------------------
728  class(crystal_t),intent(inout) :: Cryst
729  class(crystal_t),intent(out) :: new
730 
731 ! *********************************************************************
732 
733  ! Copy dimensions, scalar variables, and static arrays
734  new%natom = Cryst%natom
735  new%nsym = Cryst%nsym
736  new%ntypat = Cryst%ntypat
737  new%nirredat = Cryst%nirredat
738  new%npsp = Cryst%npsp
739  new%space_group = Cryst%space_group
740  new%timrev = Cryst%timrev
741  new%ucvol = Cryst%ucvol
742  new%use_antiferro = Cryst%use_antiferro
743  new%angdeg = Cryst%angdeg
744  new%gmet = Cryst%gmet
745  new%gprimd = Cryst%gprimd
746  new%rmet = Cryst%rmet
747  new%rprimd = Cryst%rprimd
748 
749  ! Allocate memory
750  call new%malloc()
751  if (allocated(Cryst%irredatindx)) then
752    ABI_MALLOC(new%irredatindx,(new%nirredat))
753  end if
754 
755  ! Copy dynamic arrays
756  new%indsym = Cryst%indsym
757  new%symafm = Cryst%symafm
758  new%symrec = Cryst%symrec
759  new%symrel = Cryst%symrel
760  new%symrel_cart = Cryst%symrel_cart
761  new%atindx = Cryst%atindx
762  new%atindx1 = Cryst%atindx1
763  new%typat = Cryst%typat
764  new%nattyp = Cryst%nattyp
765  new%tnons = Cryst%tnons
766  new%xcart = Cryst%xcart
767  new%xred = Cryst%xred
768  new%spinrot = Cryst%spinrot
769  new%amu = Cryst%amu
770  new%zion = Cryst%zion
771  new%znucl = Cryst%znucl
772  new%title = Cryst%title
773  if (allocated(Cryst%irredatindx)) then
774    new%irredatindx = Cryst%irredatindx
775  end if
776 
777 end subroutine crystal_copy

m_crystal/crystal_free [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

  crystal_free

FUNCTION

  Destroy the dynamic arrays in a crystal_t data type.

SOURCE

677 subroutine crystal_free(Cryst)
678 
679 !Arguments ------------------------------------
680  class(crystal_t),intent(inout) :: Cryst
681 
682 ! *********************************************************************
683 
684 !integer
685  ABI_SFREE(Cryst%indsym)
686  ABI_SFREE(Cryst%symafm)
687  ABI_SFREE(Cryst%symrec)
688  ABI_SFREE(Cryst%symrel)
689  ABI_SFREE(Cryst%symrel_cart)
690  ABI_SFREE(Cryst%atindx)
691  ABI_SFREE(Cryst%atindx1)
692  ABI_SFREE(Cryst%typat)
693  ABI_SFREE(Cryst%nattyp)
694  ABI_SFREE(Cryst%irredatindx)
695 
696 !real
697  ABI_SFREE(Cryst%tnons)
698  ABI_SFREE(Cryst%xcart)
699  ABI_SFREE(Cryst%xred)
700  ABI_SFREE(Cryst%zion)
701  ABI_SFREE(Cryst%znucl)
702  ABI_SFREE(Cryst%amu)
703  ABI_SFREE(Cryst%spinrot)
704 
705 !character
706  ABI_SFREE(Cryst%title)
707 
708 end subroutine crystal_free

m_crystal/crystal_index_atoms [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

  crystal_index_atoms

FUNCTION

  Generate index table of atoms, in order for them to be used type after type.

SOURCE

452 subroutine crystal_index_atoms(Cryst)
453 
454 !Arguments ------------------------------------
455  class(crystal_t),intent(inout) :: Cryst
456 
457 !Local variables-------------------------------
458 !scalars
459  integer :: iat,indx,itypat
460 
461 ! *********************************************************************
462 
463  !ABI_MALLOC(Cryst%nattyp,(Cryst%ntypat))
464  indx=1
465  do itypat=1,Cryst%ntypat
466    Cryst%nattyp(itypat)=0
467    do iat=1,Cryst%natom
468      if (Cryst%typat(iat)==itypat) then
469        Cryst%atindx (iat )=indx
470        Cryst%atindx1(indx)=iat
471        indx=indx+1
472        Cryst%nattyp(itypat)=Cryst%nattyp(itypat)+1
473      end if
474    end do
475  end do
476 
477 end subroutine crystal_index_atoms

m_crystal/crystal_init [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

  crystal_init

FUNCTION

  Initialize a crystal_t data type.
  Ideally the routine should work in two different modes:
  Either the symmetries are directly supplied or the space group
  is determined starting from the definition of the unit cell.
  Only the first method is implemented, the second one should be
  a wrapper for the symmetry finder library. To implement the
  second case I have to add additional entries in the object
  and I have also to pass an object describing the (optional) geometry builder.

INPUTS

  natom=number of atom
  ntypat=number of type of atoms
  nsym=number of symmetry operations
  rprimd(3,3)=dimensional lattive vector (real space)
  typat(natom)=type of each atom
  xred(3,natom)=reduced coordinates of each atom
  symrel(3,3,nsym) [optional]=symmetry operations in real space
  space_group=Space group (0 if not available)
  tnons(3,nsym) [optional]=fractional Translations
  symafm(nsym) [optional]=  ferromagnetic symmetries
  remove_inv [optional]= if .TRUE. the inversion is removed from the set of symmetries
  timrev ==2 => take advantage of time-reversal symmetry
         ==1 ==> do not use time-reversal symmetry

OUTPUT

  Cryst<crystal_t>= the object completely initialized.

TODO

  Add additional entries in the class:
  1) Info on space and point group (generators?).
  2) alchemy
  3) masses and nuclear (pseudo&AE) charge
  4) forces stresses, velocities.
  5) constraints for the relaxation
  6) Likely I will need also info on the electric field and berryopt

SOURCE

320 subroutine crystal_init(amu,Cryst,space_group,natom,npsp,ntypat,nsym,rprimd,typat,xred,&
321                         zion,znucl,timrev,use_antiferro,remove_inv,title,&
322                         symrel,tnons,symafm) ! Optional
323 
324 !Arguments ------------------------------------
325 !scalars
326  integer,intent(in) :: natom,ntypat,nsym,timrev,space_group,npsp
327  type(crystal_t),intent(inout) :: Cryst
328  logical,intent(in) :: remove_inv,use_antiferro
329 !arrays
330  integer,intent(in) :: typat(natom)
331  integer,optional,intent(in) :: symrel(3,3,nsym),symafm(nsym)
332  real(dp),intent(in) :: amu(ntypat),xred(3,natom),rprimd(3,3),zion(ntypat),znucl(npsp)
333  real(dp),optional,intent(in) :: tnons(3,nsym)
334  character(len=*),intent(in) :: title(ntypat)
335 
336 !Local variables-------------------------------
337 !scalars
338  integer :: pinv,nsym_noI
339  !character(len=500) :: msg
340 !arrays
341  integer,pointer :: symrel_noI(:,:,:)
342  real(dp),pointer :: tnons_noI(:,:)
343 ! *************************************************************************
344 
345  !@crystal_t
346  Cryst%natom = natom
347  Cryst%ntypat = ntypat
348  Cryst%npsp = npsp
349  Cryst%space_group = space_group
350  Cryst%nsym = nsym
351  Cryst%timrev = timrev
352  Cryst%use_antiferro = use_antiferro
353  Cryst%rprimd = rprimd
354 
355  call Cryst%free()
356  call Cryst%malloc()
357 
358  Cryst%amu   = amu
359  Cryst%typat = typat
360  Cryst%xred  = xred
361  Cryst%zion  = zion
362  Cryst%znucl = znucl
363  Cryst%title = title
364 
365  call Cryst%compute_geometry()
366 
367  call Cryst%index_atoms()
368 
369  ! TODO: Make this more elegant
370  if (PRESENT(symrel).and.PRESENT(tnons).and.PRESENT(symafm)) then
371    if (.not.remove_inv) then
372      ! Just a copy
373      Cryst%symrel=symrel
374      Cryst%tnons=tnons
375      Cryst%symafm=symafm
376    else
377      ! Remove inversion, just to be compatible with old GW implementation
378      ! TODO should be removed!
379      call remove_inversion(nsym,symrel,tnons,nsym_noI,symrel_noI,tnons_noI,pinv)
380      Cryst%nsym=nsym_noI
381      ABI_SFREE(Cryst%symrel)
382      ABI_SFREE(Cryst%symrec)
383      ABI_SFREE(Cryst%tnons)
384      ABI_SFREE(Cryst%symafm)
385      ABI_MALLOC(Cryst%symrel,(3,3,nsym_noI))
386      ABI_MALLOC(Cryst%symrec,(3,3,nsym_noI))
387      ABI_MALLOC(Cryst%tnons,(3,nsym_noI))
388      ABI_MALLOC(Cryst%symafm,(nsym_noI))
389      Cryst%symrel=symrel_noI
390      Cryst%tnons=tnons_noI
391      if (ANY(symafm==-1)) then
392        ABI_BUG('Solve the problem with inversion before adding ferromagnetic symmetries')
393      end if
394      Cryst%symafm=1
395      ABI_FREE(symrel_noI)
396      ABI_FREE(tnons_noI)
397    end if
398 
399  else
400    ! Find symmetries symrec,symrel,tnons,symafm
401    ! TODO This should be a wrapper around the abinit library whose usage is not so straightforward
402    ABI_BUG('NotImplememented: symrel, symrec and tnons should be specied')
403  end if
404 
405  ! Compute all symetries and construct tables.
406  call Cryst%compute_sym()
407 
408 end subroutine crystal_init

m_crystal/crystal_malloc [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

  crystal_malloc

FUNCTION

  Allocate the dynamic arrays in a crystal_t data type.

SOURCE

607 subroutine crystal_malloc(Cryst)
608 
609 !Arguments ------------------------------------
610  class(crystal_t),intent(inout) :: Cryst
611 
612 !Local variables-------------------------------
613 !scalars
614  integer :: ii
615 ! *********************************************************************
616 
617 !integer
618  ABI_MALLOC(Cryst%typat,(Cryst%natom))
619  ABI_MALLOC(Cryst%xred,(3,Cryst%natom))
620  ABI_MALLOC(Cryst%xcart,(3,Cryst%natom))
621  ABI_MALLOC(Cryst%zion,(Cryst%ntypat))
622  ABI_MALLOC(Cryst%znucl,(Cryst%npsp))
623  ABI_MALLOC(Cryst%amu, (Cryst%ntypat))
624 
625  ABI_MALLOC(Cryst%symrel,(3,3,Cryst%nsym))
626  ABI_MALLOC(Cryst%symrec,(3,3,Cryst%nsym))
627  ABI_MALLOC(Cryst%tnons,(3,Cryst%nsym))
628  ABI_MALLOC(Cryst%symafm,(Cryst%nsym))
629  ABI_MALLOC(Cryst%symrel_cart, (3, 3, Cryst%nsym))
630  ABI_MALLOC(Cryst%indsym,(4, Cryst%nsym, Cryst%natom))
631 
632  ABI_MALLOC(Cryst%atindx,(Cryst%natom))
633  ABI_MALLOC(Cryst%atindx1,(Cryst%natom))
634  ABI_MALLOC(Cryst%nattyp,(Cryst%ntypat))
635  ABI_MALLOC(Cryst%spinrot, (4, Cryst%nsym))
636 
637  ABI_MALLOC(Cryst%title,(Cryst%ntypat))
638 
639  ! nirredat must first be computed from indsym
640  !ABI_MALLOC(Cryst%irredatindx,(Cryst%nirredat))
641 
642  Cryst%typat = zero
643  Cryst%xred = zero
644  Cryst%xcart = zero
645  Cryst%zion = zero
646  Cryst%znucl = zero
647  Cryst%amu = zero
648  Cryst%symrel = zero
649  Cryst%symrec = zero
650  Cryst%tnons = zero
651  Cryst%symafm = zero
652  Cryst%symrel_cart = zero
653  Cryst%indsym = zero
654  Cryst%atindx = zero
655  Cryst%atindx1 = zero
656  Cryst%nattyp = zero
657  Cryst%spinrot = zero
658 
659  do ii=1,Cryst%ntypat
660    Cryst%title(ii) = ''
661  end do
662 
663 end subroutine crystal_malloc

m_crystal/crystal_ncread [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

 crystal_ncread

FUNCTION

 Read the crystal object from a NETCDF file.

INPUTS

  cryst<crystal_t>=Object defining the unit cell and its symmetries.
  ncid=NC file handle.

OUTPUT

  crystal

SOURCE

1657 subroutine crystal_ncread(cryst, ncid)
1658 
1659 !Arguments ------------------------------------
1660 !scalars
1661  class(crystal_t),intent(inout) :: cryst
1662  integer,intent(in) :: ncid
1663 
1664 !Local variables ------------------------------------
1665 !scalars
1666  integer :: use_antiferro
1667 
1668 #ifdef HAVE_NETCDF
1669 
1670 ! *************************************************************************
1671 
1672    !NCF_CHECK(nf90_inq_varid(ncid, varname, varid))
1673 
1674  ! ---------------
1675  ! Read dimensions
1676  ! ---------------
1677  NCF_CHECK(nctk_get_dim(ncid, "number_of_atoms", cryst%natom))
1678  NCF_CHECK(nctk_get_dim(ncid, "number_of_atom_species", cryst%ntypat))
1679  NCF_CHECK(nctk_get_dim(ncid, "number_of_atom_pseudopotentials", cryst%npsp))
1680  NCF_CHECK(nctk_get_dim(ncid, "number_of_symmetry_operations", cryst%nsym))
1681 
1682  ! ---------------
1683  ! Allocate memory
1684  ! ---------------
1685  call cryst%malloc()
1686 
1687  ! ------------
1688  ! read scalars
1689  ! ------------
1690  NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, "space_group"), cryst%space_group))
1691  NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, "time_reversal"), cryst%timrev))
1692  NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, "use_antiferromagnetic_symmetries"), use_antiferro))
1693  cryst%use_antiferro = .False.
1694  if (use_antiferro/=0) cryst%use_antiferro = .True.
1695 
1696  ! -----------
1697  ! read arrays
1698  ! -----------
1699  NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, "primitive_vectors"), cryst%rprimd))
1700  NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, "reduced_symmetry_matrices"), cryst%symrel))
1701  NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, "reduced_symmetry_translations"), cryst%tnons))
1702  NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, "atom_species"), cryst%typat))
1703  NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, "reduced_atom_positions"), cryst%xred))
1704  NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, "atomic_numbers"), cryst%znucl(1:cryst%ntypat)))
1705  NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, "atomic_mass_units"), cryst%amu))
1706 
1707  NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, "symafm"), cryst%symafm))
1708  NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, "symrel_cart"), cryst%symrel_cart))
1709  NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, "indsym"), cryst%indsym))
1710 
1711  if (cryst%npsp == cryst%ntypat) then
1712    NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, "valence_charges"), cryst%zion))
1713  end if
1714 
1715  ! Ignore those
1716  !NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, "pseudopotential_types"), psp_desc))
1717  !NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, "atom_species_names"), symbols_long))
1718  !NCF_CHECK(nf90_get_var(ncid, nctk_idname(ncid, "chemical_symbols"), symbols))
1719 
1720  ! -----------------------
1721  ! Complete initialization
1722  ! -----------------------
1723 
1724  call Cryst%compute_geometry()
1725  call Cryst%index_atoms()
1726  call Cryst%compute_sym()
1727 
1728 #else
1729  ABI_ERROR("netcdf library not available")
1730 #endif
1731 
1732 end subroutine crystal_ncread

m_crystal/crystal_ncwrite [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

 crystal_ncwrite

FUNCTION

 Output system geometry to a file, using the NETCDF file format and ETSF I/O.
 Data are taken from the crystal_t object.

INPUTS

  cryst<crystal_t>=Object defining the unit cell and its symmetries.
  ncid=NC file handle.

OUTPUT

  Only writing

NOTES

  Alchemy not treated, since crystal should be initialized at the beginning of the run.

SOURCE

1461 integer function crystal_ncwrite(cryst, ncid) result(ncerr)
1462 
1463 !Arguments ------------------------------------
1464 !scalars
1465  class(crystal_t),intent(in) :: cryst
1466  integer,intent(in) :: ncid
1467 
1468 !Local variables-------------------------------
1469 !scalars
1470  integer :: itypat
1471  character(len=500) :: msg
1472  character(len=etsfio_charlen) :: symmorphic
1473  type(atomdata_t) :: atom
1474 !arrays
1475  character(len=2) :: symbols(cryst%ntypat)
1476  character(len=80) :: psp_desc(cryst%ntypat),symbols_long(cryst%ntypat)
1477 
1478 ! *************************************************************************
1479 
1480 #ifdef HAVE_NETCDF
1481 
1482  ! TODO alchemy not treated correctly by ETSF_IO specs.
1483  if (cryst%isalchemical()) then
1484    write(msg,"(3a)")&
1485     "Alchemical crystals are not fully supported by the netcdf format",ch10,&
1486     "Important parameters (e.g. znucl, symbols) are not written with the correct value"
1487    ABI_WARNING(msg)
1488  end if
1489 
1490  symmorphic = yesno(cryst%isymmorphic())
1491 
1492  ! Define dimensions.
1493  ! npsp added in v9.
1494  ncerr = nctk_def_dims(ncid, [ &
1495    nctkdim_t("complex", 2), nctkdim_t("symbol_length", 2),&
1496    nctkdim_t("character_string_length", 80), nctkdim_t("number_of_cartesian_directions", 3),&
1497    nctkdim_t("number_of_reduced_dimensions", 3), nctkdim_t("number_of_vectors", 3),&
1498    nctkdim_t("number_of_atoms", cryst%natom), nctkdim_t("number_of_atom_species", cryst%ntypat),&
1499    nctkdim_t("number_of_atom_pseudopotentials", cryst%npsp),&
1500    nctkdim_t("number_of_symmetry_operations", cryst%nsym)], defmode=.True.)
1501  NCF_CHECK(ncerr)
1502 
1503  ! Define variables
1504  ncerr = nctk_def_iscalars(ncid, [character(len=nctk_slen) :: &
1505     "space_group", "time_reversal", "use_antiferromagnetic_symmetries"])
1506  NCF_CHECK(ncerr)
1507 
1508  ncerr = nctk_def_arrays(ncid, [ &
1509   ! Atomic structure and symmetry operations
1510   nctkarr_t("primitive_vectors", "dp", "number_of_cartesian_directions, number_of_vectors"), &
1511   nctkarr_t("reduced_symmetry_matrices", "int", &
1512     "number_of_reduced_dimensions, number_of_reduced_dimensions, number_of_symmetry_operations"), &
1513   nctkarr_t("reduced_symmetry_translations", "dp", "number_of_reduced_dimensions, number_of_symmetry_operations"), &
1514   nctkarr_t("atom_species", "int", "number_of_atoms"), &
1515   nctkarr_t("reduced_atom_positions", "dp", "number_of_reduced_dimensions, number_of_atoms"), &
1516   nctkarr_t("atomic_numbers", "dp", "number_of_atom_species"), &
1517   nctkarr_t("atom_species_names", "char", "character_string_length, number_of_atom_species"), &
1518   nctkarr_t("chemical_symbols", "char", "symbol_length, number_of_atom_species"), &
1519   nctkarr_t('atomic_mass_units', "dp", "number_of_atom_species"), &
1520   ! Atomic information.
1521   nctkarr_t("valence_charges", "dp", "number_of_atom_species"), &  ! NB: This variable is not written if alchemical
1522   nctkarr_t("pseudopotential_types", "char", "character_string_length, number_of_atom_species") &
1523  ])
1524  NCF_CHECK(ncerr)
1525 
1526  ! Some variables require the "symmorphic" attribute.
1527  NCF_CHECK(nf90_put_att(ncid, vid("reduced_symmetry_matrices"), "symmorphic", symmorphic))
1528  NCF_CHECK(nf90_put_att(ncid, vid("reduced_symmetry_translations"), "symmorphic", symmorphic))
1529 
1530  ! At this point we have an ETSF-compliant file. Add additional data for internal use in abinit.
1531  ncerr = nctk_def_arrays(ncid, [ &
1532    nctkarr_t('symafm', "int", "number_of_symmetry_operations"), &
1533    nctkarr_t('symrel_cart', "dp", "three, three, number_of_symmetry_operations"), &
1534    nctkarr_t('indsym', "int", "four, number_of_symmetry_operations, number_of_atoms") &
1535  ])
1536  NCF_CHECK(ncerr)
1537 
1538  ! Set-up atomic symbols.
1539  do itypat=1,cryst%ntypat
1540    call atomdata_from_znucl(atom, cryst%znucl(itypat))
1541    symbols(itypat) = atom%symbol
1542    write(symbols_long(itypat),'(a2,a78)') symbols(itypat),REPEAT(CHAR(0),78)
1543    write(psp_desc(itypat),'(2a)') &
1544      cryst%title(itypat)(1:MIN(80,LEN_TRIM(cryst%title(itypat)))),REPEAT(CHAR(0),MAX(0,80-LEN_TRIM(cryst%title(itypat))))
1545  end do
1546 
1547  ! Write data.
1548  NCF_CHECK(nctk_set_datamode(ncid))
1549  NCF_CHECK(nf90_put_var(ncid, vid("space_group"), cryst%space_group))
1550  NCF_CHECK(nf90_put_var(ncid, vid("primitive_vectors"), cryst%rprimd))
1551  NCF_CHECK(nf90_put_var(ncid, vid("reduced_symmetry_matrices"), cryst%symrel))
1552  NCF_CHECK(nf90_put_var(ncid, vid("reduced_symmetry_translations"), cryst%tnons))
1553  NCF_CHECK(nf90_put_var(ncid, vid("atom_species"), cryst%typat))
1554  NCF_CHECK(nf90_put_var(ncid, vid("reduced_atom_positions"), cryst%xred))
1555  NCF_CHECK(nf90_put_var(ncid, vid("atomic_numbers"), cryst%znucl(1:cryst%ntypat)))
1556  NCF_CHECK(nf90_put_var(ncid, vid("atom_species_names"), symbols_long))
1557  NCF_CHECK(nf90_put_var(ncid, vid("chemical_symbols"), symbols))
1558  NCF_CHECK(nf90_put_var(ncid, vid('atomic_mass_units'), cryst%amu))
1559  NCF_CHECK(nf90_put_var(ncid, vid("pseudopotential_types"), psp_desc))
1560  if (cryst%npsp == cryst%ntypat) then
1561    NCF_CHECK(nf90_put_var(ncid, vid("valence_charges"), cryst%zion))
1562  end if
1563 
1564  NCF_CHECK(nf90_put_var(ncid, vid("symafm"), cryst%symafm))
1565  NCF_CHECK(nf90_put_var(ncid, vid("symrel_cart"), cryst%symrel_cart))
1566  NCF_CHECK(nf90_put_var(ncid, vid("indsym"), cryst%indsym))
1567 
1568 ! Variables pertaining to the symmetry of the wavefunctions.
1569 ! Note that these variables will be used in crystal_compare
1570  NCF_CHECK(nf90_put_var(ncid, vid("time_reversal"), cryst%timrev))
1571 
1572  if (cryst%use_antiferro) then
1573     NCF_CHECK(nf90_put_var(ncid, vid("use_antiferromagnetic_symmetries"), 1))
1574  else
1575     NCF_CHECK(nf90_put_var(ncid, vid("use_antiferromagnetic_symmetries"), 0))
1576  end if
1577 
1578 
1579 #else
1580  ABI_ERROR("netcdf library not available")
1581 #endif
1582 
1583 contains
1584  integer function vid(vname)
1585    character(len=*),intent(in) :: vname
1586    vid = nctk_idname(ncid, vname)
1587  end function vid
1588 
1589 end function crystal_ncwrite

m_crystal/crystal_ncwrite_path [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

 crystal_ncwrite_path

FUNCTION

 Output system geometry to a file, using the NETCDF file format and ETSF I/O.

INPUTS

  crystal<crystal_t>=Object defining the unit cell and its symmetries.
  path=filename

OUTPUT

  Only writing

SOURCE

1610 integer function crystal_ncwrite_path(crystal, path) result(ncerr)
1611 
1612 !Arguments ------------------------------------
1613 !scalars
1614  character(len=*),intent(in) :: path
1615  class(crystal_t),intent(in) :: crystal
1616 
1617 #ifdef HAVE_NETCDF
1618 !Local variables-------------------------------
1619 !scalars
1620  integer :: ncid
1621 
1622 ! *************************************************************************
1623 
1624  ncerr = nf90_noerr
1625  if (file_exists(path)) then
1626    NCF_CHECK(nctk_open_modify(ncid, path, xmpi_comm_self))
1627  else
1628    ncerr = nctk_open_create(ncid, path, xmpi_comm_self)
1629    NCF_CHECK_MSG(ncerr, sjoin("creating:", path))
1630  end if
1631 
1632  NCF_CHECK(crystal_ncwrite(crystal, ncid))
1633  NCF_CHECK(nf90_close(ncid))
1634 #endif
1635 
1636 end function crystal_ncwrite_path

m_crystal/crystal_point_group [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

  crystal_point_group

FUNCTION

  Return the symmetries of the point group of the crystal.

INPUTS

  [include_timrev]=If True, time-reversal symmetry is included in the point group unless
    the system has spatial inversion. Default: False

OUTPUT

  ptg_nsym=Number of symmetries in the point group
  ptg_symrel(3,3,ptg_nsym)=Rotations in real space
  ptg_symrec(3,3,ptg_nsym)=Rotations in reciprocal space
  has_inversion=True if spatial inversion is present in the point group.

SOURCE

1367 subroutine crystal_point_group(cryst, ptg_nsym, ptg_symrel, ptg_symrec, has_inversion, include_timrev)
1368 
1369 !Arguments ------------------------------------
1370 !scalars
1371  class(crystal_t),intent(in) :: cryst
1372  integer,intent(out) :: ptg_nsym
1373  logical,optional,intent(in) :: include_timrev
1374  logical,intent(out) :: has_inversion
1375 !arrays
1376  integer,allocatable,intent(out) :: ptg_symrel(:,:,:),ptg_symrec(:,:,:)
1377 
1378 !Local variables-------------------------------
1379 !scalars
1380  integer :: isym, search, tmp_nsym, ierr
1381  logical :: found, my_include_timrev, debug
1382 !arrays
1383  integer :: work_symrel(3,3,cryst%nsym)
1384  integer,allocatable :: symafm(:)
1385 
1386 ! *************************************************************************
1387 
1388  my_include_timrev = .False.; if (present(include_timrev)) my_include_timrev = include_timrev
1389 
1390  tmp_nsym = 1; work_symrel(:,:,1) = cryst%symrel(:,:,1)
1391  do isym=2,cryst%nsym
1392    if (cryst%symafm(isym) == -1) cycle
1393    do search=1,tmp_nsym
1394      found = all(work_symrel(:,:,search) == cryst%symrel(:,:,isym))
1395      if (found) exit
1396    end do
1397    if (.not. found) then
1398      tmp_nsym = tmp_nsym + 1
1399      work_symrel(:,:,tmp_nsym) = cryst%symrel(:,:,isym)
1400    end if
1401  end do
1402 
1403  has_inversion = .False.
1404  do isym=1,tmp_nsym
1405    if (all(work_symrel(:,:,isym) == inversion_3d) ) then
1406      has_inversion = .True.; exit
1407    end if
1408  end do
1409 
1410  ! Now we know the symmetries of the point group.
1411  ptg_nsym = tmp_nsym; if (.not. has_inversion .and. my_include_timrev) ptg_nsym = 2 * tmp_nsym
1412  ABI_MALLOC(ptg_symrel, (3,3,ptg_nsym))
1413  ABI_MALLOC(ptg_symrec, (3,3,ptg_nsym))
1414 
1415  ptg_symrel(:,:,1:tmp_nsym) = work_symrel(:,:,1:tmp_nsym)
1416  do isym=1,tmp_nsym
1417    call mati3inv(ptg_symrel(:,:,isym), ptg_symrec(:,:,isym))
1418  end do
1419 
1420  if (.not. has_inversion .and. my_include_timrev) then
1421    ptg_symrel(:,:,tmp_nsym+1:) = -work_symrel(:,:,1:tmp_nsym)
1422    do isym=tmp_nsym+1,ptg_nsym
1423      call mati3inv(ptg_symrel(:,:,isym), ptg_symrec(:,:,isym))
1424    end do
1425  end if
1426 
1427  debug = .False.
1428  if (debug) then
1429    ABI_MALLOC(symafm, (ptg_nsym))
1430    symafm = 1
1431    call sg_multable(ptg_nsym, symafm, ptg_symrel, ierr)
1432    ABI_CHECK(ierr == 0, "point group is not a group! See messages above")
1433    ABI_FREE(symafm)
1434  end if
1435 
1436 end subroutine crystal_point_group

m_crystal/crystal_print [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

  crystal_print

FUNCTION

  Print the content of crystal_t data type

INPUTS

  Cryst<crystal_t>=The structure.
  [unit]=Unit number for output. Defaults to std_out
  [prtvol]=Verbosity level. If prtvol== -1, only lattice parameters are printed. Defaults to 0
  [mode_paral]=Either "COLL" or "PERS"
  [header]=String to be printed as header for additional info.

OUTPUT

  Only printing

SOURCE

 978 subroutine crystal_print(Cryst, header, unit, mode_paral, prtvol)
 979 
 980 !Arguments ------------------------------------
 981 !scalars
 982  integer,optional,intent(in) :: unit, prtvol
 983  character(len=*),optional,intent(in) :: mode_paral
 984  character(len=*),optional,intent(in) :: header
 985  class(crystal_t),intent(in) :: Cryst
 986 
 987 !Local variables-------------------------------
 988  integer :: my_unt,my_prtvol,nu,iatom, isym, ii, nsym
 989  character(len=4) :: my_mode
 990  character(len=500) :: msg
 991 ! *********************************************************************
 992 
 993  my_unt   =std_out; if (PRESENT(unit      )) my_unt   =unit
 994  my_prtvol=0      ; if (PRESENT(prtvol    )) my_prtvol=prtvol
 995  my_mode  ='COLL' ; if (PRESENT(mode_paral)) my_mode  =mode_paral
 996 
 997  msg=' ==== Info on the Cryst% object ==== '
 998  if (PRESENT(header)) msg=' ==== '//TRIM(ADJUSTL(header))//' ==== '
 999  call wrtout(my_unt, sjoin(ch10, msg), my_mode)
1000 
1001  write(msg,'(a)')' Real(R)+Recip(G) space primitive vectors, cartesian coordinates (Bohr,Bohr^-1):'
1002  call wrtout(my_unt,msg,my_mode)
1003  do nu=1,3
1004    write(msg,'(1x,a,i1,a,3f11.7,2x,a,i1,a,3f11.7)')&
1005     'R(',nu,')=',Cryst%rprimd(:,nu)+tol10, &
1006     'G(',nu,')=',Cryst%gprimd(:,nu)+tol10  ! tol10 is used to be consistent with metric.F90
1007    call wrtout(my_unt,msg,my_mode)
1008  end do
1009 
1010  write(msg,'(a,1p,e15.7,a)')' Unit cell volume ucvol=',Cryst%ucvol+tol10,' bohr^3'
1011  call wrtout(my_unt,msg,my_mode)
1012 
1013  write(msg,'(a,3es16.8,a)')' Angles (23,13,12)=',Cryst%angdeg(1:3),' degrees'
1014  call wrtout(my_unt,msg,my_mode)
1015 
1016  if (Cryst%timrev==1) then
1017    msg = ' Time-reversal symmetry is not present '
1018  else if (Cryst%timrev==2) then
1019    msg = ' Time-reversal symmetry is present '
1020  else
1021    ABI_BUG(sjoin('Wrong value for timrev:', itoa(cryst%timrev)))
1022  end if
1023  call wrtout(my_unt,msg,my_mode)
1024  if (my_prtvol == -1) return
1025 
1026  if (my_prtvol > 0) then
1027    call print_symmetries(Cryst%nsym, Cryst%symrel, Cryst%tnons, Cryst%symafm, unit=my_unt, mode_paral=my_mode)
1028    if (Cryst%use_antiferro) call wrtout(my_unt,' System has magnetic symmetries ',my_mode)
1029 
1030    ! Print indsym using the same format as in symatm
1031    nsym = cryst%nsym
1032    do iatom=1,cryst%natom
1033      write(msg, '(a,i0,a)' )' symatm: atom number ',iatom,' is reached starting at atom'
1034      call wrtout(std_out, msg)
1035      do ii=1,(nsym-1)/24+1
1036        if (cryst%natom<100) then
1037          write(msg, '(1x,24i3)' ) (cryst%indsym(4,isym,iatom),isym=1+(ii-1)*24,min(nsym,ii*24))
1038        else
1039          write(msg, '(1x,24i6)' ) (cryst%indsym(4,isym,iatom),isym=1+(ii-1)*24,min(nsym,ii*24))
1040        end if
1041        call wrtout(std_out, msg)
1042      end do
1043    end do
1044 
1045  end if
1046 
1047  call wrtout(my_unt, " Reduced atomic positions [iatom, xred, symbol]:", my_mode)
1048  do iatom=1,cryst%natom
1049    write(msg,"(i5,a,2x,3f11.7,2x,a)")iatom,")",cryst%xred(:,iatom), cryst%symbol_type(cryst%typat(iatom))
1050    call wrtout(my_unt,msg,my_mode)
1051  end do
1052 
1053 end subroutine crystal_print

m_crystal/crystal_print_abivars [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

  crystal_print_abivars

FUNCTION

   Print unit cell info in Abinit/abivars format

INPUTS

  unit=Output unit

OUTPUT

  Only printing

SOURCE

1071 subroutine crystal_print_abivars(cryst, unit)
1072 
1073 !Arguments ------------------------------------
1074 !scalars
1075  integer,intent(in) :: unit
1076  class(crystal_t),intent(in) :: cryst
1077 
1078 !Local variables-------------------------------
1079  integer :: iatom, ii
1080  !character(len=500) :: fmt
1081 ! *********************************************************************
1082 
1083  if (unit == dev_null) return
1084 
1085  ! Write variables using standard Abinit input format.
1086  write(unit, "(/,/,a)")" # Abinit variables"
1087  write(unit, "(a)")" acell 1.0 1.0 1.0"
1088  write(unit, "(a)")" rprimd"
1089  do ii=1,3
1090     write(unit, "(3(f11.7,1x))")cryst%rprimd(:, ii)
1091  end do
1092  write(unit, "(a, i0)")" natom ", cryst%natom
1093  write(unit, "(a, i0)")" ntypat ", cryst%ntypat
1094  write(unit, strcat("(a, ", itoa(cryst%natom), "(i0,1x))")) " typat ", cryst%typat
1095  write(unit, strcat("(a, ", itoa(cryst%npsp), "(f5.1,1x))")) " znucl ", cryst%znucl
1096  write(unit, "(a)")" xred"
1097  do iatom=1,cryst%natom
1098    write(unit,"(1x, 3f11.7,2x,2a)")cryst%xred(:,iatom), " # ", cryst%symbol_type(cryst%typat(iatom))
1099  end do
1100 
1101  ! Write variables using the abivars format supported by structure variable.
1102  !write(unit, "(/,/,a)")" # Abivars format (external file with structure variable)"
1103  !write(unit, "(a)")" acell 1.0 1.0 1.0"
1104  !write(unit, "(a)")" rprimd"
1105  !do ii=1,3
1106  !   write(unit, "(1x, 3(f11.7,1x))")cryst%rprimd(:, ii)
1107  !end do
1108  !write(unit, "(a, i0)")" natom ", cryst%natom
1109  !write(unit, "(a)")" xred_symbols"
1110  !do iatom=1,cryst%natom
1111  !  write(unit,"(1x, 3f11.7,2x,a)")cryst%xred(:,iatom), cryst%symbol_type(cryst%typat(iatom))
1112  !end do
1113 
1114 end subroutine crystal_print_abivars

m_crystal/crystal_symmetrize_cart_tens33 [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

  crystal_symmetrize_cart_tens33

FUNCTION

  Symmetrize a cartesian 3x3 tensor
  Use spatial and time-reversal symmetry (TR) to symmetrize a Cartesian 3x3 tensor of real elements.

INPUTS

  v(3)=Vector in Cartesian coordinates.
  time_opt=Prefacator that defines how the vectors transforms under TR. Usually +1 or -1
      Note that TR is used only if cryst%timrev == 2. time_opt = 0 disables TR for testing purposes.

SOURCE

2123 function crystal_symmetrize_cart_tens33(cryst, t, time_opt) result(tsum)
2124 
2125 !Arguments ------------------------------------
2126  class(crystal_t),intent(in) :: cryst
2127  real(dp),intent(in) :: t(3,3)
2128  integer,intent(in) :: time_opt
2129  real(dp) :: tsum(3,3)
2130 
2131 !Local variables-------------------------------
2132  integer :: isym, itime, nsym_sum
2133  real(dp) :: tsym(3,3), tsign
2134 ! *************************************************************************
2135 
2136  tsum = zero; nsym_sum = 0
2137 
2138  do itime=1,cryst%timrev
2139    tsign = 1
2140    if (itime == cryst%timrev) then
2141      if (time_opt == 0) cycle
2142      tsign = time_opt
2143    end if
2144    do isym=1, cryst%nsym
2145      nsym_sum = nsym_sum + 1
2146      tsym = tsign * matmul((cryst%symrel_cart(:,:,isym)), matmul(t, transpose(cryst%symrel_cart(:,:,isym))))
2147      tsum = tsum + tsym
2148    end do
2149  end do
2150 
2151  tsum = tsum / nsym_sum
2152 
2153 end function crystal_symmetrize_cart_tens33

m_crystal/crystal_symmetrize_cart_vec3 [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

  crystal_symmetrize_cart_vec3

FUNCTION

  Use spatial and time-reversal symmetry (TR) to symmetrize a Cartesian vector of real elements.

INPUTS

  v(3)=Vector in Cartesian coordinates.
  time_opt=Prefactor that defines how the vectors transforms under TR. Usually +1 or -1
      Note that TR is used only if cryst%timrev == 2. time_opt = 0 disables TR for testing purposes.

SOURCE

2075 function crystal_symmetrize_cart_vec3(cryst, v, time_opt) result(vsum)
2076 
2077 !Arguments ------------------------------------
2078  class(crystal_t),intent(in) :: cryst
2079  real(dp),intent(in) :: v(3)
2080  integer,intent(in) :: time_opt
2081  real(dp) :: vsum(3)
2082 
2083 !Local variables-------------------------------
2084  integer :: isym, itime, nsym_sum
2085  real(dp) :: vsym(3), tsign
2086 ! *************************************************************************
2087 
2088  vsum = zero; nsym_sum = 0
2089  do itime=1,cryst%timrev
2090    tsign = 1
2091    if (itime == cryst%timrev) then
2092      if (time_opt == 0) cycle
2093      tsign = time_opt
2094    end if
2095    do isym=1, cryst%nsym
2096      nsym_sum = nsym_sum + 1
2097      vsym = matmul(cryst%symrel_cart(:,:,isym), v) * tsign
2098      vsum = vsum + vsym
2099    end do
2100  end do
2101  vsum = vsum / nsym_sum
2102 
2103 end function crystal_symmetrize_cart_vec3

m_crystal/crystal_t [ Types ]

[ Top ] [ m_crystal ] [ Types ]

NAME

 crystal_t

FUNCTION

 Structure defining the unit cell (geometry, atomic positions and symmetry operations in real and reciprocal space)

SOURCE

 60  type,public :: crystal_t
 61 
 62 !scalars
 63   !integer :: point_group                    ! Point group
 64   !integer :: bravais,crystsys               ! Bravais lattice, Crystal system
 65   !integer :: nptsym                         ! No of point symmetries of the Bravais lattice
 66   !integer :: bravais(11)                    ! bravais(1)=iholohedry, bravais(2)=center
 67                                              ! bravais(3:11)=coordinates of rprim in the axes of the conventional
 68                                              ! bravais lattice (*2 if center/=0)
 69   !integer,pointer ptsymrel(:,:,:)
 70   !ptsymrel(3,3,nptsym)
 71   ! nptsym point-symmetry operations of the Bravais lattice in real space in terms of primitive translations.
 72 
 73   integer :: natom
 74   ! Number of atoms
 75 
 76   integer :: nsym
 77   ! Number of symmetry operations
 78 
 79   integer :: ntypat
 80   ! Number of type of atoms
 81 
 82   integer :: nirredat
 83   ! Number of irreducibel atoms
 84 
 85   integer :: npsp
 86   ! No. of pseudopotentials
 87 
 88   integer :: space_group
 89   ! Space group
 90 
 91   integer :: timrev
 92   ! TODO BE CAREFUL here, as the convention used in abinit is different.
 93   ! 1 => do not use time-reversal symmetry.
 94   ! 2 => take advantage of time-reversal symmetry.
 95 
 96   real(dp) :: ucvol
 97   ! Real space unit cell volume.
 98 
 99   !real(dp) :: bzvol
100   ! Real space unit cell volume.
101 
102   logical :: use_antiferro
103   ! .TRUE. if AFM symmetries are present and used.
104 
105 !arrays
106   real(dp) :: angdeg(3)
107   ! Angles among rprim (degree).
108 
109   real(dp) :: gmet(3,3)
110   ! Reciprocal space metric ($\textrm{bohr}^{-2}$).
111 
112   real(dp) :: gprimd(3,3)
113   ! Dimensional primitive translations for reciprocal space ($\textrm{bohr}^{-1}$)
114 
115   real(dp) :: rmet(3,3)
116   ! Metric in real space.
117 
118   real(dp) :: rprimd(3,3)
119   ! Direct lattice vectors, Bohr units.
120 
121   integer,allocatable :: indsym(:,:,:)
122   ! indsym(4,nsym,natom)
123   ! indirect indexing array for atoms, see symatm.F90.
124 
125   integer,allocatable :: symafm(:)
126   ! symafm(nsym)
127   ! (Anti)Ferromagnetic symmetries. +1/-1
128 
129   integer,allocatable :: symrec(:,:,:)
130   ! symrec(3,3,nsym)
131   ! Symmetry operation in reciprocal space (reduced coordinates)
132 
133   integer,allocatable :: symrel(:,:,:)
134   ! symrel(3,3,nsym)
135   ! Symmetry operations in direct space (reduced coordinates).
136 
137   real(dp),allocatable :: symrel_cart(:,:,:)
138   ! symrel_cart(3,3,nsym)
139   ! Symmetry operations in cartesian coordinates (same order as symrel)
140 
141   integer,allocatable :: atindx(:)
142   integer,allocatable :: atindx1(:)
143   ! atindx(natom), atindx1(natom)
144   ! Index tables for atoms useful to treat atoms type after type.
145 
146   integer,allocatable :: typat(:)
147   integer,allocatable :: nattyp(:)
148   ! typat(natom), nattyp(ntypat)
149   ! Type of each natom and number of atoms of each type.
150 
151   integer,allocatable :: irredatindx(:)
152   ! Index of irreducible atoms
153 
154   real(dp),allocatable :: tnons(:,:)
155   ! tnons(3,nsym)
156   ! Fractional translations (reduced coordinates)
157 
158   real(dp),allocatable :: xcart(:,:)
159   ! xcart(3,natom)
160   ! Cartesian coordinates.
161 
162   real(dp),allocatable :: xred(:,:)
163   ! xred(3,natom)
164   ! Reduced coordinates.
165 
166   real(dp),allocatable :: spinrot(:,:)
167   ! spinrot(4,nsym)
168   ! spinor rotation matrices.
169 
170   real(dp),allocatable :: amu(:)
171   !  amu(ntypat)
172   !  mass of the atoms (atomic mass unit)
173 
174   real(dp),allocatable :: zion(:)
175   ! zion(ntypat)
176   ! Charge of the pseudo-ion
177   ! (No of valence electrons needed to screen exactly the pseudopotential).
178 
179   !real(dp),allocatable :: znucltypat(:)
180    ! znucltypat(ntypat)
181    ! The atomic number of each type of atom (might be alchemy wrt psps)
182 
183   real(dp),allocatable :: znucl(:)
184   ! znucl(npsp)
185   ! Nuclear charge for each type of pseudopotential
186 
187   character(len=132),allocatable :: title(:)
188    ! title(ntypat)
189    ! The content of first line read from the psp file
190 
191  contains
192 
193    procedure :: ncwrite => crystal_ncwrite
194    ! Write the object in netcdf format
195 
196    procedure :: ncwrite_path => crystal_ncwrite_path
197    ! Dump the object to netcdf file.
198 
199    procedure :: ncread => crystal_ncread
200    ! Read the object from a netcdf file.
201 
202    procedure :: isymmorphic
203    ! True if space group is symmorphic.
204 
205    procedure :: idx_spatial_inversion
206    ! Return the index of the spatial inversion, 0 if not present.
207 
208    procedure :: isalchemical
209    ! True if we are using alchemical pseudopotentials.
210 
211    procedure :: malloc => crystal_malloc
212    ! Allocate memory.
213 
214    procedure :: free => crystal_free
215    ! Free memory.
216 
217    procedure :: copy => crystal_copy
218    ! Copy object.
219 
220    procedure :: bcast => crystal_bcast
221    ! Master broadcasts data and others allocate their arrays.
222 
223    procedure :: new_without_symmetries => crystal_without_symmetries
224    ! Return new object without symmetries (actually nsym = 1 and identity operation)
225 
226    procedure :: get_point_group => crystal_point_group
227    ! Return the symmetries of the point group of the crystal.
228 
229    procedure :: index_atoms => crystal_index_atoms
230    ! Generate index table of atoms.
231 
232    procedure :: compute_sym => crystal_compute_sym
233    ! Compute all symetries and construct tables.
234 
235    procedure :: compute_geometry => crystal_compute_geometry
236    ! Compute the different metrics and the angle between primitive vectors.
237 
238    procedure :: symbol_type
239    ! Return the atomic symbol from the itypat index.
240 
241    procedure :: symbol_iatom
242    ! Return the atomic symbol from the iatom index.
243 
244    procedure :: adata_type
245    ! Return atomic data from the itypat index.
246 
247    procedure :: compare => crystal_compare
248    ! Compare two crystalline structures, write warning messages if they differ, return exit status
249 
250    procedure :: print => crystal_print
251    ! Print dimensions and basic info stored in the object
252 
253    procedure :: print_abivars => crystal_print_abivars
254    ! Print unit cell info in Abinit/abivars format
255 
256    procedure :: symmetrize_cart_vec3 => crystal_symmetrize_cart_vec3
257    ! Symmetrize a 3d cartesian vector
258 
259    procedure :: symmetrize_cart_tens33 => crystal_symmetrize_cart_tens33
260    ! Symmetrize a cartesian 3x3 tensor
261 
262    procedure :: get_redcart_qdirs => get_redcart_qdirs
263    ! Return predefined list of 6 q-versors in reciprocal space reduced coordinates.
264 
265  end type crystal_t
266 
267  public :: crystal_init            ! Main Creation method.
268  public :: crystal_free            ! Main Destruction method.
269  public :: symbols_crystal         ! Return an array with the atomic symbol:["Sr","Ru","O1","O2","O3"]
270  public :: prt_cif                 ! Print CIF file.
271  public :: prtposcar               ! output VASP style POSCAR and FORCES files.

m_crystal/crystal_without_symmetries [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

  crystal_without_symmetries

FUNCTION

  Return new crystal_t object without symmetries (actually nsym = 1 and identity operation)

INPUTS

OUTPUT

SOURCE

424 type(crystal_t) function crystal_without_symmetries(self) result(new)
425 
426 !Arguments ------------------------------------
427  class(crystal_t), intent(in) :: self
428 
429 !Local variables-------------------------------
430  integer,parameter :: timrev1 = 1, new_symafm(1) = 1
431  real(dp),parameter :: new_tnons(3,1) = zero
432 ! *************************************************************************
433 
434  call crystal_init(self%amu, new, 1, self%natom, self%npsp, self%ntypat, 1, self%rprimd, self%typat, &
435   self%xred, self%zion, self%znucl, timrev1, .False., .False., self%title, &
436   symrel=identity_3d, tnons=new_tnons, symafm=new_symafm)
437 
438 end function crystal_without_symmetries

m_crystal/get_recart_qdirs [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

  get_recart_qdirs

FUNCTION

  Return predefined list of 6 q-versors in reciprocal space reduced coordinates.
  First 3 entries are along the recip. space lattice vectors, then along the Cartesian axis x,y,z.
  The optional qlen argument, can be used to rescale the vectors. Default: 1

INPUTS

SOURCE

2171 subroutine get_redcart_qdirs(cryst, nq, qdirs, qlen)
2172 
2173  class(crystal_t),intent(in) :: cryst
2174  integer,intent(out) :: nq
2175  real(dp),allocatable,intent(out) :: qdirs(:,:)
2176  real(dp),optional,intent(in) :: qlen
2177 
2178 !Local variables-------------------------------
2179  integer :: iq
2180  real(dp) :: qred2cart(3,3), qcart2red(3,3)
2181 
2182 ! *************************************************************************
2183 
2184  qred2cart = two_pi * cryst%gprimd
2185  qcart2red = qred2cart
2186  call matrginv(qcart2red, 3, 3)
2187 
2188  nq = 6
2189  ABI_MALLOC(qdirs, (3, nq))
2190  qdirs(:,1) = [one, zero, zero]  ! (100)
2191  qdirs(:,2) = [zero, one, zero]  ! (010)
2192  qdirs(:,3) = [zero, zero, one]  ! (001)
2193  qdirs(:,4) = matmul(qcart2red, [one, zero, zero]) ! (x)
2194  qdirs(:,5) = matmul(qcart2red, [zero, one, zero]) ! (y)
2195  qdirs(:,6) = matmul(qcart2red, [zero, zero, one]) ! (z)
2196 
2197  ! normalization
2198  do iq=1,nq
2199    qdirs(:,iq) = qdirs(:,iq) / normv(qdirs(:,iq), cryst%gmet, "G")
2200  end do
2201 
2202  if (present(qlen)) qdirs = qlen * qdirs
2203 
2204 end subroutine get_redcart_qdirs

m_crystal/idx_spatial_inversion [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

  idx_spatial_inversion

FUNCTION

  Return the index of the spatial inversion, 0 if not present

SOURCE

1193 pure function idx_spatial_inversion(Cryst) result(inv_idx)
1194 
1195 !Arguments ------------------------------------
1196 !scalars
1197  integer :: inv_idx
1198  class(crystal_t),intent(in) :: Cryst
1199 
1200 !Local variables-------------------------------
1201 !scalars
1202  integer :: isym
1203 
1204 ! *************************************************************************
1205 
1206  inv_idx=0
1207  do isym=1,cryst%nsym
1208    if (all(cryst%symrel(:,:,isym) == inversion_3d)) then
1209     inv_idx=isym; return
1210    end if
1211  end do
1212 
1213 end function idx_spatial_inversion

m_crystal/isalchemical [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

  isalchemical

FUNCTION

  Returns .TRUE. if we are using alchemical pseudopotentials

SOURCE

1252 pure logical function isalchemical(Cryst) result(ans)
1253 
1254 !Arguments ------------------------------------
1255  class(crystal_t),intent(in) :: Cryst
1256 
1257 ! *************************************************************************
1258 
1259  ans = (Cryst%npsp /= Cryst%ntypat)
1260 
1261 end function isalchemical

m_crystal/isymmorphic [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

  isymmorphic

FUNCTION

  Returns .TRUE. if space group is symmorphic, i.e. all fractional translations are zero.

SOURCE

1227 pure function isymmorphic(Cryst) result(ans)
1228 
1229 !Arguments ------------------------------------
1230 !scalars
1231  logical :: ans
1232  class(crystal_t),intent(in) :: Cryst
1233 
1234 ! *************************************************************************
1235 
1236  ans = ALL(ABS(Cryst%tnons) < tol6)
1237 
1238 end function isymmorphic

m_crystal/prt_cif [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

 prt_cif

FUNCTION

   print out CIF format file

INPUTS

OUTPUT

NOTES

SOURCE

1752 subroutine prt_cif(brvltt, ciffname, natom, nsym, ntypat, rprimd, &
1753                    spgaxor, spgroup, spgorig, symrel, tnon, typat, xred, znucl)
1754 
1755 !Arguments ------------------------------------
1756 !scalars
1757  integer,intent(in) :: natom, ntypat, nsym
1758  integer, intent(in) :: brvltt, spgaxor, spgroup, spgorig
1759 !arrays
1760  integer, intent(in) :: typat(natom)
1761  integer, intent(in) :: symrel(3,3,nsym)
1762  character(len=*), intent(in) :: ciffname
1763  real(dp), intent(in) :: tnon(3,nsym)
1764  real(dp), intent(in) :: rprimd(3,3)
1765  real(dp), intent(in) :: xred(3,natom)
1766  real(dp), intent(in) :: znucl(ntypat)
1767 
1768 !Local variables -------------------------------
1769 !scalars
1770  integer :: unitcif, iatom, isym, sporder, itypat, nat_this_type
1771  real(dp) :: ucvol
1772  type(atomdata_t) :: atom
1773 !arrays
1774  character(len=80) :: tmpstring
1775  character(len=1) :: brvsb
1776  character(len=15) :: intsb,ptintsb,ptschsb,schsb
1777  character(len=35) :: intsbl
1778  character(len=10) :: str_nat_type
1779  character(len=100) :: chemformula
1780  character(len=500) :: msg
1781  real(dp) :: angle(3), gprimd(3,3), rmet(3,3), gmet(3,3)
1782 
1783 !*************************************************************************
1784 
1785  ! open file in append mode xlf and other compilers refuse append mode
1786  if (open_file(ciffname,msg,newunit=unitcif) /=0) then
1787    ABI_WARNING(msg)
1788    return
1789  end if
1790 
1791  ! print title for dataset
1792  write (unitcif,'(a)') 'data_set'
1793 
1794  ! print cell parameters a,b,c, angles, volume
1795  call metric(gmet,gprimd,-1,rmet,rprimd,ucvol)
1796  angle(1)=acos(rmet(2,3)/sqrt(rmet(2,2)*rmet(3,3)))/two_pi*360.0_dp
1797  angle(2)=acos(rmet(1,3)/sqrt(rmet(1,1)*rmet(3,3)))/two_pi*360.0_dp
1798  angle(3)=acos(rmet(1,2)/sqrt(rmet(1,1)*rmet(2,2)))/two_pi*360.0_dp
1799 
1800  write (unitcif,'(a,E20.10)') '_cell_length_a                     ', sqrt(rmet(1,1))*Bohr_Ang
1801  write (unitcif,'(a,E20.10)') '_cell_length_b                     ', sqrt(rmet(2,2))*Bohr_Ang
1802  write (unitcif,'(a,E20.10)') '_cell_length_c                     ', sqrt(rmet(3,3))*Bohr_Ang
1803  write (unitcif,'(a,E20.10)') '_cell_angle_alpha                  ', angle(1)
1804  write (unitcif,'(a,E20.10)') '_cell_angle_beta                   ', angle(2)
1805  write (unitcif,'(a,E20.10)') '_cell_angle_gamma                  ', angle(3)
1806  write (unitcif,'(a,E20.10)') '_cell_volume                       ', ucvol*(Bohr_Ang)**3
1807 
1808  ! print reduced positions
1809  write (unitcif,'(a)') 'loop_'
1810  write (unitcif,'(a,E20.10)') '  _atom_site_label                   '
1811  write (unitcif,'(a,E20.10)') '  _atom_site_fract_x                 '
1812  write (unitcif,'(a,E20.10)') '  _atom_site_fract_y                 '
1813  write (unitcif,'(a,E20.10)') '  _atom_site_fract_z                 '
1814  do iatom = 1, natom
1815    call atomdata_from_znucl(atom,znucl(typat(iatom)))
1816    write (unitcif,'(2a,3E20.10)') '  ', atom%symbol, xred(:,iatom)
1817  end do
1818 
1819 !other specs in CIF dictionary which may be useful:
1820 !GEOM_BOND GEOM_ANGLE GEOM_TORSION
1821 
1822  ! print chemical composition in simplest form
1823  chemformula = "'"
1824  do itypat = 1, ntypat
1825    nat_this_type = 0
1826    do iatom = 1, natom
1827      if (typat(iatom) == itypat) nat_this_type = nat_this_type+1
1828    end do
1829    call atomdata_from_znucl(atom,znucl(itypat))
1830    call int2char10(nat_this_type, str_nat_type)
1831    chemformula = trim(chemformula) // atom%symbol // trim(str_nat_type) // "  "
1832  end do
1833  chemformula = trim(chemformula) // "'"
1834  write (unitcif,'(2a)') '_chemical_formula_analytical              ', chemformula
1835 
1836  !FIXME: check that brvltt is correctly used here - is it equal to bravais(1) in the invars routines?
1837  if (brvltt==1) then
1838    write (unitcif,'(a)') '_symmetry_cell_setting             triclinic'
1839  else if(brvltt==2)then
1840    write (unitcif,'(a)') '_symmetry_cell_setting             monoclinic'
1841  else if(brvltt==3)then
1842    write (unitcif,'(a)') '_symmetry_cell_setting             orthorhombic'
1843  else if(brvltt==4)then
1844    write (unitcif,'(a)') '_symmetry_cell_setting             tetragonal'
1845  else if(brvltt==5)then
1846    write (unitcif,'(a)') '_symmetry_cell_setting             rhombohedral'
1847  else if(brvltt==6)then
1848    write (unitcif,'(a)') '_symmetry_cell_setting             hexagonal'
1849  else if(brvltt==7)then
1850    write (unitcif,'(a)') '_symmetry_cell_setting             cubic'
1851  end if
1852 
1853  call spgdata(brvsb,intsb,intsbl,ptintsb,ptschsb,schsb,spgaxor,spgroup,sporder,spgorig)
1854 
1855  ! print symmetry operations
1856  write (unitcif,'(a,I6)') "_symmetry_Int_Tables_number          ", spgroup
1857  write (unitcif,'(5a)') "_symmetry_space_group_name_H-M        '", brvsb, " ", trim(intsb), "'"
1858  write (unitcif,'(a)') ''
1859  write (unitcif,'(a)') 'loop_'
1860  write (unitcif,'(a)') '  _symmetry_equiv_pos_as_xyz           '
1861  do isym = 1, nsym
1862    call  symrel2string(symrel(:,:,isym), tnon(:,isym), tmpstring)
1863    write (unitcif,'(2a)') '  ', trim(tmpstring)
1864  end do
1865 
1866  close(unitcif)
1867 
1868 end subroutine prt_cif

m_crystal/prtposcar [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

  prtposcar

FUNCTION

  output VASP style POSCAR and FORCES files for use with frozen phonon codes, like
  PHON from Dario Alfe' or frophon
  IMPORTANT: the order of atoms is fixed such that typat is re-grouped.
  First typat=1 then typat=2, etc...
  Only master should call this routine in MPI-mode.

INPUTS

  fcart = forces on atoms in cartesian coordinates
  natom = number of atoms
  ntypat = number of types of atoms
  rprimd = lattice vectors for the primitive cell
  typat = type for each of the natom atoms
  ucvol = unit cell volume
  xred = reduced positions of the atoms
  znucl = nuclear charge of each atomic type

 OUTPUTS
   Only files written

SOURCE

1958 subroutine prtposcar(fcart, fnameradix, natom, ntypat, rprimd, typat, ucvol, xred, znucl)
1959 
1960 !Arguments ------------------------------------
1961 !scalars
1962  integer, intent(in) :: natom, ntypat
1963  real(dp), intent(in) :: ucvol
1964 !arrays
1965  integer, intent(in) :: typat(natom)
1966  real(dp), intent(in) :: fcart(3,natom)
1967  real(dp), intent(in) :: rprimd(3,3)
1968  real(dp), intent(in) :: xred(3,natom)
1969  real(dp), intent(in) :: znucl(ntypat)
1970  character(len=fnlen), intent(in) :: fnameradix
1971 
1972 !Local variables-------------------------------
1973 !scalars
1974  integer :: iatom, itypat, iout
1975  type(atomdata_t) :: atom
1976 ! arrays
1977  integer :: natoms_this_type(ntypat)
1978  character(len=2) :: symbol
1979  character(len=7) :: natoms_this_type_str
1980  character(len=100) :: chem_formula, natoms_all_types
1981  character(len=500) :: msg
1982 
1983 !************************************************************************
1984 
1985  ! Output POSCAR file for positions, atom types etc
1986  if (open_file(trim(fnameradix)//"_POSCAR", msg, newunit=iout) /= 0) then
1987    ABI_ERROR(msg)
1988  end if
1989 
1990  natoms_this_type = 0
1991  do itypat=1,ntypat
1992    do iatom=1,natom
1993      if (typat(iatom) == itypat) natoms_this_type(itypat) = natoms_this_type(itypat) + 1
1994    end do
1995  end do
1996 
1997  chem_formula = ""
1998  do itypat=1, ntypat
1999    call atomdata_from_znucl(atom, znucl(itypat))
2000    symbol = atom%symbol
2001    if (natoms_this_type(itypat) < 10) then
2002      write (natoms_this_type_str, '(I1)') natoms_this_type(itypat)
2003    else if (natoms_this_type(itypat) < 100) then
2004      write (natoms_this_type_str, '(I2)') natoms_this_type(itypat)
2005    else if (natoms_this_type(itypat) < 1000) then
2006      write (natoms_this_type_str, '(I3)') natoms_this_type(itypat)
2007    end if
2008    chem_formula = trim(chem_formula) // symbol // trim(natoms_this_type_str)
2009  end do
2010 
2011  write (iout,'(2a)') "ABINIT generated POSCAR file. Title string - should be chemical formula... ",trim(chem_formula)
2012 
2013  write (iout,'(E24.14)') -ucvol*Bohr_Ang*Bohr_Ang*Bohr_Ang
2014  write (iout,'(3E24.14,1x)') Bohr_Ang*rprimd(:,1) ! (angstr? bohr?)
2015  write (iout,'(3E24.14,1x)') Bohr_Ang*rprimd(:,2)
2016  write (iout,'(3E24.14,1x)') Bohr_Ang*rprimd(:,3)
2017 
2018  natoms_all_types = "   "
2019  do itypat=1, ntypat
2020    write (natoms_this_type_str, '(I7)') natoms_this_type(itypat)
2021    natoms_all_types = trim(natoms_all_types) // "   " // trim(natoms_this_type_str)
2022  end do
2023 
2024  write (iout,'(a)') trim(natoms_all_types)
2025  write (iout,'(a)') "Direct"
2026 
2027  do itypat=1, ntypat
2028    do iatom=1,natom
2029      if (typat(iatom) /= itypat) cycle
2030      write (iout,'(3(E24.14,1x))') xred(:,iatom)
2031    end do
2032  end do
2033  close (iout)
2034 
2035  ! output FORCES file for forces in same order as positions above
2036  if (open_file(trim(fnameradix)//"_FORCES", msg, newunit=iout) /= 0 ) then
2037    ABI_ERROR(msg)
2038  end if
2039 
2040  !ndisplacements
2041  !iatom_displaced displacement_red_coord(3)
2042  !forces_cart_ev_Angstr(3)
2043  !...
2044  !<repeat for other displaced atoms>
2045  write (iout,'(I7)') 1
2046  write (iout,'(a)') '1 0 0 0        ! TO BE FILLED IN '
2047  do itypat=1, ntypat
2048    do iatom=1,natom
2049      if (typat(iatom) /= itypat) cycle
2050      write (iout,'(3(E24.14,1x))') Ha_eV/Bohr_Ang*fcart(:,iatom)
2051    end do
2052  end do
2053 
2054  close (iout)
2055 
2056 end subroutine prtposcar

m_crystal/symbol_iatom [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

  symbol_iatom

FUNCTION

  Return the atomic symbol from the iatom index

SOURCE

1331 function symbol_iatom(crystal, iatom) result(symbol)
1332 
1333 !Arguments ------------------------------------
1334 !scalars
1335  integer,intent(in) :: iatom
1336  character(len=2) :: symbol
1337  class(crystal_t),intent(in) :: crystal
1338 
1339 ! *************************************************************************
1340 
1341  symbol = crystal%symbol_type(crystal%typat(iatom))
1342 
1343 end function symbol_iatom

m_crystal/symbol_type [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

  symbol_type

FUNCTION

  Return the atomic symbol from the itypat index

SOURCE

1300 function symbol_type(crystal, itypat) result(symbol)
1301 
1302 !Arguments ------------------------------------
1303 !scalars
1304  integer,intent(in) :: itypat
1305  character(len=2) :: symbol
1306  class(crystal_t),intent(in) :: crystal
1307 
1308 !Local variables-------------------------------
1309 !scalars
1310  type(atomdata_t) :: atom
1311 
1312 ! *************************************************************************
1313 
1314  atom = crystal%adata_type(itypat)
1315  symbol = atom%symbol
1316 
1317 end function symbol_type

m_crystal/symbols_crystal [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

 symbols_crystal

FUNCTION

 Return a array with the symbol of each atoms with indexation e.g.
 ["Sr","Ru","O1","O2","O3"]

INPUTS

 natom = number of atoms
 ntypat = number of typat
 npsp =  number of pseudopotentials
 znucl = Nuclear charge for each type of pseudopotential

OUTPUT

 symbols = array with the symbol of each atoms

SOURCE

1139 subroutine symbols_crystal(natom,ntypat,npsp,symbols,typat,znucl)
1140 
1141 !Arguments ------------------------------------
1142 !scalars
1143  integer,intent(in) :: natom,ntypat,npsp
1144 !arrays
1145  real(dp),intent(in):: znucl(npsp)
1146  integer,intent(in) :: typat(natom)
1147  character(len=5),intent(out) :: symbols(natom)
1148  character(len=3) :: powerchar
1149 
1150 !Local variables-------------------------------
1151 !scalar
1152  integer :: ia,ii,itypat,jj
1153 ! *************************************************************************
1154 
1155  ! Fill the symbols array
1156  do ia=1,natom
1157    symbols(ia) = adjustl(znucl2symbol(znucl(typat(ia))))
1158  end do
1159  itypat = 0
1160  do itypat =1,ntypat
1161    ii = 0
1162    do ia=1,natom
1163      if(typat(ia)==itypat) then
1164        ii = ii + 1
1165      end if
1166    end do
1167    if(ii>1)then
1168      jj=1
1169      do ia=1,natom
1170        if(typat(ia)==itypat) then
1171          write(powerchar,'(I0)') jj
1172          symbols(ia) = trim(symbols(ia))//trim(powerchar)
1173          jj=jj+1
1174        end if
1175      end do
1176    end if
1177  end do
1178 
1179 end subroutine symbols_crystal

m_crystal/symrel2string [ Functions ]

[ Top ] [ m_crystal ] [ Functions ]

NAME

 symrel2string

FUNCTION

INPUTS

OUTPUT

NOTES

SOURCE

1885 subroutine symrel2string(symrel1, tnon, string)
1886 
1887 !Arguments ------------------------------------
1888 !scalars
1889  integer, intent(in) :: symrel1(3,3)
1890  real(dp), intent(in) :: tnon(3)
1891  character(len=80), intent(out) :: string
1892 
1893 !Local variables -------------------------
1894 !scalars
1895  integer :: i1,i2
1896  character(len=1) :: xyz(3)
1897 
1898 ! *********************************************************************
1899 
1900  xyz(1) = 'x'
1901  xyz(2) = 'y'
1902  xyz(3) = 'z'
1903 
1904  string = ''
1905  do i1=1,3
1906    if (abs(tnon(i1)) > tol10) then
1907      ! find fraction 1/n for tnon, otherwise do not know what to print
1908      if (abs(one-two*tnon(i1)) < tol10) string = trim(string)//'1/2'
1909      if (abs(one+two*tnon(i1)) < tol10) string = trim(string)//'-1/2'
1910 
1911      if (abs(one-three*tnon(i1)) < tol10) string = trim(string)//'1/3'
1912      if (abs(one+three*tnon(i1)) < tol10) string = trim(string)//'-1/3'
1913      if (abs(two-three*tnon(i1)) < tol10) string = trim(string)//'2/3'
1914      if (abs(two+three*tnon(i1)) < tol10) string = trim(string)//'-2/3'
1915 
1916      if (abs(one-six*tnon(i1)) < tol10) string = trim(string)//'1/6'
1917      if (abs(one+six*tnon(i1)) < tol10) string = trim(string)//'-1/6'
1918      if (abs(five-six*tnon(i1)) < tol10) string = trim(string)//'5/6'
1919      if (abs(five+six*tnon(i1)) < tol10) string = trim(string)//'-5/6'
1920    end if
1921    do i2=1,3
1922      ! FIXME: check if this is correct ordering for symrel(i1,i2) looks ok
1923      if (symrel1(i1,i2) == 1)  string = trim(string)//'+'//xyz(i2)
1924      if (symrel1(i1,i2) == -1) string = trim(string)//'-'//xyz(i2)
1925    end do
1926    if (i1 /= 3) string = trim(string)//','
1927  end do
1928 
1929 end subroutine symrel2string