TABLE OF CONTENTS


ABINIT/m_ab7_symmetry [ Modules ]

[ Top ] [ Modules ]

NAME

 m_ab7_symmetry

FUNCTION

COPYRIGHT

  Copyright (C) 2008-2018 ABINIT group (DC)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

PARENTS

SOURCE

  18 #if defined HAVE_CONFIG_H
  19 #include "config.h"
  20 #endif
  21 
  22 #include "abi_common.h"
  23 
  24 module m_ab7_symmetry
  25 
  26   use defs_basis
  27   use m_abicore
  28 
  29   use m_symtk,     only : mati3inv, mati3det, symatm, symcharac
  30   use m_symfind,   only : symfind, symanal, symlatt
  31   use m_geometry,  only : metric
  32   use m_spgdata,   only : spgdata
  33 
  34   implicit none
  35 
  36   private
  37 
  38   integer, parameter, public :: AB7_MAX_SYMMETRIES = 384
  39 
  40   type, public :: symmetry_type
  41      ! The input characteristics
  42      real(dp) :: tolsym
  43      real(dp) :: rprimd(3,3), gprimd(3,3), rmet(3,3)
  44      integer :: nAtoms
  45      integer, pointer :: typeAt(:)
  46      real(dp), pointer :: xRed(:,:)
  47 
  48      logical :: withField
  49      real(dp) :: field(3)
  50 
  51      logical :: withJellium
  52 
  53      integer :: nzchempot
  54 
  55      integer :: withSpin
  56      real(dp), pointer :: spinAt(:,:)
  57 
  58      logical :: withSpinOrbit
  59 
  60      integer :: vacuum(3)
  61 
  62      ! The output characteristics
  63      ! The bravais parameters
  64      integer :: nBravSym
  65      integer :: bravais(11), bravSym(3, 3, AB7_MAX_SYMMETRIES)
  66      ! The symmetry matrices
  67      logical  :: auto
  68      integer  :: nSym
  69      integer, pointer  :: sym(:,:,:)
  70      real(dp), pointer :: transNon(:,:)
  71      integer, pointer  :: symAfm(:)
  72      ! Some additional information
  73      integer          :: multiplicity
  74      real(dp)         :: genAfm(3)
  75      integer          :: spaceGroup, pointGroupMagn
  76      integer, pointer :: indexingAtoms(:,:,:)
  77   end type symmetry_type
  78 
  79   ! We store here a list of symmetry objects to be able to
  80   ! call several symmetry operations on different objects.
  81   ! The simplest portable way to do it, is to create
  82   ! a list of Fortran structure and to use the list index
  83   ! as an identifier that can be given to the other languages.
  84   type, private :: symmetry_list
  85      integer                       :: id
  86      type(symmetry_list),  pointer :: next
  87      type(symmetry_type)           :: data
  88   end type symmetry_list
  89   type(symmetry_list), pointer :: my_symmetries
  90   integer :: n_symmetries = 0
  91 
  92   logical, private, parameter :: AB_DBG = .false.
  93 
  94   public :: symmetry_new
  95   public :: symmetry_free
  96   public :: symmetry_set_tolerance
  97   public :: symmetry_set_lattice
  98   public :: symmetry_set_structure
  99   public :: symmetry_set_collinear_spin
 100   public :: symmetry_set_spin
 101   public :: symmetry_set_spin_orbit
 102   public :: symmetry_set_field
 103   public :: symmetry_set_jellium
 104   public :: symmetry_set_periodicity
 105   public :: symmetry_set_n_sym
 106 
 107   public :: symmetry_get_from_id
 108   public :: symmetry_get_n_atoms
 109   public :: symmetry_get_n_sym
 110   public :: symmetry_get_multiplicity
 111   public :: symmetry_get_bravais
 112   public :: symmetry_get_matrices
 113   public :: symmetry_get_matrices_p
 114   public :: symmetry_get_group
 115   public :: symmetry_get_equivalent_atom
 116   public :: symmetry_get_type
 117 
 118 contains
 119 
 120   subroutine new_item(token)
 121 
 122 
 123 !This section has been created automatically by the script Abilint (TD).
 124 !Do not modify the following lines by hand.
 125 #undef ABI_FUNC
 126 #define ABI_FUNC 'new_item'
 127 !End of the abilint section
 128 
 129     type(symmetry_list), pointer :: token
 130 
 131     ! We allocate a new list token and prepend it.
 132     if (AB_DBG) write(std_err,*) "AB symmetry: create a new token."
 133 
 134     ! Init case, very first call.
 135     if (n_symmetries == 0) then
 136        nullify(my_symmetries)
 137     end if
 138 
 139     ! Normal treatment.
 140     n_symmetries = n_symmetries + 1
 141 
 142     ABI_DATATYPE_ALLOCATE(token,)
 143     token%id = n_symmetries
 144     call new_symmetry(token%data)
 145     token%next => my_symmetries
 146 
 147     my_symmetries => token
 148     if (AB_DBG) write(std_err,*) "AB symmetry: creation OK with id ", token%id
 149   end subroutine new_item
 150 
 151   subroutine free_item(token)
 152 
 153 
 154 !This section has been created automatically by the script Abilint (TD).
 155 !Do not modify the following lines by hand.
 156 #undef ABI_FUNC
 157 #define ABI_FUNC 'free_item'
 158 !End of the abilint section
 159 
 160     type(symmetry_list), pointer :: token
 161 
 162     type(symmetry_list), pointer :: tmp
 163 
 164     if (.not. associated(token)) then
 165        return
 166     end if
 167 
 168     call free_symmetry(token%data)
 169 
 170     if (AB_DBG) write(std_err,*) "AB symmetry: free request on token ", token%id
 171     ! We remove token from the list.
 172     if (my_symmetries%id == token%id) then
 173        my_symmetries => token%next
 174     else
 175        tmp => my_symmetries
 176        do
 177           if (.not.associated(tmp)) then
 178              return
 179           end if
 180           if (associated(tmp%next) .and. tmp%next%id == token%id) then
 181              exit
 182           end if
 183           tmp => tmp%next
 184        end do
 185        tmp%next => token%next
 186     end if
 187     ABI_DATATYPE_DEALLOCATE(token)
 188     if (AB_DBG) write(std_err,*) "AB symmetry: free done"
 189   end subroutine free_item
 190 
 191   subroutine get_item(token, id)
 192 
 193 
 194 !This section has been created automatically by the script Abilint (TD).
 195 !Do not modify the following lines by hand.
 196 #undef ABI_FUNC
 197 #define ABI_FUNC 'get_item'
 198 !End of the abilint section
 199 
 200     type(symmetry_list), pointer :: token
 201     integer, intent(in) :: id
 202 
 203     type(symmetry_list), pointer :: tmp
 204 
 205     if (AB_DBG) write(std_err,*) "AB symmetry: request list element ", id
 206     nullify(token)
 207 
 208     tmp => my_symmetries
 209     do
 210        if (.not. associated(tmp)) then
 211           exit
 212        end if
 213        if (tmp%id == id) then
 214           token => tmp
 215           return
 216        end if
 217        tmp => tmp%next
 218     end do
 219   end subroutine get_item
 220 
 221   subroutine symmetry_get_from_id(sym, id, errno)
 222 
 223 
 224 !This section has been created automatically by the script Abilint (TD).
 225 !Do not modify the following lines by hand.
 226 #undef ABI_FUNC
 227 #define ABI_FUNC 'symmetry_get_from_id'
 228 !End of the abilint section
 229 
 230     type(symmetry_type), pointer :: sym
 231     integer, intent(in) :: id
 232     integer, intent(out) :: errno
 233 
 234     type(symmetry_list), pointer :: token
 235 
 236     errno = AB7_NO_ERROR
 237     call get_item(token, id)
 238     if (associated(token)) then
 239        sym => token%data
 240        if (sym%nSym <= 0) then
 241           ! We do the computation of the matrix part.
 242           call compute_matrices(sym, errno)
 243        end if
 244     else
 245        errno = AB7_ERROR_OBJ
 246        nullify(sym)
 247     end if
 248   end subroutine symmetry_get_from_id
 249 
 250   subroutine new_symmetry(sym)
 251 
 252 
 253 !This section has been created automatically by the script Abilint (TD).
 254 !Do not modify the following lines by hand.
 255 #undef ABI_FUNC
 256 #define ABI_FUNC 'new_symmetry'
 257 !End of the abilint section
 258 
 259     type(symmetry_type), intent(out) :: sym
 260 
 261     if (AB_DBG) write(std_err,*) "AB symmetry: create a new symmetry object."
 262     nullify(sym%xRed)
 263     nullify(sym%spinAt)
 264     nullify(sym%typeAt)
 265     sym%tolsym   = tol8
 266     sym%auto     = .true.
 267     sym%nSym     = 0
 268     nullify(sym%sym)
 269     nullify(sym%symAfm)
 270     nullify(sym%transNon)
 271     sym%nBravSym = -1
 272     sym%withField   = .false.
 273     sym%withJellium = .false.
 274     sym%nzchempot = 0
 275     sym%withSpin = 1
 276     sym%withSpinOrbit = .false.
 277     sym%multiplicity = -1
 278     nullify(sym%indexingAtoms)
 279     sym%vacuum = 0
 280   end subroutine new_symmetry
 281 
 282   subroutine free_symmetry(sym)
 283 
 284 
 285 !This section has been created automatically by the script Abilint (TD).
 286 !Do not modify the following lines by hand.
 287 #undef ABI_FUNC
 288 #define ABI_FUNC 'free_symmetry'
 289 !End of the abilint section
 290 
 291     type(symmetry_type), intent(inout) :: sym
 292 
 293     if (AB_DBG) write(std_err,*) "AB symmetry: free a symmetry."
 294 
 295     if (associated(sym%xRed))  then
 296       ABI_DEALLOCATE(sym%xRed)
 297     end if
 298     if (associated(sym%spinAt))  then
 299       ABI_DEALLOCATE(sym%spinAt)
 300     end if
 301     if (associated(sym%typeAt))  then
 302       ABI_DEALLOCATE(sym%typeAt)
 303     end if
 304     if (associated(sym%indexingAtoms))  then
 305       ABI_DEALLOCATE(sym%indexingAtoms)
 306     end if
 307     if (associated(sym%sym))  then
 308       ABI_DEALLOCATE(sym%sym)
 309     end if
 310     if (associated(sym%symAfm))  then
 311       ABI_DEALLOCATE(sym%symAfm)
 312     end if
 313     if (associated(sym%transNon))  then
 314       ABI_DEALLOCATE(sym%transNon)
 315     end if
 316   end subroutine free_symmetry
 317 
 318 
 319 
 320 
 321 
 322   subroutine symmetry_new(id)
 323 
 324 
 325 !This section has been created automatically by the script Abilint (TD).
 326 !Do not modify the following lines by hand.
 327 #undef ABI_FUNC
 328 #define ABI_FUNC 'symmetry_new'
 329 !End of the abilint section
 330 
 331     integer, intent(out) :: id
 332 
 333     type(symmetry_list), pointer :: token
 334 
 335     if (AB_DBG) write(std_err,*) "AB symmetry: call new symmetry."
 336     call new_item(token)
 337     id = token%id
 338   end subroutine symmetry_new
 339 
 340   subroutine symmetry_free(id)
 341 
 342 
 343 !This section has been created automatically by the script Abilint (TD).
 344 !Do not modify the following lines by hand.
 345 #undef ABI_FUNC
 346 #define ABI_FUNC 'symmetry_free'
 347 !End of the abilint section
 348 
 349     integer, intent(in) :: id
 350 
 351     type(symmetry_list), pointer :: token
 352 
 353     if (AB_DBG) write(std_err,*) "AB symmetry: call free symmetry."
 354 
 355     call get_item(token, id)
 356     if (associated(token)) then
 357       call free_item(token)
 358     end if
 359   end subroutine symmetry_free
 360 
 361   subroutine symmetry_set_tolerance(id, tolsym, errno)
 362 
 363 
 364 !This section has been created automatically by the script Abilint (TD).
 365 !Do not modify the following lines by hand.
 366 #undef ABI_FUNC
 367 #define ABI_FUNC 'symmetry_set_tolerance'
 368 !End of the abilint section
 369 
 370     integer, intent(in) :: id
 371     real(dp), intent(in) :: tolsym
 372     integer, intent(out) :: errno
 373 
 374     type(symmetry_list), pointer :: token
 375 
 376     if (AB_DBG) write(std_err,*) "AB symmetry: call set tolerance."
 377 
 378     errno = AB7_NO_ERROR
 379     call get_item(token, id)
 380     if (.not. associated(token)) then
 381        errno = AB7_ERROR_OBJ
 382        return
 383     end if
 384 
 385     token%data%tolsym = tolsym
 386 
 387     ! We unset all the computed symmetries
 388     token%data%nBravSym = -1
 389     if (token%data%auto) then
 390        token%data%nSym  = 0
 391     end if
 392   end subroutine symmetry_set_tolerance
 393 
 394   subroutine symmetry_set_lattice(id, rprimd, errno)
 395 
 396 
 397 !This section has been created automatically by the script Abilint (TD).
 398 !Do not modify the following lines by hand.
 399 #undef ABI_FUNC
 400 #define ABI_FUNC 'symmetry_set_lattice'
 401 !End of the abilint section
 402 
 403     integer, intent(in) :: id
 404     real(dp), intent(in) :: rprimd(3,3)
 405     integer, intent(out) :: errno
 406 
 407     type(symmetry_list), pointer :: token
 408     real(dp) :: ucvol
 409     real(dp) :: gmet(3,3)
 410 
 411     if (AB_DBG) write(std_err,*) "AB symmetry: call set lattice."
 412     if (AB_DBG) write(std_err, "(A,3F12.6,A)") "  (", rprimd(:,1), ")"
 413     if (AB_DBG) write(std_err, "(A,3F12.6,A)") "  (", rprimd(:,2), ")"
 414     if (AB_DBG) write(std_err, "(A,3F12.6,A)") "  (", rprimd(:,3), ")"
 415 
 416     errno = AB7_NO_ERROR
 417     call get_item(token, id)
 418     if (.not. associated(token)) then
 419        errno = AB7_ERROR_OBJ
 420        return
 421     end if
 422 
 423     token%data%rprimd = rprimd
 424     call metric(gmet, token%data%gprimd, -1, token%data%rmet, rprimd, ucvol)
 425 
 426     ! We unset all the computed symmetries
 427     token%data%nBravSym = -1
 428     if (token%data%auto) then
 429        token%data%nSym  = 0
 430     end if
 431   end subroutine symmetry_set_lattice
 432 
 433   subroutine symmetry_set_structure(id, nAtoms, typeAt, xRed, errno)
 434 
 435 
 436 !This section has been created automatically by the script Abilint (TD).
 437 !Do not modify the following lines by hand.
 438 #undef ABI_FUNC
 439 #define ABI_FUNC 'symmetry_set_structure'
 440 !End of the abilint section
 441 
 442     integer, intent(in) :: id
 443     integer, intent(in) :: nAtoms
 444     integer, intent(in) :: typeAt(nAtoms)
 445     real(dp), intent(in) :: xRed(3,nAtoms)
 446     integer, intent(out) :: errno
 447 
 448     type(symmetry_list), pointer :: token
 449     integer :: i
 450 
 451     if (AB_DBG) write(std_err,*) "AB symmetry: call set structure."
 452     if (AB_DBG) write(std_err, "(A,I3,A)") "  ", nAtoms, " atoms"
 453     if (AB_DBG) then
 454        do i = 1, nAtoms, 1
 455           if (AB_DBG) write(std_err,"(A,3F12.6,I3)") "  ", xRed(:, i), typeAt(i)
 456        end do
 457     end if
 458 
 459     errno = AB7_NO_ERROR
 460     call get_item(token, id)
 461     if (.not. associated(token)) then
 462        errno = AB7_ERROR_OBJ
 463        return
 464     end if
 465 
 466     token%data%nAtoms =  nAtoms
 467     ABI_ALLOCATE(token%data%typeAt,(nAtoms))
 468     token%data%typeAt = typeAt
 469     ABI_ALLOCATE(token%data%xRed,(3, nAtoms))
 470     token%data%xRed   = xRed
 471 
 472     ! We unset only the symmetries
 473     if (token%data%auto) then
 474        token%data%nSym = 0
 475     end if
 476     if (associated(token%data%indexingAtoms)) then
 477       ABI_DEALLOCATE(token%data%indexingAtoms)
 478     end if
 479   end subroutine symmetry_set_structure
 480 
 481   subroutine symmetry_set_spin(id, nAtoms, spinAt, errno)
 482 
 483 
 484 !This section has been created automatically by the script Abilint (TD).
 485 !Do not modify the following lines by hand.
 486 #undef ABI_FUNC
 487 #define ABI_FUNC 'symmetry_set_spin'
 488 !End of the abilint section
 489 
 490     integer, intent(in) :: id
 491     integer, intent(in) :: nAtoms
 492     real(dp), intent(in) :: spinAt(3,nAtoms)
 493     integer, intent(out) :: errno
 494 
 495     type(symmetry_list), pointer :: token
 496     integer :: i
 497 
 498     if (AB_DBG) write(std_err,*) "AB symmetry: call set spin."
 499     if (AB_DBG) then
 500        do i = 1, nAtoms, 1
 501           if (AB_DBG) write(std_err,"(A,3F12.6)") "  ", spinAt(:, i)
 502        end do
 503     end if
 504 
 505     errno = AB7_NO_ERROR
 506     call get_item(token, id)
 507     if (.not. associated(token)) then
 508        errno = AB7_ERROR_OBJ
 509        return
 510     end if
 511     if (token%data%nAtoms /= nAtoms) then
 512        errno = AB7_ERROR_ARG
 513        return
 514     end if
 515 
 516     token%data%withSpin = 4
 517     ABI_ALLOCATE(token%data%spinAt,(3, nAtoms))
 518     token%data%spinAt = spinAt
 519 
 520     ! We unset only the symmetries
 521     if (token%data%auto) then
 522        token%data%nSym  = 0
 523     end if
 524   end subroutine symmetry_set_spin
 525 
 526   subroutine symmetry_set_collinear_spin(id, nAtoms, spinAt, errno)
 527 
 528 
 529 !This section has been created automatically by the script Abilint (TD).
 530 !Do not modify the following lines by hand.
 531 #undef ABI_FUNC
 532 #define ABI_FUNC 'symmetry_set_collinear_spin'
 533 !End of the abilint section
 534 
 535     integer, intent(in) :: id
 536     integer, intent(in) :: nAtoms
 537     integer, intent(in) :: spinAt(nAtoms)
 538     integer, intent(out) :: errno
 539 
 540     type(symmetry_list), pointer :: token
 541     integer :: i
 542 
 543     if (AB_DBG) write(std_err,*) "AB symmetry: call set collinear spin."
 544     if (AB_DBG) then
 545        do i = 1, nAtoms, 1
 546           if (AB_DBG) write(std_err,"(A,I3)") "  ", spinAt(i)
 547        end do
 548     end if
 549 
 550     errno = AB7_NO_ERROR
 551     call get_item(token, id)
 552     if (.not. associated(token)) then
 553        errno = AB7_ERROR_OBJ
 554        return
 555     end if
 556     if (token%data%nAtoms /= nAtoms) then
 557        errno = AB7_ERROR_ARG
 558        return
 559     end if
 560 
 561     token%data%withSpin = 2
 562     ABI_ALLOCATE(token%data%spinAt,(3,nAtoms))
 563     token%data%spinAt = 0._dp
 564     token%data%spinAt(3, :) = real(spinAt)
 565 
 566     ! We unset only the symmetries
 567     if (token%data%auto) then
 568        token%data%nSym  = 0
 569     end if
 570   end subroutine symmetry_set_collinear_spin
 571 
 572   subroutine symmetry_set_spin_orbit(id, withSpinOrbit, errno)
 573 
 574 
 575 !This section has been created automatically by the script Abilint (TD).
 576 !Do not modify the following lines by hand.
 577 #undef ABI_FUNC
 578 #define ABI_FUNC 'symmetry_set_spin_orbit'
 579 !End of the abilint section
 580 
 581     integer, intent(in) :: id
 582     logical, intent(in) :: withSpinOrbit
 583     integer, intent(out) :: errno
 584 
 585     type(symmetry_list), pointer :: token
 586 
 587     if (AB_DBG) write(std_err,*) "AB symmetry: call set spin orbit."
 588 
 589     errno = AB7_NO_ERROR
 590     call get_item(token, id)
 591     if (.not. associated(token)) then
 592        errno = AB7_ERROR_OBJ
 593        return
 594     end if
 595 
 596     token%data%withSpinOrbit = withSpinOrbit
 597 
 598     ! We unset only the symmetries
 599     if (token%data%auto) then
 600        token%data%nSym  = 0
 601     end if
 602   end subroutine symmetry_set_spin_orbit
 603 
 604   subroutine symmetry_set_field(id, field, errno)
 605 
 606 
 607 !This section has been created automatically by the script Abilint (TD).
 608 !Do not modify the following lines by hand.
 609 #undef ABI_FUNC
 610 #define ABI_FUNC 'symmetry_set_field'
 611 !End of the abilint section
 612 
 613     integer, intent(in) :: id
 614     real(dp), intent(in) :: field(3)
 615     integer, intent(out) :: errno
 616 
 617     type(symmetry_list), pointer :: token
 618 
 619     if (AB_DBG) write(std_err,*) "AB symmetry: call set field."
 620 
 621     errno = AB7_NO_ERROR
 622     call get_item(token, id)
 623     if (.not. associated(token)) then
 624        errno = AB7_ERROR_OBJ
 625        return
 626     end if
 627 
 628     token%data%withField = .true.
 629     token%data%field = field
 630 
 631     ! We unset all the computed symmetries
 632     token%data%nBravSym = -1
 633     if (token%data%auto) then
 634        token%data%nSym  = 0
 635     end if
 636   end subroutine symmetry_set_field
 637 
 638   subroutine symmetry_set_jellium(id, jellium, errno)
 639 
 640 
 641 !This section has been created automatically by the script Abilint (TD).
 642 !Do not modify the following lines by hand.
 643 #undef ABI_FUNC
 644 #define ABI_FUNC 'symmetry_set_jellium'
 645 !End of the abilint section
 646 
 647     integer, intent(in) :: id
 648     logical, intent(in) :: jellium
 649     integer, intent(out) :: errno
 650 
 651     type(symmetry_list), pointer :: token
 652 
 653     if (AB_DBG) write(std_err,*) "AB symmetry: call set jellium."
 654 
 655     errno = AB7_NO_ERROR
 656     call get_item(token, id)
 657     if (.not. associated(token)) then
 658        errno = AB7_ERROR_OBJ
 659        return
 660     end if
 661 
 662     token%data%withJellium = jellium
 663 
 664     ! We unset only the symmetries
 665     if (token%data%auto) then
 666        token%data%nSym  = 0
 667     end if
 668   end subroutine symmetry_set_jellium
 669 
 670 
 671   subroutine symmetry_set_nzchempot(id, nzchempot, errno)
 672 
 673 
 674 !This section has been created automatically by the script Abilint (TD).
 675 !Do not modify the following lines by hand.
 676 #undef ABI_FUNC
 677 #define ABI_FUNC 'symmetry_set_nzchempot'
 678 !End of the abilint section
 679 
 680     integer, intent(in) :: id
 681     integer, intent(in) :: nzchempot
 682     integer, intent(out) :: errno
 683 
 684     type(symmetry_list), pointer :: token
 685 
 686     if (AB_DBG) write(std_err,*) "AB symmetry: call set nzchempot."
 687 
 688     errno = AB7_NO_ERROR
 689     call get_item(token, id)
 690     if (.not. associated(token)) then
 691        errno = AB7_ERROR_OBJ
 692        return
 693     end if
 694 
 695     token%data%nzchempot = nzchempot
 696 
 697     ! We unset only the symmetries
 698     if (token%data%auto) then
 699        token%data%nSym  = 0
 700     end if
 701   end subroutine symmetry_set_nzchempot
 702 
 703   subroutine symmetry_set_periodicity(id, periodic, errno)
 704 
 705 
 706 !This section has been created automatically by the script Abilint (TD).
 707 !Do not modify the following lines by hand.
 708 #undef ABI_FUNC
 709 #define ABI_FUNC 'symmetry_set_periodicity'
 710 !End of the abilint section
 711 
 712     integer, intent(in) :: id
 713     logical, intent(in) :: periodic(3)
 714     integer, intent(out) :: errno
 715 
 716     type(symmetry_list), pointer :: token
 717 
 718     if (AB_DBG) write(std_err,*) "AB symmetry: call set periodicity."
 719     if (AB_DBG) write(std_err, "(A,3L1,A)") "  (", periodic, ")"
 720 
 721     errno = AB7_NO_ERROR
 722     call get_item(token, id)
 723     if (.not. associated(token)) then
 724        errno = AB7_ERROR_OBJ
 725        return
 726     end if
 727 
 728     token%data%vacuum = 0
 729     if (.not. periodic(1)) token%data%vacuum(1) = 1
 730     if (.not. periodic(2)) token%data%vacuum(2) = 1
 731     if (.not. periodic(3)) token%data%vacuum(3) = 1
 732   end subroutine symmetry_set_periodicity
 733 
 734 
 735 
 736 
 737 
 738   subroutine symmetry_get_n_atoms(id, nAtoms, errno)
 739     !scalars
 740 
 741 !This section has been created automatically by the script Abilint (TD).
 742 !Do not modify the following lines by hand.
 743 #undef ABI_FUNC
 744 #define ABI_FUNC 'symmetry_get_n_atoms'
 745 !End of the abilint section
 746 
 747     integer, intent(in) :: id
 748     integer, intent(out) :: errno
 749     integer, intent(out) :: nAtoms
 750 
 751     type(symmetry_list), pointer :: token
 752 
 753     if (AB_DBG) write(std_err,*) "AB symmetry: call get nAtoms."
 754 
 755     errno = AB7_NO_ERROR
 756     call get_item(token, id)
 757     if (.not. associated(token)) then
 758        errno = AB7_ERROR_OBJ
 759        return
 760     end if
 761 
 762     nAtoms = token%data%nAtoms
 763   end subroutine symmetry_get_n_atoms
 764 
 765   subroutine compute_bravais(sym)
 766 
 767 
 768 !This section has been created automatically by the script Abilint (TD).
 769 !Do not modify the following lines by hand.
 770 #undef ABI_FUNC
 771 #define ABI_FUNC 'compute_bravais'
 772 !End of the abilint section
 773 
 774     type(symmetry_type), intent(inout) :: sym
 775 
 776     integer :: berryopt
 777 
 778     ! We do the computation
 779     if (sym%withField) then
 780        berryopt = 4
 781     else
 782        berryopt = 0
 783     end if
 784     if (AB_DBG) write(std_err,*) "AB symmetry: call ABINIT symlatt."
 785     call symlatt(sym%bravais, AB7_MAX_SYMMETRIES, &
 786          & sym%nBravSym, sym%bravSym, sym%rprimd, sym%tolsym)
 787     if (AB_DBG) write(std_err,*) "AB symmetry: call ABINIT OK."
 788     if (AB_DBG) write(std_err, "(A,I3)") "  nSymBrav :", sym%nBravSym
 789     if (AB_DBG) write(std_err, "(A,I3)") "  holohedry:", sym%bravais(1)
 790     if (AB_DBG) write(std_err, "(A,I3)") "  center   :", sym%bravais(2)
 791   end subroutine compute_bravais
 792 
 793   subroutine symmetry_get_bravais(id, bravais, holohedry, center, &
 794        & nBravSym, bravSym, errno)
 795     !scalars
 796 
 797 !This section has been created automatically by the script Abilint (TD).
 798 !Do not modify the following lines by hand.
 799 #undef ABI_FUNC
 800 #define ABI_FUNC 'symmetry_get_bravais'
 801 !End of the abilint section
 802 
 803     integer, intent(in) :: id
 804     integer, intent(out) :: errno
 805     integer, intent(out) :: nBravSym, holohedry, center
 806     !arrays
 807     integer, intent(out) :: bravais(3,3), bravSym(3, 3, AB7_MAX_SYMMETRIES)
 808 
 809     type(symmetry_list), pointer :: token
 810 
 811     if (AB_DBG) write(std_err,*) "AB symmetry: call get bravais."
 812 
 813     errno = AB7_NO_ERROR
 814     call get_item(token, id)
 815     if (.not. associated(token)) then
 816        errno = AB7_ERROR_OBJ
 817        return
 818     end if
 819 
 820     if (token%data%nBravSym < 0) then
 821        ! We do the computation
 822        call compute_bravais(token%data)
 823     end if
 824 
 825     holohedry = token%data%bravais(1)
 826     center    = token%data%bravais(2)
 827     bravais   = reshape(token%data%bravais(3:11), (/ 3,3 /))
 828     nBravSym  = token%data%nBravSym
 829     bravSym(:, :, 1:nBravSym) = token%data%bravSym(:, :, 1:nBravSym)
 830   end subroutine symmetry_get_bravais
 831 
 832   subroutine compute_matrices(sym, errno)
 833 
 834 
 835 !This section has been created automatically by the script Abilint (TD).
 836 !Do not modify the following lines by hand.
 837 #undef ABI_FUNC
 838 #define ABI_FUNC 'compute_matrices'
 839 !End of the abilint section
 840 
 841     type(symmetry_type), intent(inout) :: sym
 842     integer, intent(out) :: errno
 843 
 844     integer :: berryopt, jellslab, noncol
 845     integer :: use_inversion
 846     real(dp), pointer :: spinAt_(:,:)
 847     integer  :: sym_(3, 3, AB7_MAX_SYMMETRIES)
 848     real(dp) :: transNon_(3, AB7_MAX_SYMMETRIES)
 849     integer  :: symAfm_(AB7_MAX_SYMMETRIES)
 850 
 851     errno = AB7_NO_ERROR
 852 
 853     if (sym%nBravSym < 0) then
 854        ! We do the computation of the Bravais part.
 855        call compute_bravais(sym)
 856     end if
 857 
 858     if (sym%withField) then
 859        berryopt = 4
 860     else
 861        berryopt = 0
 862     end if
 863     if (sym%withJellium) then
 864        jellslab = 1
 865     else
 866        jellslab = 0
 867     end if
 868     if (sym%withSpin == 4) then
 869        noncol = 1
 870        spinAt_ => sym%spinAt
 871     else if (sym%withSpin == 2) then
 872        noncol = 0
 873        spinAt_ => sym%spinAt
 874     else
 875        noncol = 0
 876        ABI_ALLOCATE(spinAt_,(3, sym%nAtoms))
 877        spinAt_ = 0
 878     end if
 879     if (sym%withSpinOrbit) then
 880        use_inversion = 0
 881     else
 882        use_inversion = 1
 883     end if
 884 
 885     if (sym%nsym == 0) then
 886        if (AB_DBG) write(std_err,*) "AB symmetry: call ABINIT symfind."
 887        call symfind(berryopt, sym%field, sym%gprimd, jellslab, AB7_MAX_SYMMETRIES, &
 888             & sym%nAtoms, noncol, sym%nBravSym, sym%nSym, sym%nzchempot, 0, sym%bravSym, spinAt_, &
 889             & symAfm_, sym_, transNon_, sym%tolsym, sym%typeAt, &
 890             & use_inversion, sym%xRed)
 891        if (AB_DBG) write(std_err,*) "AB symmetry: call ABINIT OK."
 892        if (AB_DBG) write(std_err, "(A,I3)") "  nSym:", sym%nSym
 893        if (associated(sym%sym))  then
 894           ABI_DEALLOCATE(sym%sym)
 895        end if
 896        if (associated(sym%symAfm))  then
 897           ABI_DEALLOCATE(sym%symAfm)
 898        end if
 899        if (associated(sym%transNon))  then
 900           ABI_DEALLOCATE(sym%transNon)
 901        end if
 902        ABI_ALLOCATE(sym%sym,(3, 3, sym%nSym))
 903        sym%sym(:,:,:) = sym_(:,:, 1:sym%nSym)
 904        ABI_ALLOCATE(sym%symAfm,(sym%nSym))
 905        sym%symAfm(:) = symAfm_(1:sym%nSym)
 906        ABI_ALLOCATE(sym%transNon,(3, sym%nSym))
 907        sym%transNon(:,:) = transNon_(:, 1:sym%nSym)
 908     else if (sym%nsym < 0) then
 909        sym%nsym = -sym%nsym
 910        sym_(:,:, 1:sym%nSym) = sym%sym(:,:,:)
 911        transNon_(:, 1:sym%nSym) = sym%transNon(:,:)
 912        symAfm_(1:sym%nSym) = sym%symAfm(:)
 913     end if
 914 
 915     if (sym%withSpin == 1) then
 916        ABI_DEALLOCATE(spinAt_)
 917     end if
 918 
 919     if (AB_DBG) write(std_err,*) "AB symmetry: call ABINIT symanal."
 920     call symanal(sym%bravais, 0, sym%genAfm, AB7_MAX_SYMMETRIES, sym%nSym, &
 921          & sym%pointGroupMagn, sym%rprimd, sym%spaceGroup, symAfm_, &
 922          & sym_, transNon_, sym%tolsym)
 923     if (AB_DBG) write(std_err,*) "AB symmetry: call ABINIT OK."
 924     sym%transNon(:,:) = transNon_(:, 1:sym%nSym)
 925 
 926     if (sym%bravais(1) < 0) then
 927        sym%multiplicity = 2
 928     else
 929        sym%multiplicity = 1
 930     end if
 931     if (AB_DBG) write(std_err, "(A,I3)") "  multi:", sym%multiplicity
 932     if (AB_DBG) write(std_err, "(A,I3)") "  space:", sym%spaceGroup
 933   end subroutine compute_matrices
 934 
 935   subroutine symmetry_get_n_sym(id, nSym, errno)
 936     !scalars
 937 
 938 !This section has been created automatically by the script Abilint (TD).
 939 !Do not modify the following lines by hand.
 940 #undef ABI_FUNC
 941 #define ABI_FUNC 'symmetry_get_n_sym'
 942 !End of the abilint section
 943 
 944     integer, intent(in) :: id
 945     integer, intent(out) :: errno
 946     integer, intent(out) :: nSym
 947 
 948     type(symmetry_list), pointer :: token
 949 
 950     if (AB_DBG) write(std_err,*) "AB symmetry: call get nSym."
 951 
 952     errno = AB7_NO_ERROR
 953     call get_item(token, id)
 954     if (.not. associated(token)) then
 955        errno = AB7_ERROR_OBJ
 956        return
 957     end if
 958 
 959     if (token%data%nSym <= 0) then
 960        ! We do the computation of the matrix part.
 961        call compute_matrices(token%data, errno)
 962     end if
 963 
 964     nSym = token%data%nSym
 965   end subroutine symmetry_get_n_sym
 966 
 967   subroutine symmetry_set_n_sym(id, nSym, sym, transNon, symAfm, errno)
 968     !scalars
 969 
 970 !This section has been created automatically by the script Abilint (TD).
 971 !Do not modify the following lines by hand.
 972 #undef ABI_FUNC
 973 #define ABI_FUNC 'symmetry_set_n_sym'
 974 !End of the abilint section
 975 
 976     integer, intent(in)  :: id
 977     integer, intent(in)  :: nSym
 978     integer, intent(in)  :: sym(3, 3, nSym)
 979     real(dp), intent(in) :: transNon(3, nSym)
 980     integer, intent(in)  :: symAfm(nSym)
 981     integer, intent(out) :: errno
 982 
 983     type(symmetry_list), pointer :: token
 984 
 985     if (AB_DBG) write(std_err,*) "AB symmetry: call get nSym."
 986 
 987     errno = AB7_NO_ERROR
 988     call get_item(token, id)
 989     if (.not. associated(token)) then
 990        errno = AB7_ERROR_OBJ
 991        return
 992     end if
 993 
 994     if (nSym <= 0) then
 995        errno = AB7_ERROR_ARG
 996        return
 997     else
 998        ABI_ALLOCATE(token%data%sym, (3, 3, nSym))
 999        token%data%sym(:,:,:) = sym(:,:,:)
1000        ABI_ALLOCATE(token%data%symAfm, (nSym))
1001        token%data%symAfm(:) = symAfm(:)
1002        ABI_ALLOCATE(token%data%transNon, (3, nSym))
1003        token%data%transNon(:,:) = transNon(:,:)
1004 
1005        token%data%auto = .false.
1006        token%data%nsym = -nSym
1007     end if
1008 
1009     ! We do the computation of the matrix part.
1010     call compute_matrices(token%data, errno)
1011   end subroutine symmetry_set_n_sym
1012 
1013   subroutine symmetry_get_matrices(id, nSym, sym, transNon, symAfm, errno)
1014 
1015 
1016 !This section has been created automatically by the script Abilint (TD).
1017 !Do not modify the following lines by hand.
1018 #undef ABI_FUNC
1019 #define ABI_FUNC 'symmetry_get_matrices'
1020 !End of the abilint section
1021 
1022     integer, intent(in) :: id
1023     integer, intent(out) :: errno
1024     integer, intent(out) :: nSym
1025     integer, intent(out)  :: sym(3, 3, AB7_MAX_SYMMETRIES)
1026     integer, intent(out)  :: symAfm(AB7_MAX_SYMMETRIES)
1027     real(dp), intent(out) :: transNon(3, AB7_MAX_SYMMETRIES)
1028 
1029     type(symmetry_list), pointer :: token
1030 
1031     if (AB_DBG) write(std_err,*) "AB symmetry: call get matrices."
1032 
1033     errno = AB7_NO_ERROR
1034     call get_item(token, id)
1035     if (.not. associated(token)) then
1036        errno = AB7_ERROR_OBJ
1037        return
1038     end if
1039 
1040     if (token%data%nSym <= 0) then
1041        ! We do the computation of the matrix part.
1042        call compute_matrices(token%data, errno)
1043     end if
1044 
1045     nSym                = token%data%nSym
1046     sym(:, :, 1:nSym)   = token%data%sym(:, :,:)
1047     symAfm(1:nSym)      = token%data%symAfm(:)
1048     transNon(:, 1:nSym) = token%data%transNon(:,:)
1049   end subroutine symmetry_get_matrices
1050 
1051   subroutine symmetry_get_matrices_p(id, nSym, sym, transNon, symAfm, errno)
1052 
1053 
1054 !This section has been created automatically by the script Abilint (TD).
1055 !Do not modify the following lines by hand.
1056 #undef ABI_FUNC
1057 #define ABI_FUNC 'symmetry_get_matrices_p'
1058 !End of the abilint section
1059 
1060     integer, intent(in) :: id
1061     integer, intent(out) :: errno
1062     integer, intent(out) :: nSym
1063     integer, pointer  :: sym(:,:,:)
1064     integer, pointer  :: symAfm(:)
1065     real(dp), pointer :: transNon(:,:)
1066 
1067     type(symmetry_list), pointer :: token
1068 
1069     if (AB_DBG) write(std_err,*) "AB symmetry: call get matrices as pointers."
1070 
1071     errno = AB7_NO_ERROR
1072     call get_item(token, id)
1073     if (.not. associated(token)) then
1074        errno = AB7_ERROR_OBJ
1075        return
1076     end if
1077 
1078     if (token%data%nSym <= 0) then
1079        ! We do the computation of the matrix part.
1080        call compute_matrices(token%data, errno)
1081     end if
1082 
1083     nSym     =  token%data%nSym
1084     sym      => token%data%sym
1085     symAfm   => token%data%symAfm
1086     transNon => token%data%transNon
1087   end subroutine symmetry_get_matrices_p
1088 
1089   subroutine symmetry_get_multiplicity(id, multiplicity, errno)
1090 
1091 
1092 !This section has been created automatically by the script Abilint (TD).
1093 !Do not modify the following lines by hand.
1094 #undef ABI_FUNC
1095 #define ABI_FUNC 'symmetry_get_multiplicity'
1096 !End of the abilint section
1097 
1098     integer, intent(in) :: id
1099     integer, intent(out) :: multiplicity, errno
1100 
1101     type(symmetry_list), pointer :: token
1102 
1103     if (AB_DBG) write(std_err,*) "AB symmetry: call get multiplicity."
1104 
1105     errno = AB7_NO_ERROR
1106     call get_item(token, id)
1107     if (.not. associated(token)) then
1108        errno = AB7_ERROR_OBJ
1109        return
1110     end if
1111 
1112     if (token%data%multiplicity < 0) then
1113        ! We do the computation of the matrix part.
1114        call compute_matrices(token%data, errno)
1115     end if
1116     multiplicity = token%data%multiplicity
1117   end subroutine symmetry_get_multiplicity
1118 
1119   subroutine symmetry_get_group(id, spaceGroup, spaceGroupId, &
1120        & pointGroupMagn, genAfm, errno)
1121 
1122 
1123 !This section has been created automatically by the script Abilint (TD).
1124 !Do not modify the following lines by hand.
1125 #undef ABI_FUNC
1126 #define ABI_FUNC 'symmetry_get_group'
1127 !End of the abilint section
1128 
1129     integer, intent(in)            :: id
1130     integer, intent(out)           :: errno
1131     real(dp), intent(out)          :: genAfm(3)
1132     character(len=15), intent(out) :: spaceGroup
1133     integer, intent(out)           :: spaceGroupId, pointGroupMagn
1134 
1135     type(symmetry_list), pointer  :: token
1136     integer :: sporder
1137     character(len=1)  :: brvLattice
1138     character(len=15) :: ptintsb,ptschsb,schsb,spgrp
1139     character(len=35) :: intsbl
1140 
1141     if (AB_DBG) write(std_err,*) "AB symmetry: call get group."
1142 
1143     errno = AB7_NO_ERROR
1144     call get_item(token, id)
1145     if (.not. associated(token)) then
1146        errno = AB7_ERROR_OBJ
1147        return
1148     end if
1149 
1150     if (token%data%multiplicity < 0) then
1151        ! We do the computation of the matrix part.
1152        call compute_matrices(token%data, errno)
1153     end if
1154 
1155     if (token%data%multiplicity /= 1) then
1156        errno = AB7_ERROR_SYM_NOT_PRIMITIVE
1157        return
1158     end if
1159 
1160     call spgdata(brvLattice,spgrp,intsbl,ptintsb,ptschsb,&
1161          &  schsb,1,token%data%spaceGroup,sporder,1)
1162 
1163     write(spaceGroup, "(3A)") brvLattice, " ", trim(spgrp(1:13))
1164     pointGroupMagn = token%data%pointGroupMagn
1165     spaceGroupId   = token%data%spaceGroup
1166     genAfm         = token%data%genAfm
1167   end subroutine symmetry_get_group
1168 
1169   subroutine compute_equivalent_atoms(sym)
1170 
1171 
1172 !This section has been created automatically by the script Abilint (TD).
1173 !Do not modify the following lines by hand.
1174 #undef ABI_FUNC
1175 #define ABI_FUNC 'compute_equivalent_atoms'
1176 !End of the abilint section
1177 
1178     type(symmetry_type), intent(inout) :: sym
1179 
1180     integer, allocatable :: symrec(:,:,:)
1181     integer :: isym
1182 
1183     if (.not. associated(sym%indexingAtoms)) then
1184       ABI_ALLOCATE(sym%indexingAtoms,(4, sym%nSym, sym%nAtoms))
1185     end if
1186 
1187     !Get the symmetry matrices in terms of reciprocal basis
1188     ABI_ALLOCATE(symrec,(3, 3, sym%nSym))
1189     do isym = 1, sym%nSym, 1
1190        call mati3inv(sym%sym(:,:,isym), symrec(:,:,isym))
1191     end do
1192 
1193     !Obtain a list of rotated atom labels:
1194     call symatm(sym%indexingAtoms, sym%nAtoms, sym%nSym, symrec, &
1195          & sym%transNon, sym%tolsym, sym%typeAt, sym%xRed)
1196 
1197     ABI_DEALLOCATE(symrec)
1198   end subroutine compute_equivalent_atoms
1199 
1200   subroutine symmetry_get_equivalent_atom(id, equiv, iAtom, errno)
1201 
1202 
1203 !This section has been created automatically by the script Abilint (TD).
1204 !Do not modify the following lines by hand.
1205 #undef ABI_FUNC
1206 #define ABI_FUNC 'symmetry_get_equivalent_atom'
1207 !End of the abilint section
1208 
1209     integer, intent(in)  :: id
1210     integer, intent(in)  :: iAtom
1211     integer, intent(out) :: equiv(4, AB7_MAX_SYMMETRIES)
1212     integer, intent(out) :: errno
1213 
1214     type(symmetry_list), pointer  :: token
1215 
1216     if (AB_DBG) write(std_err,*) "AB symmetry: call get equivalent."
1217 
1218     errno = AB7_NO_ERROR
1219     call get_item(token, id)
1220     if (.not. associated(token)) then
1221        errno = AB7_ERROR_OBJ
1222        return
1223     end if
1224 
1225     if (iAtom < 1 .or. iAtom > token%data%nAtoms) then
1226        errno = AB7_ERROR_ARG
1227        return
1228     end if
1229 
1230     if (.not. associated(token%data%indexingAtoms)) then
1231        ! We do the computation of the matrix part.
1232        call compute_equivalent_atoms(token%data)
1233     end if
1234 
1235     equiv(:, 1:token%data%nSym) = token%data%indexingAtoms(:,:,iAtom)
1236   end subroutine symmetry_get_equivalent_atom
1237 
1238   subroutine symmetry_get_type(id, iSym, label, type, errno)
1239 
1240 
1241 !This section has been created automatically by the script Abilint (TD).
1242 !Do not modify the following lines by hand.
1243 #undef ABI_FUNC
1244 #define ABI_FUNC 'symmetry_get_type'
1245 !End of the abilint section
1246 
1247     integer, intent(in)  :: id, iSym
1248     character(len = 128), intent(out) :: label
1249     integer, intent(out) :: errno, type
1250 
1251     integer :: det
1252     type(symmetry_list), pointer :: token
1253     type(symmetry_type), pointer :: sym
1254 
1255     if (AB_DBG) write(std_err,*) "AB symmetry: call get type."
1256 
1257     errno = AB7_NO_ERROR
1258     call get_item(token, id)
1259     if (.not. associated(token)) then
1260        errno = AB7_ERROR_OBJ
1261        return
1262     end if
1263     sym => token%data
1264 
1265     if (iSym < 1 .or. iSym > sym%nSym) then
1266        errno = AB7_ERROR_ARG
1267        return
1268     end if
1269 
1270     if (sym%multiplicity < 0) then
1271        ! We do the computation of the matrix part.
1272        call compute_matrices(sym, errno)
1273     end if
1274 
1275     ! Calculate determinant.
1276     call mati3det(sym%sym(:,:,iSym),det)
1277     call symcharac(sym%bravais(2), det, sym%bravais(1), iSym, label, &
1278          sym%sym(:,:,iSym), sym%transNon(:, iSym), type)
1279   end subroutine symmetry_get_type
1280 
1281 end module m_ab7_symmetry