TABLE OF CONTENTS
- ABINIT/m_crystal
- m_crystal/adata_type
- m_crystal/crystal_bcast
- m_crystal/crystal_compare
- m_crystal/crystal_compute_geometry
- m_crystal/crystal_compute_sym
- m_crystal/crystal_copy
- m_crystal/crystal_free
- m_crystal/crystal_index_atoms
- m_crystal/crystal_init
- m_crystal/crystal_malloc
- m_crystal/crystal_ncread
- m_crystal/crystal_ncwrite
- m_crystal/crystal_ncwrite_path
- m_crystal/crystal_point_group
- m_crystal/crystal_print
- m_crystal/crystal_print_abivars
- m_crystal/crystal_symmetrize_cart_tens33
- m_crystal/crystal_symmetrize_cart_vec3
- m_crystal/crystal_t
- m_crystal/crystal_without_symmetries
- m_crystal/get_recart_qdirs
- m_crystal/idx_spatial_inversion
- m_crystal/isalchemical
- m_crystal/isymmorphic
- m_crystal/prt_cif
- m_crystal/prtposcar
- m_crystal/symbol_iatom
- m_crystal/symbol_type
- m_crystal/symbols_crystal
- m_crystal/symrel2string
ABINIT/m_crystal [ 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