TABLE OF CONTENTS


m_abi_linalg/abi_cheev [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

 abi_cheev

FUNCTION

INPUTS

SOURCE

243 subroutine abi_cheev(jobz,uplo,n,a,lda,w,work,lwork,rwork,info)
244 
245 !This section has been created automatically by the script Abilint (TD).
246 !Do not modify the following lines by hand.
247 #undef ABI_FUNC
248 #define ABI_FUNC 'abi_cheev'
249 !End of the abilint section
250 
251  implicit none
252 
253 !Arguments ------------------------------------
254  character(len=1), intent(in) :: jobz
255  character(len=1), intent(in) :: uplo
256  integer, intent(in) :: n
257  integer, intent(in) :: lda
258  integer, intent(in) :: lwork
259  complex(spc), target, intent(inout) :: a(lda,*)
260  complex(spc), target, intent(inout) :: work(*)
261  real(sp),target, intent(inout) :: rwork(lwork)
262  real(sp),target, intent(out) :: w(n)
263  integer, intent(out) :: info
264 
265 #ifdef HAVE_LINALG_PLASMA
266  integer :: jobz_plasma_a
267  type(c_ptr) :: plasma_work
268 #endif
269     ! *********************************************************************
270 
271 #ifdef HAVE_LINALG_PLASMA
272  ! FDahm & LNGuyen  (November 2012) :
273  !  In Plasma v 2.4.6, eigen routines support only
274  !the eigenvalues computation (jobz=N) and not the
275  !full eigenvectors bases determination (jobz=V)
276  if (LSAME(jobz,'N')) then
277     jobz_plasma_a = jobz_plasma(jobz)
278 
279     call PLASMA_Alloc_Workspace_cheev(n,n,plasma_work,info)
280 
281     info = PLASMA_cheev_c(jobz_plasma_a,uplo_plasma(uplo),&
282 &     n,c_loc(a),lda,c_loc(w),plasma_work,c_loc(rwork),lwork)
283 
284     call PLASMA_Dealloc_handle(plasma_work,info)
285  else
286 #endif
287    call cheev(jobz,uplo,n,a,lda,w,work,lwork,rwork,info)
288 #ifdef HAVE_LINALG_PLASMA
289  end if
290 #endif
291 
292 end subroutine abi_cheev

m_abi_linalg/abi_cheev_alloc [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

 abi_cheev_alloc

FUNCTION

INPUTS

PARENTS

SOURCE

309   subroutine abi_cheev_alloc(jobz,uplo,n,a,w)
310 
311 !This section has been created automatically by the script Abilint (TD).
312 !Do not modify the following lines by hand.
313 #undef ABI_FUNC
314 #define ABI_FUNC 'abi_cheev_alloc'
315 
316 !End of the abilint section
317 
318  implicit none
319  !Arguments ------------------------------------
320  character(len=1), intent(in) :: jobz
321  character(len=1), intent(in) :: uplo
322  integer, intent(in) :: n
323  complex(spc), DEV_CONTARRD intent(inout) :: a(:,:) ! Here I need lda to be consistent
324  real(sp), intent(out) :: w(n)
325 
326  integer :: info
327 
328 ! *********************************************************************
329 
330  call abi_cheev(jobz,uplo,n,a,n,w,eigen_c_work,eigen_c_lwork,eigen_c_rwork,info)
331 
332 end subroutine abi_cheev_alloc

m_abi_linalg/abi_dheev [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

 abi_dheev

FUNCTION

INPUTS

PARENTS

SOURCE

 33   subroutine abi_dheev(jobz,uplo,n,a,lda,w,work,lwork,rwork,info, &
 34 &       x_cplx,istwf_k,timopt,tim_xeigen,use_slk,use_gpu)
 35 
 36 !This section has been created automatically by the script Abilint (TD).
 37 !Do not modify the following lines by hand.
 38 #undef ABI_FUNC
 39 #define ABI_FUNC 'abi_dheev'
 40 !End of the abilint section
 41 
 42  implicit none
 43 
 44 !Arguments ------------------------------------
 45  character(len=1), intent(in) :: jobz
 46  character(len=1), intent(in) :: uplo
 47  integer, intent(in) :: n
 48  integer, intent(in) :: lda
 49  integer, intent(in) :: lwork
 50  integer, intent(out) :: info
 51  real(dp),target, intent(inout) :: a(lda,*)  ! FIXME should be cplex * lda
 52  real(dp),target,intent(inout) :: work(*)
 53  real(dp),optional,target,intent(inout) :: rwork(lwork)
 54  real(dp), target,intent(out) :: w(n)
 55  integer, optional, intent(in) :: x_cplx
 56  integer, optional, intent(in) :: timopt,tim_xeigen,use_gpu
 57  integer, optional, intent(in) :: use_slk
 58  integer, optional, intent(in) :: istwf_k
 59 
 60 !Local variables-------------------------------
 61  character(len=500) :: msg
 62  integer :: istwf_k_ , usegpu_,use_slk_
 63  real(dp) :: tsec(2)
 64  integer  :: cplx_
 65 #ifdef HAVE_LINALG_MAGMA
 66  integer :: liwork_,lwork_, lrwork_,iwk_(1)
 67  integer, dimension(:), allocatable :: iwork
 68  real(dp) :: rwk_(1)
 69  complex(dpc) :: cwk_(1)
 70 #endif
 71 #ifdef HAVE_LINALG_PLASMA
 72  integer :: jobz_plasma_a
 73  type(c_ptr) :: plasma_work
 74 #endif
 75 
 76 ! *********************************************************************
 77 
 78  if (present(tim_xeigen).and.present(timopt)) then
 79    if(abs(timopt)==3) then
 80      call timab(tim_xeigen,1,tsec)
 81    end if
 82  end if
 83 
 84  cplx_=1 ; if(PRESENT(x_cplx)) cplx_ = x_cplx
 85  usegpu_=0;if (present(use_gpu)) usegpu_=use_gpu
 86  istwf_k_=1;if (present(istwf_k)) istwf_k_=istwf_k
 87  use_slk_ = 0; if(present(use_slk)) use_slk_ = 1
 88 
 89  if( n > eigen_d_maxsize ) then
 90     write(msg,'(a,2i3)')' Eigen size higher than max size set!!',n,eigen_d_maxsize
 91     MSG_ERROR(msg)
 92  endif
 93 
 94 #ifdef HAVE_LINALG_MAGMA
 95  if (usegpu_==1) then
 96     !work only if  lwork=n**2 + 33*n?
 97     if (cplx_ == 2 .and. present(rwork)) then
 98        !old
 99        !lwork_=33*n + (n**2)
100        !lrwork_= 1 + 5*n + 2*(n**2)
101        !liwork_ = 3 + 5*n
102        !new
103        call magmaf_zheevd(jobz,uplo,n,a,lda,w,cwk_(1),-1,rwk_(1),-1,iwk_(1),-1,info)
104        lwork_=int(real(cwk_(1))) ; lrwork_=int(rwk_(1)) ; liwork_=iwk_(1)
105        ABI_ALLOCATE(iwork,(liwork_))
106        call magmaf_zheevd(jobz,uplo,n,a,lda,w,work(1:2*lwork_),lwork_,rwork(1:lrwork_),lrwork_,iwork,liwork_,info)
107        ABI_DEALLOCATE(iwork)
108     else
109        !old
110        !lwork_= 1+n*(6 + 2*n)
111        !lrwork_= 1 + 5*n + 2*(n**2)
112        !liwork_ = 3 + 5*n
113        !new
114        call magmaf_dsyevd(jobz,uplo,n,a,lda,w,rwk_(1),-1,iwk_(1),-1,info)
115        lwork_=int(rwk_(1)) ; liwork_=iwk_(1)
116        ABI_ALLOCATE(iwork,(liwork_))
117        call magmaf_dsyevd(jobz,uplo,n,a,lda,w,work(1:lwork_),lwork_,iwork(1:liwork_),liwork_,info)
118        ABI_DEALLOCATE(iwork)
119     endif
120  else
121 #endif
122 
123 #ifdef HAVE_LINALG_SCALAPACK
124  if( use_slk_ == 1.and.( n > maxval(abi_processor%grid%dims(1:2))) )  then
125     ABI_CHECK(present(x_cplx),"x_cplx must be present")
126     call compute_eigen1(abi_communicator,abi_processor,cplx_,n,n,a,w,istwf_k_)
127     info = 0 ! This is to avoid unwanted warning but it's not clean
128  else
129 #endif
130 
131 #ifdef HAVE_LINALG_PLASMA
132  ! FDahm & LNGuyen  (November 2012) :
133  !  In Plasma v 2.4.6, eigen routines support only
134  !the eigenvalues computation (jobz=N) and not the
135  !full eigenvectors bases determination (jobz=V)
136  if (LSAME(jobz,'N')) then
137     jobz_plasma_a = jobz_plasma(jobz)
138 
139     if ( cplx_ == 2 .and. present(rwork)) then
140        call PLASMA_Alloc_Workspace_zheev(n,n,plasma_work,info)
141        info = PLASMA_zheev_c(jobz_plasma_a,uplo_plasma(uplo),n,c_loc(a),lda,c_loc(w),plasma_work,c_loc(rwork),lwork)
142        ABI_CHECK(info==0,"PLASMA_zheev_c returned info !=0")
143        call PLASMA_Dealloc_handle(plasma_work,info)
144     else
145        call PLASMA_Alloc_Workspace_dsyev(n,n,plasma_work,info)
146        info = PLASMA_dsyev_c(jobz_plasma_a,uplo_plasma(uplo),n,c_loc(a),lda,c_loc(w),plasma_work,c_loc(rwork),lwork)
147        ABI_CHECK(info==0,"PLASMA_dsyev_c returned info !=0")
148        call PLASMA_Dealloc_handle(plasma_work,info)
149     endif
150  else
151 #endif
152    if ( cplx_ == 2 .and. present(rwork)) then
153       call zheev(jobz,uplo,n,a,lda,w,work,lwork/2,rwork,info)
154    else
155 #ifdef FC_NAG
156       ! MG: This hack needed to pass paral[25] paral[29] and mpiio on nag@petrus with np=4
157       if (n < 0) write(std_out, *)"work: ",work(1:3)
158 #endif
159       call dsyev(jobz,uplo,n,a,lda,w,work,lwork,info)
160    endif
161 #ifdef HAVE_LINALG_PLASMA
162  end if
163 #endif
164 #ifdef HAVE_LINALG_SCALAPACK
165  end if
166 #endif
167 #ifdef HAVE_LINALG_MAGMA
168  end if
169 #endif
170 
171  if(info/=0) then
172     write(msg,'(a,i0)')' Problem in abi_xheev, info= ',info
173     !MSG_WARNING(msg)
174     MSG_ERROR(msg)
175  endif
176 
177  if (present(tim_xeigen).and.present(timopt)) then
178    if(abs(timopt)==3) then
179      call timab(tim_xeigen,2,tsec)
180    end if
181  end if
182 
183 end subroutine abi_dheev

m_abi_linalg/abi_dheev_alloc [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

 abi_dheev_alloc

FUNCTION

INPUTS

PARENTS

SOURCE

200   subroutine abi_dheev_alloc(jobz,uplo,n,a,w, &
201 &       x_cplx,istwf_k,timopt,tim_xeigen,use_slk,use_gpu)
202 
203 !This section has been created automatically by the script Abilint (TD).
204 !Do not modify the following lines by hand.
205 #undef ABI_FUNC
206 #define ABI_FUNC 'abi_dheev_alloc'
207 
208 !End of the abilint section
209 
210  implicit none
211  !Arguments ------------------------------------
212  character(len=1), intent(in) :: jobz
213  character(len=1), intent(in) :: uplo
214  integer, intent(in) :: n
215  real(dp), DEV_CONTARRD intent(inout) :: a(:,:)  ! Here I need lda to be consistent
216  real(dp), intent(out) :: w(n)
217  integer, optional, intent(in) :: x_cplx
218  integer, optional, intent(in) :: timopt,tim_xeigen,use_slk,use_gpu
219  integer, optional, intent(in) :: istwf_k
220 
221  integer :: info
222 
223 ! *********************************************************************
224 
225  call abi_dheev(jobz,uplo,n,a,n,w,eigen_d_work,eigen_d_lwork,eigen_z_rwork,info, &
226 &         x_cplx,istwf_k,timopt,tim_xeigen,use_slk,use_gpu)
227 
228 end subroutine abi_dheev_alloc

m_abi_linalg/abi_xheev [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

  abi_xheev

FUNCTION

  abi_xheev is the generic function that compute
  all eigenvalues and, optionally, eigenvectors of a
  symmetric or hermitian matrix A.

COPYRIGHT

  Copyright (C) 2001-2018 ABINIT group (LNguyen,FDahm (CS))
  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


m_abi_linalg/abi_zheev [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

 abi_zheev

FUNCTION

INPUTS

PARENTS

SOURCE

349 subroutine abi_zheev(jobz,uplo,n,a,lda,w,work,lwork,rwork,info)
350 
351 !This section has been created automatically by the script Abilint (TD).
352 !Do not modify the following lines by hand.
353 #undef ABI_FUNC
354 #define ABI_FUNC 'abi_zheev'
355 !End of the abilint section
356 
357  implicit none
358 
359 !Arguments ------------------------------------
360  character(len=1), intent(in) :: jobz
361  character(len=1), intent(in) :: uplo
362  integer, intent(in) :: n
363  integer, intent(in) :: lda
364  integer, intent(in) :: lwork
365  complex(dpc),target,intent(inout) :: a(lda,*)
366  complex(dpc), target, intent(inout) :: work(*)
367  real(dp),target, intent(inout) :: rwork(lwork)        ! TODO Check this!
368  real(dp),target,intent(out) :: w(n)
369  integer, intent(out) :: info
370 
371 #ifdef HAVE_LINALG_PLASMA
372  !Optional Arguments ------------------------------------
373  integer :: jobz_plasma_a
374  type(c_ptr) :: plasma_work
375 #endif
376 
377  ! *********************************************************************
378 
379 #ifdef HAVE_LINALG_PLASMA
380  ! FDahm & LNGuyen  (November 2012) :
381  !  In Plasma v 2.4.6, eigen routines support only
382  !the eigenvalues computation (jobz=N) and not the
383  !full eigenvectors bases determination (jobz=V)
384  if (LSAME(jobz,'N')) then
385     jobz_plasma_a = jobz_plasma(jobz)
386 
387     call PLASMA_Alloc_Workspace_zheev(n,n,plasma_work,info)
388 
389     info = PLASMA_zheev_c(jobz_plasma_a,uplo_plasma(uplo),&
390 &     n,c_loc(a),lda,c_loc(w),plasma_work,c_loc(rwork),lwork)
391 
392     call PLASMA_Dealloc_handle(plasma_work,info)
393  else
394 #endif
395    call zheev(jobz,uplo,n,a,lda,w,work,lwork,rwork,info)
396 #ifdef HAVE_LINALG_PLASMA
397  end if
398 #endif
399 
400 end subroutine abi_zheev

m_abi_linalg/abi_zheev_alloc [ Functions ]

[ Top ] [ m_abi_linalg ] [ Functions ]

NAME

 abi_zheev_alloc

FUNCTION

INPUTS

PARENTS

SOURCE

417 subroutine abi_zheev_alloc(jobz,uplo,n,a,w)
418 
419 !This section has been created automatically by the script Abilint (TD).
420 !Do not modify the following lines by hand.
421 #undef ABI_FUNC
422 #define ABI_FUNC 'abi_zheev_alloc'
423 
424 !End of the abilint section
425 
426  implicit none
427 
428 !Arguments ------------------------------------
429  character(len=1), intent(in) :: jobz
430  character(len=1), intent(in) :: uplo
431  integer, intent(in) :: n
432  complex(dpc), DEV_CONTARRD intent(inout) :: a(:,:)   ! Here I need lda to be consistent
433  real(dp), intent(out) :: w(n)
434 
435  integer :: info
436 
437 ! *********************************************************************
438 
439  call abi_zheev(jobz,uplo,n,a,n,w,eigen_z_work,eigen_z_lwork,eigen_z_rwork,info)
440 
441 end subroutine abi_zheev_alloc