TABLE OF CONTENTS


ABINIT/m_abi_linalg [ Modules ]

[ Top ] [ Modules ]

NAME

  m_abi_linalg

FUNCTION

  management of Linear Algebra wrappers routines
 with support of different external library (scalapack, elpa, plasma, magma, ... )

COPYRIGHT

  Copyright (C) 2012-2018 ABINIT group (LNguyen,FDahm)
  This file is distributed under the terms of the
  GNU General Public License, see ~ABINIT/Infos/copyright
  or http://www.gnu.org/copyleft/gpl.txt .

SOURCE

18 #if defined HAVE_CONFIG_H
19 #include "config.h"
20 #endif
21 
22 #include "abi_common.h"
23 
24 module m_abi_linalg
25 
26   use defs_basis
27   use m_errors
28   use m_abicore
29   use m_xmpi
30   use m_xomp
31   use m_slk
32   use iso_c_binding
33 #ifdef HAVE_LINALG_ELPA
34  use m_elpa
35 #endif
36 #ifdef HAVE_LINALG_PLASMA
37  use plasma, except_dp => dp, except_sp => sp
38 #endif
39  use m_time,  only : timab
40  
41 
42  implicit none
43 
44  private

m_abi_linalg/abi_linalg_eigen_setmaxsize [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

 abi_linalg_eigen_setmaxsize

FUNCTION

INPUTS

PARENTS

      m_abi_linalg

CHILDREN

SOURCE

415  subroutine abi_linalg_eigen_setmaxsize(max_eigen_pb_size)
416 
417 
418 !This section has been created automatically by the script Abilint (TD).
419 !Do not modify the following lines by hand.
420 #undef ABI_FUNC
421 #define ABI_FUNC 'abi_linalg_eigen_setmaxsize'
422 !End of the abilint section
423 
424  implicit none
425 
426 !Arguments ------------------------------------
427  integer,intent(in) :: max_eigen_pb_size
428 
429 !******************************************************************
430 
431 !Single precision work arrays (real and complex)
432  eigen_s_maxsize = max_eigen_pb_size
433  eigen_s_lwork = max(1,3*eigen_s_maxsize)
434  if(allocated(eigen_s_work)) then
435    ABI_DEALLOCATE(eigen_s_work)
436  end if
437  ABI_ALLOCATE(eigen_s_work,(max(1,eigen_s_lwork)))
438 
439  eigen_c_maxsize = max_eigen_pb_size
440  eigen_c_lwork = max(1,2*eigen_c_maxsize-1)
441  if(allocated(eigen_c_work)) then
442    ABI_DEALLOCATE(eigen_c_work)
443  end if
444 
445  ABI_ALLOCATE(eigen_c_work,(max(1,eigen_c_lwork)))
446  if(allocated(eigen_c_rwork)) then
447    ABI_DEALLOCATE(eigen_c_rwork)
448  end if
449  ABI_ALLOCATE(eigen_c_rwork,(max(1, 3*eigen_c_maxsize-2)))
450 
451 !Double precision work arrays (real and complex)
452  eigen_d_maxsize = max_eigen_pb_size
453  eigen_d_lwork = max(eigen_d_lwork,2*(2*eigen_d_maxsize - 1)) !for zh[ep][eg]v
454  eigen_d_lwork = max(eigen_d_lwork,(3*eigen_d_maxsize))       !for ds[ep][eg]v
455 #ifdef HAVE_LINALG_MAGMA
456  eigen_d_lwork = max(eigen_d_lwork,2*(2*(eigen_d_maxsize**2) + 6*eigen_d_maxsize +1))!for magma_zhe[eg]v
457  eigen_d_lwork = max(eigen_d_lwork, eigen_d_maxsize**2 + 33*eigen_d_maxsize)         !for magma_dsy[eg]v
458 #endif
459  if(allocated(eigen_d_work)) then
460    ABI_DEALLOCATE(eigen_d_work)
461  end if
462  ABI_ALLOCATE(eigen_d_work,(max(1,eigen_d_lwork)))
463 
464  eigen_z_maxsize = max_eigen_pb_size
465  eigen_z_lwork = max(1,2*eigen_z_maxsize-1)
466  if(allocated(eigen_z_work)) then
467    ABI_DEALLOCATE(eigen_z_work)
468  end if
469  ABI_ALLOCATE(eigen_z_work,(max(1,eigen_z_lwork)))
470  if(allocated(eigen_z_rwork)) then
471    ABI_DEALLOCATE(eigen_z_rwork)
472  end if
473 
474 #ifdef HAVE_LINALG_MAGMA
475  ABI_ALLOCATE(eigen_z_rwork,(max(1, 2*(eigen_z_maxsize**2) + 5*eigen_z_maxsize + 1)))
476 #else
477  ABI_ALLOCATE(eigen_z_rwork,(max(1, 3*eigen_z_maxsize-2)))
478 #endif
479 
480  end subroutine abi_linalg_eigen_setmaxsize

m_abi_linalg/abi_linalg_finalize [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

 abi_linalg_finalize

FUNCTION

INPUTS

PARENTS

      compute_kgb_indicator,driver

CHILDREN

SOURCE

498  subroutine abi_linalg_finalize(only_scalapack) ! optional argument
499 
500 
501 !This section has been created automatically by the script Abilint (TD).
502 !Do not modify the following lines by hand.
503 #undef ABI_FUNC
504 #define ABI_FUNC 'abi_linalg_finalize'
505 !End of the abilint section
506 
507  implicit none
508 
509 !Arguments ------------------------------------
510  logical,intent(in),optional :: only_scalapack
511 
512 !Local variables ------------------------------
513 #ifdef HAVE_LINALG_PLASMA
514  integer :: info
515 #endif
516 
517 !******************************************************************
518 
519 #ifdef HAVE_LINALG_SCALAPACK
520  if (abi_communicator/=xmpi_comm_null) then
521    call xmpi_comm_free(abi_communicator)
522    call end_scalapack(abi_processor)
523  end if
524 #endif
525  if (present(only_scalapack)) then
526    if (only_scalapack) return
527  end if
528 
529 #ifdef HAVE_LINALG_ELPA
530  call elpa_func_uninit()
531 #endif
532 
533 #ifdef HAVE_LINALG_PLASMA
534  call PLASMA_Finalize(info)
535 #endif
536 
537 #ifdef HAVE_LINALG_MAGMA
538 #ifdef HAVE_LINALG_MAGMA_15
539  call magmaf_finalize()
540 #endif
541 #endif
542 
543 #ifdef HAVE_GPU_CUDA
544  call gpu_linalg_shutdown()
545 #endif
546 
547 !Memory freeing
548  if(allocated(eigen_s_work)) then
549    ABI_DEALLOCATE(eigen_s_work)
550  end if
551  if(allocated(eigen_d_work)) then
552    ABI_DEALLOCATE(eigen_d_work)
553  end if
554  if(allocated(eigen_c_work)) then
555    ABI_DEALLOCATE(eigen_c_work)
556  end if
557  if(allocated(eigen_c_rwork)) then
558    ABI_DEALLOCATE(eigen_c_rwork)
559  end if
560  if(allocated(eigen_z_work)) then
561    ABI_DEALLOCATE(eigen_z_work)
562  end if
563  if(allocated(eigen_z_rwork)) then
564    ABI_DEALLOCATE(eigen_z_rwork)
565  end if
566 
567  eigen_s_maxsize = 0
568  eigen_d_maxsize = 0
569  eigen_c_maxsize = 0
570  eigen_z_maxsize = 0
571  eigen_s_lwork   = 0
572  eigen_d_lwork   = 0
573  eigen_c_lwork   = 0
574  eigen_z_lwork   = 0
575 
576  end subroutine abi_linalg_finalize

m_abi_linalg/abi_linalg_init [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

 abi_linalg_init

FUNCTION

 Initalization of linear algebra environnement

INPUTS

PARENTS

      compute_kgb_indicator,driver

CHILDREN

SOURCE

278  subroutine abi_linalg_init(comm_scalapack,eigen_group_size,max_eigen_pb_size,my_rank,&
279 &                           only_scalapack) ! optional parameter
280 
281 #if defined HAVE_MPI2
282   use mpi
283 #endif
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 'abi_linalg_init'
289 !End of the abilint section
290 
291  implicit none
292 
293 #if defined HAVE_MPI1
294  include 'mpif.h'
295 #endif
296 
297 !Arguments ------------------------------------
298  integer,intent(in) :: comm_scalapack,eigen_group_size,max_eigen_pb_size,my_rank
299  logical,intent(in),optional :: only_scalapack
300 
301 !Local variables ------------------------------
302 #ifdef HAVE_LINALG_SCALAPACK
303  integer :: abi_info1,rank,commsize
304  integer :: sizecart(2)
305  logical :: periodic(2), reorder, keepdim(2)
306  integer :: commcart
307 #endif
308 #ifdef HAVE_LINALG_PLASMA
309  integer :: abi_info2,core_id,num_cores,num_cores_node
310  integer,allocatable :: affinity(:)
311 #endif
312 
313 !******************************************************************
314 
315 #ifdef HAVE_LINALG_SCALAPACK
316 !Scalapack initalization
317  if (eigen_group_size>0) then
318    abi_communicator = comm_scalapack
319    rank=xmpi_comm_rank(comm_scalapack)
320 
321    ! We create abi_communicator using a cartesian grid, and store its complement
322    commsize = MIN(eigen_group_size, xmpi_comm_size(comm_scalapack))
323    sizecart = (/commsize, xmpi_comm_size(comm_scalapack)/commsize/)
324    periodic = (/.true.,.true./)
325    reorder = .false.
326    call MPI_CART_CREATE(comm_scalapack,2,sizecart,periodic,reorder,commcart,abi_info1)
327    keepdim = (/.true., .false./)
328    call MPI_CART_SUB(commcart, keepdim, abi_communicator,abi_info1)
329    keepdim = (/.false., .true./)
330    call MPI_CART_SUB(commcart, keepdim, abi_complement_communicator,abi_info1)
331 
332    call init_scalapack(abi_processor,abi_communicator)
333  else
334    abi_communicator=xmpi_comm_null
335    abi_complement_communicator = xmpi_comm_null
336  end if
337 #endif
338  if (present(only_scalapack)) then
339    if (only_scalapack) return
340  end if
341 
342 #ifdef HAVE_LINALG_ELPA
343  call elpa_func_init()
344 #endif
345 
346 #ifdef HAVE_LINALG_PLASMA
347 !Plasma Initialization
348 !Because use of hybrid use of mpi+openmp+plasma,
349 !  we need to set manually the thread bindings policy
350 !  to avoid conflicts between mpi process due to plasma
351  num_cores=xomp_get_max_threads()
352  num_cores_node=xomp_get_num_cores_node()
353 
354  if (num_cores_node == 0) then ! This means that OMP is not enabled.
355    num_cores_node = 1
356    MSG_WARNING("You are using PLASMA but OpenMP is not enabled in Abinit!")
357  end if
358 
359  ABI_ALLOCATE(affinity,(num_cores))
360  do core_id =1,num_cores
361    affinity(core_id) = MOD(my_rank*num_cores + (core_id-1), num_cores_node)
362  end do
363  call PLASMA_Init_Affinity(num_cores,affinity(1),abi_info2)
364  ABI_DEALLOCATE(affinity)
365  XPLASMA_ISON = .True.
366 #endif
367 
368 #ifdef HAVE_LINALG_MAGMA
369 #ifdef HAVE_LINALG_MAGMA_15
370  call magmaf_init()
371 #endif
372 #endif
373 
374 #ifdef HAVE_GPU_CUDA
375 !Cublas initialization
376  call gpu_linalg_init()
377 #endif
378 
379 !Allocation of work sizes
380  eigen_s_maxsize = 0
381  eigen_d_maxsize = 0
382  eigen_c_maxsize = 0
383  eigen_z_maxsize = 0
384  eigen_s_lwork   = 0
385  eigen_d_lwork   = 0
386  eigen_c_lwork   = 0
387  eigen_z_lwork   = 0
388 
389  call abi_linalg_eigen_setmaxsize(max_eigen_pb_size)
390 
391  return
392 
393  ABI_UNUSED(comm_scalapack)
394  ABI_UNUSED(eigen_group_size)
395  ABI_UNUSED(my_rank)
396 
397  end subroutine abi_linalg_init

m_abi_linalg/diag_plasma [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

FUNCTION

  Convert diag character to PLASMA integer

PARENTS

SOURCE

873 integer function diag_plasma(diag)
874 
875 
876 !This section has been created automatically by the script Abilint (TD).
877 !Do not modify the following lines by hand.
878 #undef ABI_FUNC
879 #define ABI_FUNC 'diag_plasma'
880 !End of the abilint section
881 
882  implicit none
883 
884 !Arguments ------------------------------------
885 !scalars
886  character(len=1),intent(in) :: diag
887 
888 ! *************************************************************************
889 
890  if (LSAME(diag,'U')) then
891    diag_plasma = PlasmaUnit
892  else
893    diag_plasma = PlasmaNonUnit
894  end if
895 
896 end function diag_plasma

m_abi_linalg/jobz_plasma [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

FUNCTION

  Convert jobz character to PLASMA integer

PARENTS

SOURCE

911 integer function jobz_plasma(jobz)
912 
913 
914 !This section has been created automatically by the script Abilint (TD).
915 !Do not modify the following lines by hand.
916 #undef ABI_FUNC
917 #define ABI_FUNC 'jobz_plasma'
918 !End of the abilint section
919 
920  implicit none
921 
922 !Arguments ------------------------------------
923 !scalars
924  character(len=1),intent(in) :: jobz
925 
926 ! *************************************************************************
927 
928  if (LSAME(jobz,'N')) then
929    jobz_plasma = PlasmaNoVec
930  else
931    jobz_plasma = PlasmaVec
932  end if
933 
934 end function jobz_plasma

m_abi_linalg/laparams_t [ Types ]

[ Top ] [ m_abi_linalg ] [ Types ]

NAME

 laparams_t

FUNCTION

  Gather the parameters and the options to be passed to the
  linear algebra routines.

SOURCE

 92  integer,public,parameter :: LIB_LINALG=1
 93  integer,public,parameter :: LIB_MAGMA=2
 94  integer,public,parameter :: LIB_PLASMA=3
 95  integer,public,parameter :: LIB_SLK=4
 96 
 97  type,public :: laparams_t
 98 
 99    integer :: lib = LIB_LINALG
100      ! Flag defining the library to use.
101 
102    integer :: comm = xmpi_comm_null
103      ! MPI communicator (used if lib == LIB_SLK)
104 
105    !integer :: timopt, tim
106    ! Options for timing.
107  end type laparams_t

m_abi_linalg/linalg_allow_gemm3m [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

FUNCTION

  Programmatic interface to enable the use of [Z,C]GEMM3M calls

SOURCE

589 subroutine linalg_allow_gemm3m(bool)
590 
591 
592 !This section has been created automatically by the script Abilint (TD).
593 !Do not modify the following lines by hand.
594 #undef ABI_FUNC
595 #define ABI_FUNC 'linalg_allow_gemm3m'
596 !End of the abilint section
597 
598  implicit none
599 
600 !Arguments ------------------------------------
601 !scalars
602  logical,intent(in) :: bool
603 
604 ! *************************************************************************
605 
606  XGEMM3M_ISON = bool
607 
608 end subroutine linalg_allow_gemm3m

m_abi_linalg/linalg_allow_plasma [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

FUNCTION

  Programmatic interface to enable the use of PLASMA
  False to disable PLASMA version.

SOURCE

719 subroutine linalg_allow_plasma(bool)
720 
721 
722 !This section has been created automatically by the script Abilint (TD).
723 !Do not modify the following lines by hand.
724 #undef ABI_FUNC
725 #define ABI_FUNC 'linalg_allow_plasma'
726 !End of the abilint section
727 
728  implicit none
729 
730 !Arguments ------------------------------------
731 !scalars
732  logical,intent(in) :: bool
733 
734 ! *************************************************************************
735 
736  XPLASMA_ISON = bool
737 #ifndef HAVE_LINALG_PLASMA
738  ! Just to be on the safe-side.
739  ! I have to use a weird set of branches to make abirules happy in the BLAS/LAPACK
740  ! wrappers, and one cannot set XPLASMA_MODE to .True. if PLASMA is not available.
741  XPLASMA_ISON = .False.
742 #endif
743 
744 end subroutine linalg_allow_plasma

m_abi_linalg/side_plasma [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

FUNCTION

  Convert side character to PLASMA integer

PARENTS

SOURCE

835 integer function side_plasma(side)
836 
837 
838 !This section has been created automatically by the script Abilint (TD).
839 !Do not modify the following lines by hand.
840 #undef ABI_FUNC
841 #define ABI_FUNC 'side_plasma'
842 !End of the abilint section
843 
844  implicit none
845 
846 !Arguments ------------------------------------
847 !scalars
848  character(len=1),intent(in) :: side
849 
850 ! *************************************************************************
851 
852  if(LSAME(side,'L')) then
853     side_plasma = PlasmaLeft
854  else
855     side_plasma = PlasmaRight
856  end if
857 
858 end function side_plasma

m_abi_linalg/trans_plasma [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

FUNCTION

  Convert trans character to PLASMA integer

PARENTS

SOURCE

795 integer function trans_plasma(trans)
796 
797 
798 !This section has been created automatically by the script Abilint (TD).
799 !Do not modify the following lines by hand.
800 #undef ABI_FUNC
801 #define ABI_FUNC 'trans_plasma'
802 !End of the abilint section
803 
804  implicit none
805 
806 !Arguments ------------------------------------
807 !scalars
808  character(len=1),intent(in) :: trans
809 
810 ! *************************************************************************
811 
812  if (LSAME(trans,'C')) then
813    trans_plasma = PlasmaConjTrans
814  else if (LSAME(trans,'T')) then
815    trans_plasma = PlasmaTrans
816  else
817    trans_plasma = PlasmaNoTrans
818  end if
819 
820 end function trans_plasma

m_abi_linalg/uplo_plasma [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

FUNCTION

  Convert uplo character to PLASMA integer

SOURCE

757 integer function uplo_plasma(uplo)
758 
759 
760 !This section has been created automatically by the script Abilint (TD).
761 !Do not modify the following lines by hand.
762 #undef ABI_FUNC
763 #define ABI_FUNC 'uplo_plasma'
764 !End of the abilint section
765 
766  implicit none
767 
768 !Arguments ------------------------------------
769 !scalars
770  character(len=1),intent(in) :: uplo
771 
772 ! *************************************************************************
773 
774  if (LSAME(uplo,'U')) then
775     uplo_plasma = PlasmaUpper
776  else
777     uplo_plasma = PlasmaLower
778  end if
779 
780 end function uplo_plasma

m_abi_linalg/use_cgemm3m [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

  use_cgemm3m

FUNCTION

  Enable the use of CGEMM3M

PARENTS

NOTES

  See use_zgemm3m

SOURCE

682 pure logical function use_cgemm3m(m,n,k)
683 
684 
685 !This section has been created automatically by the script Abilint (TD).
686 !Do not modify the following lines by hand.
687 #undef ABI_FUNC
688 #define ABI_FUNC 'use_cgemm3m'
689 !End of the abilint section
690 
691  implicit none
692 
693 !Arguments ------------------------------------
694 !scalars
695  integer,intent(in) :: m,n,k
696 
697 ! *************************************************************************
698 
699  use_cgemm3m = .False.
700  if (XGEMM3M_ISON) use_cgemm3m = ((m * n * k) > CGEMM3M_LIMIT)
701 #ifndef HAVE_LINALG_GEMM3M
702  use_cgemm3m = .False.
703 #endif
704 
705 end function use_cgemm3m

m_abi_linalg/use_zgemm3m [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

  use_zgemm3m

FUNCTION

  Enable the use of ZGEMM3M

PARENTS

NOTES

  The CGEMM3M and ZGEMM3M routines use an algorithm requiring 3 real matrix
  multiplications and 5 real matrix additions to compute the complex matrix
  product; CGEMM(3S) and ZGEMM(3S) use 4 real matrix multiplications and 2
  real matrix additions. Because the matrix multiplication time is usually
  the limiting performance factor in these routines, CGEMM3M and ZGEMM3M
  may run up to 33 percent faster than CGEMM and ZGEMM.  Because of other
  overhead associated with the 3M routines, however, these performance
  improvements may not always be realized.  For example, on one processor
  the 3M routines will generally run more slowly than the standard complex
  matrix multiplication routines when m * n * k < FACTOR, where m, n, and k
  are the input matrix dimensions and FACTOR is approximately 200000 for
  CGEMM3M and 325000 for ZGEMM3M.
  from: http://techpubs.sgi.com/library/tpl/cgi-bin/getdoc.cgi?coll=0650&db=man&raw=1&fname=/usr/share/catman/p_man/cat3/SCSL/ZGEMM3M.z

SOURCE

639 pure logical function use_zgemm3m(m,n,k)
640 
641 
642 !This section has been created automatically by the script Abilint (TD).
643 !Do not modify the following lines by hand.
644 #undef ABI_FUNC
645 #define ABI_FUNC 'use_zgemm3m'
646 !End of the abilint section
647 
648  implicit none
649 
650 !Arguments ------------------------------------
651 !scalars
652  integer,intent(in) :: m,n,k
653 
654 ! *************************************************************************
655 
656  use_zgemm3m = .False.
657  if (XGEMM3M_ISON) use_zgemm3m = ((m * n * k) > ZGEMM3M_LIMIT)
658 
659 #ifndef HAVE_LINALG_GEMM3M
660  use_zgemm3m = .False.
661 #endif
662 
663 end function use_zgemm3m