TABLE OF CONTENTS


ABINIT/m_abi2big [ Modules ]

[ Top ] [ Modules ]

NAME

  m_abi2big

FUNCTION

  Module to copy objects from ABINIT to BigDFT and viceversa.

COPYRIGHT

  Copyright (C) 2012-2018 ABINIT group (TR,DC,MT)
  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

CHILDREN

SOURCE

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 module m_abi2big
28 
29  use defs_basis
30  use defs_wvltypes
31  use m_errors
32  use m_xmpi
33 
34  use m_geometry,      only : mkrdim, xcart2xred, xred2xcart
35 
36  implicit none
37 
38  private
39 
40  public :: wvl_vtrial_abi2big
41  !to copy vtrial to wvl_den%rhov and viceversa.
42 
43  public :: wvl_rho_abi2big
44  !to copy a density from ABINIT to BigDFT or viceversa
45 
46  public :: wvl_vhartr_abi2big
47  ! to copy V_hartree from ABINIT to BigDFT and viceversa
48 
49  public :: wvl_vxc_abi2big
50  ! to copy Vxc from ABINIT to BigDFT and viceversa
51 
52  public :: wvl_occ_abi2big
53  ! to copy occupations from/to ABINIT to/from BigDFT
54 
55  public :: wvl_eigen_abi2big
56  ! to copy eigenvalues from/to ABINIT to/from BigDFT
57 
58  public :: wvl_occopt_abi2big
59  ! maps occupation method in ABINIT and BigDFT
60 
61  public :: wvl_rhov_abi2big
62  !generic routine to copy a density or potential from/to
63  !ABINIT to/from BigDFT.
64  interface wvl_rhov_abi2big
65    module procedure wvl_rhov_abi2big_2D_4D
66    module procedure wvl_rhov_abi2big_1D_4D
67    module procedure wvl_rhov_abi2big_2D_2D
68    module procedure wvl_rhov_abi2big_1D_2D
69    module procedure wvl_rhov_abi2big_2D_1D
70    module procedure wvl_rhov_abi2big_1D_1D
71  end interface wvl_rhov_abi2big
72 
73  logical,parameter :: hmem=.false. !high memory
74 !!  Set hmem=.false. if memory is limited. It will copy element by element.
75 !!  If  hmem=.true. all elements are copied at once.
76 
77  public :: wvl_setngfft
78 
79  public :: wvl_setBoxGeometry
80 
81 contains

ABINIT/wvl_setBoxGeometry [ Functions ]

[ Top ] [ Functions ]

NAME

 wvl_setBoxGeometry

FUNCTION

 When wavelets are used, the box definition needs to be changed.
 The box size is recomputed knowing some psp informations such as
 the radius for coarse and fine grid. Then, the atoms are translated
 to be included in the new box. Finally the FFT grid is computed using
 the fine wavelet mesh and a buffer characteristic of used wavelets plus
 a buffer used to be multiple of 2, 3 or 5.

INPUTS

  psps <type(pseudopotential_type)>=variables related to pseudopotentials
  radii= the radii for each type of atoms, giving the fine and the coarse grid.

OUTPUT

  rprimd(3,3)=dimensional primitive translations in real space (bohr)

SIDE EFFECTS

  wvl <type(wvl_internal_type)>=internal variables used by wavelets, describing
                             the box are set.
  xred(3,natom)=reduced dimensionless atomic coordinates

PARENTS

      gstate,wvl_memory,wvl_wfsinp_reformat

CHILDREN

      mkrdim,nullify_locreg_descriptors,system_size,wrtout,xcart2xred
      xred2xcart

SOURCE

1404 subroutine wvl_setBoxGeometry(prtvol, radii, rprimd, xred, wvl, wvl_crmult, wvl_frmult)
1405 
1406 #if defined HAVE_BIGDFT
1407  use BigDFT_API, only: system_size,nullify_locreg_descriptors
1408 #endif
1409 
1410 !This section has been created automatically by the script Abilint (TD).
1411 !Do not modify the following lines by hand.
1412 #undef ABI_FUNC
1413 #define ABI_FUNC 'wvl_setBoxGeometry'
1414 !End of the abilint section
1415 
1416  implicit none
1417 
1418 !Arguments ------------------------------------
1419 !scalars
1420  integer,intent(in) :: prtvol
1421  real(dp), intent(in) :: wvl_crmult, wvl_frmult
1422  type(wvl_internal_type), intent(inout) :: wvl
1423 !arrays
1424  real(dp),intent(in) :: radii(:,:)
1425  real(dp),intent(inout) :: rprimd(3,3),xred(:,:)
1426 
1427 !Local variables-------------------------------
1428 #if defined HAVE_BIGDFT
1429 !scalars
1430  integer :: ii
1431  logical,parameter :: OCLconv=.false.
1432  character(len=500) :: message
1433 !arrays
1434  real(dp) :: rprim(3,3),acell(3)
1435  real(dp),allocatable :: xcart(:,:)
1436 #endif
1437 
1438 ! *********************************************************************
1439 
1440 #if defined HAVE_BIGDFT
1441  if (prtvol == 0) then
1442    write(message, '(a,a,a,a)' ) ch10,&
1443 &   ' wvl_setBoxGeometry : Changing the box for wavelets computation.'
1444    call wrtout(std_out,message,'COLL')
1445  end if
1446 
1447 !Store xcart for each atom
1448  ABI_ALLOCATE(xcart,(3, wvl%atoms%astruct%nat))
1449  call xred2xcart(wvl%atoms%astruct%nat, rprimd, xcart, xred)
1450 
1451  call nullify_locreg_descriptors(wvl%Glr)
1452  call system_size(wvl%atoms, xcart, radii, wvl_crmult, &
1453 & wvl_frmult, wvl%h(1), wvl%h(2), wvl%h(3), OCLconv, wvl%Glr, wvl%shift)
1454 
1455  acell(:) = wvl%atoms%astruct%cell_dim(:)
1456 
1457  if (prtvol == 0) then
1458    write(message, '(a,3F12.6)' ) &
1459 &   '  | acell is now:         ', acell
1460    call wrtout(std_out,message,'COLL')
1461    write(message, '(a,2I5,a,a,2I5,a,a,2I5)' ) &
1462 &   '  | nfl1, nfu1:           ', wvl%Glr%d%nfl1, wvl%Glr%d%nfu1, ch10, &
1463 &   '  | nfl2, nfu2:           ', wvl%Glr%d%nfl2, wvl%Glr%d%nfu2, ch10, &
1464 &   '  | nfl3, nfu3:           ', wvl%Glr%d%nfl3, wvl%Glr%d%nfu3
1465    call wrtout(std_out,message,'COLL')
1466  end if
1467 
1468 !Change the metric to orthogonal one
1469  rprim(:, :) = real(0., dp)
1470  do ii = 1, 3, 1
1471    rprim(ii,ii) = real(1., dp)
1472  end do
1473  call mkrdim(acell, rprim, rprimd)
1474 
1475 !Save shifted atom positions into xred
1476  call xcart2xred(wvl%atoms%astruct%nat, rprimd, xcart, xred)
1477  ABI_DEALLOCATE(xcart)
1478 
1479  if (prtvol == 0) then
1480    write(message, '(a,3I12)' ) &
1481 &   '  | box size for datas:   ', wvl%Glr%d%n1i, wvl%Glr%d%n2i, wvl%Glr%d%n3i
1482    call wrtout(std_out,message,'COLL')
1483    write(message, '(a,3I12)' ) &
1484 &   '  | box size for wavelets:', wvl%Glr%d%n1, wvl%Glr%d%n2, wvl%Glr%d%n3
1485    call wrtout(std_out,message,'COLL')
1486  end if
1487 
1488 #else
1489  BIGDFT_NOTENABLED_ERROR()
1490  if (.false.) write(std_out,*)  prtvol,wvl_crmult,wvl_frmult,wvl%h(1),&
1491 & radii(1,1),rprimd(1,1),xred(1,1)
1492 #endif
1493 
1494 end subroutine wvl_setBoxGeometry

ABINIT/wvl_setngfft [ Functions ]

[ Top ] [ Functions ]

NAME

 wvl_setngfft

FUNCTION

 When wavelets are used, the FFT grid is used to store potentials and
 density. The size of the grid takes into account the two resolution in wavelet
 description and also the distribution over processor in the parallel case.

 The FFT grid is not in strict terms an FFT grid but rather a real space grid.
 Its dimensions are not directly compatible with FFTs. This is not relevant
 when using the wavelet part of the code and in the Poisson solver the arrays
 are extended to match FFT dimensions internally. But for other parts of the
 code, this must be taken into account.

 see doc/variables/vargs.html#ngfft for details about ngfft

SIDE EFFECTS

  mpi_enreg=informations about MPI parallelization (description of the
            density and potentials scatterring is allocated and updated).
  dtset <type(dataset_type)>=the FFT grid is changed.

PARENTS

      gstate,wvl_wfsinp_reformat

CHILDREN

      wrtout

SOURCE

1294 subroutine wvl_setngfft(me_wvl, mgfft, nfft, ngfft, nproc_wvl, n1i, n2i, n3i,n3d)
1295 
1296 
1297 !This section has been created automatically by the script Abilint (TD).
1298 !Do not modify the following lines by hand.
1299 #undef ABI_FUNC
1300 #define ABI_FUNC 'wvl_setngfft'
1301 !End of the abilint section
1302 
1303  implicit none
1304 
1305 !Arguments ------------------------------------
1306 !scalars
1307  integer, intent(out) :: mgfft, nfft
1308  integer, intent(in)  :: n1i, n2i, n3i,n3d, nproc_wvl, me_wvl
1309 !arrays
1310  integer, intent(out) :: ngfft(18)
1311 
1312 !Local variables-------------------------------
1313 !scalars
1314 #if defined HAVE_BIGDFT
1315  character(len=500) :: message
1316 #endif
1317 
1318 ! *************************************************************************
1319 
1320 #if defined HAVE_BIGDFT
1321  write(message, '(a,a,a,a)' ) ch10,&
1322 & ' wvl_setngfft : Changing the FFT grid definition.'
1323  call wrtout(std_out,message,'COLL')
1324 
1325 !Change nfft and ngfft
1326 !Now ngfft will use the density definition (since the potential size
1327 !is always smaller than the density one). ????
1328  ngfft(1) = n1i
1329  ngfft(2) = n2i
1330  ngfft(3) = n3i
1331 
1332  nfft = n1i*n2i*n3d
1333 !Set up fft array dimensions ngfft(4,5,6) to avoid cache conflicts
1334 !Code paste from getng()
1335  ngfft(4) = 2 * (ngfft(1) / 2) + 1
1336  ngfft(5) = 2 * (ngfft(2) / 2) + 1
1337  ngfft(6) = ngfft(3)
1338  if (nproc_wvl == 0) then
1339    ngfft(9)  = 0    ! paral_fft
1340    ngfft(10) = 1    ! nproc_fft
1341    ngfft(11) = 0    ! me_fft
1342    ngfft(12) = 0    ! n2proc
1343    ngfft(13) = 0    ! n3proc
1344  else
1345    ngfft(9)  = 1    ! paral_fft
1346    ngfft(10) = nproc_wvl
1347    ngfft(11) = me_wvl
1348    ngfft(12) = ngfft(2)
1349    ngfft(13) = n3d
1350  end if
1351 
1352  write(message, '(a,3I12)' ) &
1353 & '  | ngfft(1:3) is now:    ', ngfft(1:3)
1354  call wrtout(std_out,message,'COLL')
1355  write(message, '(a,3I12)' ) &
1356 & '  | ngfft(4:6) is now:    ', ngfft(4:6)
1357  call wrtout(std_out,message,'COLL')
1358 
1359 !Set mgfft
1360  mgfft= max(ngfft(1), ngfft(2), ngfft(3))
1361 
1362 #else
1363  BIGDFT_NOTENABLED_ERROR()
1364  if (.false.) write(std_out,*) mgfft,nfft,n1i,n2i,n3i,n3d,nproc_wvl,me_wvl,ngfft(1)
1365 #endif
1366 
1367 end subroutine wvl_setngfft

m_abi2big/wvl_eigen_abi2big [ Functions ]

[ Top ] [ m_abi2big ] [ Functions ]

NAME

  wvl_eigen_abi2big

FUNCTION

  Copies eigenvalues in ABINIT to BigDFT objects and viceversa.

INPUTS

  opt= 1) copy from ABINIT to BigDFT
       2) copy from BigDFT to ABINIT
  nsppol= number of spin polarization

OUTPUT

SIDE EFFECTS

 occ is copied to wfs%ks%orbs%occup, or viceversa, depending on "opt" (see above).

PARENTS

      afterscfloop,vtorho

CHILDREN

SOURCE

547 subroutine wvl_eigen_abi2big(mband,nkpt,nsppol,eigen,opt,wvl_wfs)
548 
549 
550 !This section has been created automatically by the script Abilint (TD).
551 !Do not modify the following lines by hand.
552 #undef ABI_FUNC
553 #define ABI_FUNC 'wvl_eigen_abi2big'
554 !End of the abilint section
555 
556  implicit none
557 
558 !Arguments ------------------------------------
559  integer , intent(in)  :: mband,nkpt,nsppol,opt
560  real(dp) , intent(inout)  :: eigen(mband*nkpt*nsppol)
561  type(wvl_wf_type), intent(inout) :: wvl_wfs
562 
563 !Local variables-------------------------------
564 #if defined HAVE_BIGDFT
565  integer :: ii,norb,norbd,norbu
566  character(len=100) :: message
567 #endif
568 
569 ! *************************************************************************
570 
571  DBG_ENTER("COLL")
572 
573 #if defined HAVE_BIGDFT
574 !PENDING: I am not sure this will work for nsppol==2
575 !check also the parallel case.
576 
577  norbu=wvl_wfs%ks%orbs%norbu
578  norbd=wvl_wfs%ks%orbs%norbd
579  norb =wvl_wfs%ks%orbs%norb
580  if(opt==1) then !ABINIT -> BIGDFT
581    if (nsppol == 1) then
582     wvl_wfs%ks%orbs%eval=eigen
583    else
584      wvl_wfs%ks%orbs%eval(1:norbu)=eigen(1:norbu)
585      wvl_wfs%ks%orbs%eval(norbu + 1:norb)= &
586 &       eigen(mband + 1:mband + norbd)
587    end if
588  elseif(opt==2) then !BigDFT -> ABINIT
589    if (nsppol == 1) then
590      do ii=1,norb
591        eigen(ii)=wvl_wfs%ks%orbs%eval(ii)
592      end do
593    else
594      eigen(1:norbu) = wvl_wfs%ks%orbs%eval(1:norbu)
595      eigen(mband + 1:mband + norbd) = &
596 &     wvl_wfs%ks%orbs%eval(norbu + 1:norb)
597    end if
598  else
599    message='wvl_eigen_abi2big: wrong option'
600    MSG_BUG(message)
601  end if
602 
603 #else
604  BIGDFT_NOTENABLED_ERROR()
605  if (.false.) write(std_out,*) mband,nkpt,nsppol,opt,eigen(1),wvl_wfs%ks
606 #endif
607 
608  DBG_EXIT("COLL")
609 
610 end subroutine wvl_eigen_abi2big

m_abi2big/wvl_occ_abi2big [ Functions ]

[ Top ] [ m_abi2big ] [ Functions ]

NAME

  wvl_occ_abi2big

FUNCTION

  Copies occupations in ABINIT to BigDFT objects and viceversa.

INPUTS

  opt= 1) copy from ABINIT to BigDFT
       2) copy from BigDFT to ABINIT
  nsppol= number of spin polarization

OUTPUT

SIDE EFFECTS

 occ is copied to wfs%ks%orbs%occup, or viceversa, depending on "opt" (see above).

PARENTS

      afterscfloop,gstate,vtorho,wvl_wfsinp_disk,wvl_wfsinp_scratch

CHILDREN

SOURCE

455 subroutine wvl_occ_abi2big(mband,nkpt,nsppol,occ,opt,wvl_wfs)
456 
457 
458 !This section has been created automatically by the script Abilint (TD).
459 !Do not modify the following lines by hand.
460 #undef ABI_FUNC
461 #define ABI_FUNC 'wvl_occ_abi2big'
462 !End of the abilint section
463 
464  implicit none
465 
466 !Arguments ------------------------------------
467  integer , intent(in)  :: mband,nkpt,nsppol,opt
468  real(dp) , intent(inout)  :: occ(mband*nkpt*nsppol)
469  type(wvl_wf_type), intent(inout) :: wvl_wfs
470 
471 !Local variables-------------------------------
472 #if defined HAVE_BIGDFT
473  integer :: norb,norbd,norbu,ii
474  character(len=100) :: message
475 #endif
476 
477 ! *************************************************************************
478 
479  DBG_ENTER("COLL")
480 
481 #if defined HAVE_BIGDFT
482 !PENDING: I am not sure this will work for nsppol==2
483 !check also the parallel case.
484 
485  norbu=wvl_wfs%ks%orbs%norbu
486  norbd=wvl_wfs%ks%orbs%norbd
487  norb =wvl_wfs%ks%orbs%norb
488  if(opt==1) then !ABINIT -> BIGDFT
489    if (nsppol == 1) then
490     do ii=1,norb
491       wvl_wfs%ks%orbs%occup(ii)=occ(ii)
492     end do
493    else
494      wvl_wfs%ks%orbs%occup(1:norbu)=occ(1:norbu)
495      wvl_wfs%ks%orbs%occup(norbu + 1:norb)= &
496 &       occ(mband + 1:mband + norbd)
497    end if
498  elseif(opt==2) then !BigDFT -> ABINIT
499    if (nsppol == 1) then
500      do ii=1,norb
501        occ=wvl_wfs%ks%orbs%occup
502      end do
503    else
504      occ(1:norbu) = wvl_wfs%ks%orbs%occup(1:norbu)
505      occ(mband + 1:mband + norbd) = &
506 &     wvl_wfs%ks%orbs%occup(norbu + 1:norb)
507    end if
508  else
509    message='wvl_occ_abi2big: wrong option'
510    MSG_BUG(message)
511  end if
512 
513 #else
514  BIGDFT_NOTENABLED_ERROR()
515  if (.false.) write(std_out,*) mband,nkpt,nsppol,opt,occ(1),wvl_wfs%ks
516 #endif
517 
518  DBG_EXIT("COLL")
519 
520 end subroutine wvl_occ_abi2big

m_abi2big/wvl_occopt_abi2big [ Functions ]

[ Top ] [ m_abi2big ] [ Functions ]

NAME

  wvl_occopt_abi2big

FUNCTION

  Copies occopt in ABINIT to BigDFT objects and viceversa.

INPUTS

  opt= 1) copy from ABINIT to BigDFT
       2) copy from BigDFT to ABINIT

OUTPUT

SIDE EFFECTS

NOTES

 Several smearing schemes do not exists in both codes such
 as the SMEARING_DIST_ERF in BigDFT.

PARENTS

      vtorho,wvl_wfsinp_scratch

CHILDREN

SOURCE

639 subroutine wvl_occopt_abi2big(occopt_abi,occopt_big,opt)
640 
641 #if defined HAVE_BIGDFT
642  use BigDFT_API, only : &
643 &  SMEARING_DIST_FERMI, SMEARING_DIST_COLD1, SMEARING_DIST_COLD2,&
644 &  SMEARING_DIST_METPX
645 #endif
646 
647 !This section has been created automatically by the script Abilint (TD).
648 !Do not modify the following lines by hand.
649 #undef ABI_FUNC
650 #define ABI_FUNC 'wvl_occopt_abi2big'
651 !End of the abilint section
652 
653  implicit none
654 
655 !Arguments ------------------------------------
656  integer , intent(inout)  :: occopt_abi,occopt_big
657  integer , intent(in)     :: opt
658 
659 !Local variables-------------------------------
660 #if defined HAVE_BIGDFT
661  character(len=500) :: message
662 #endif
663 
664 ! *************************************************************************
665 
666  DBG_ENTER("COLL")
667 
668 #if defined HAVE_BIGDFT
669 
670  if(opt==1) then !ABINIT -> BIGDFT
671    if(occopt_abi==3) then
672      occopt_big=SMEARING_DIST_FERMI
673    elseif(occopt_abi==4) then
674      occopt_big=SMEARING_DIST_COLD1
675    elseif(occopt_abi==5) then
676      occopt_big=SMEARING_DIST_COLD2
677    elseif(occopt_abi==6) then
678      occopt_big=SMEARING_DIST_METPX
679    else
680      write(message,'(4a)') ch10,&
681 &     ' wvl_occopt_abi2big: occopt does not have a corresponding option in BigDFT.',ch10,&
682 &     ' Action: change the value of occopt to a number between 3 and 6'
683      MSG_ERROR(message)
684    end if
685  elseif(opt==2) then !BigDFT -> ABINIT
686    if(occopt_big==SMEARING_DIST_FERMI) then
687      occopt_abi=3
688    elseif(occopt_big==SMEARING_DIST_COLD1) then
689      occopt_abi=4
690    elseif(occopt_big==SMEARING_DIST_COLD2) then
691      occopt_abi=5
692    elseif(occopt_big==SMEARING_DIST_METPX) then
693      occopt_abi=6
694    else
695 !    One should never get here.
696      write(message,'(4a)') ch10,&
697 &     ' wvl_occopt_abi2big: occopt in BigDFT does not have a corresponding option in ABINIT.',ch10,&
698 &     ' Action: contact the ABINIT group'
699      MSG_ERROR(message)
700    end if
701  else
702    message='wvl_occopt_abi2big: wrong option'
703    MSG_BUG(message)
704  end if
705 
706 #else
707  BIGDFT_NOTENABLED_ERROR()
708  if (.false.) write(std_out,*) occopt_abi,occopt_big,opt
709 #endif
710 
711  DBG_EXIT("COLL")
712 
713 end subroutine wvl_occopt_abi2big

m_abi2big/wvl_rho_abi2big [ Functions ]

[ Top ] [ m_abi2big ] [ Functions ]

NAME

  wvl_rho_abi2big

FUNCTION

  Copies the density from ABINIT to BigDFT, or viceversa.

INPUTS

  opt= 1) copy from ABINIT to BigDFT
       2) copy from BigDFT to ABINIT
  rhor(nfft,nspden)= trial potential (ABINIT array)
  wvl_den= density-potential BigDFT object

OUTPUT

SIDE EFFECTS

 Density copied from ABINIT to BigDFT or viceversa.
 It verifies that (or sets) wvl_den%rhov_is= ELECTRONIC_DENSITY.

NOTES

 It uses the generic routine wvl_rhov_abi2big.

PARENTS

      afterscfloop,mklocl,newrho,vtorho,wvl_mkrho

CHILDREN

SOURCE

207 subroutine wvl_rho_abi2big(opt,rhor,wvl_den)
208 
209 #if defined HAVE_BIGDFT
210   use BigDFT_API, only : ELECTRONIC_DENSITY
211 #endif
212 
213 !This section has been created automatically by the script Abilint (TD).
214 !Do not modify the following lines by hand.
215 #undef ABI_FUNC
216 #define ABI_FUNC 'wvl_rho_abi2big'
217 !End of the abilint section
218 
219  implicit none
220 
221 !Arguments ------------------------------------
222  integer , intent(in)  :: opt
223  real(dp) , intent(inout)  :: rhor(:,:)
224  type(wvl_denspot_type), intent(inout) :: wvl_den
225 
226 !Local variables-------------------------------
227 #if defined HAVE_BIGDFT
228  character(len=100) :: message
229 #endif
230 
231 ! *************************************************************************
232 
233  DBG_ENTER("COLL")
234 
235 #if defined HAVE_BIGDFT
236 ! write(message,'(2a)') ch10,'wvl_rho_abi2big: but why are you copying me :..o('
237 ! call wrtout(std_out,message,'COLL')
238 
239  if(opt==1) then !ABINIT -> BIGDFT
240 
241    call wvl_rhov_abi2big(opt,rhor,wvl_den%denspot%rhov)
242    wvl_den%denspot%rhov_is = ELECTRONIC_DENSITY
243 
244  elseif(opt==2) then !BigDFT -> ABINIT
245 
246    if(wvl_den%denspot%rhov_is .ne. ELECTRONIC_DENSITY) then
247      message='wvl_rho_abi2big: rhov should contain the ELECTRONIC_DENSITY'
248      MSG_BUG(message)
249    end if
250    call wvl_rhov_abi2big(opt,rhor,wvl_den%denspot%rhov)
251 
252  else
253    message='wvl_rho_abi2big: wrong option'
254    MSG_BUG(message)
255  end if
256 
257 #else
258  BIGDFT_NOTENABLED_ERROR()
259  if (.false.) write(std_out,*)  opt,rhor(1,1),wvl_den%symObj
260 #endif
261 
262  DBG_EXIT("COLL")
263 
264 end subroutine wvl_rho_abi2big

m_abi2big/wvl_rhov_abi2big_1D_1D [ Functions ]

[ Top ] [ m_abi2big ] [ Functions ]

NAME

  wvl_rhov_abi2big_1D_1D

FUNCTION

  Copies a density/potential from ABINIT to BigDFT or viceversa.
  Target : ABINIT 2D arrays (without spin), BigDFT 1D arrays (with or without spin)

INPUTS

  opt= 1: copy from ABINIT to BigDFT, 2: copy from BigDFT to ABINIT
  rhov_abi(:) = density/potential array in ABINIT
  rhov_big(:) = density/potential array in BigDFT
  [shift] = shift to be applied in rhov_abi array (parallelism)

OUTPUT

SIDE EFFECTS

  At output rhov_abi is copied to rhov_abinit, or viceversa

PARENTS

CHILDREN

SOURCE

1064 subroutine wvl_rhov_abi2big_1D_1D(opt,rhov_abi,rhov_big,shift)
1065 
1066 
1067 !This section has been created automatically by the script Abilint (TD).
1068 !Do not modify the following lines by hand.
1069 #undef ABI_FUNC
1070 #define ABI_FUNC 'wvl_rhov_abi2big_1D_1D'
1071 !End of the abilint section
1072 
1073  implicit none
1074 
1075 !Arguments ------------------------------------
1076  integer,intent(in) :: opt
1077  integer,intent(in),optional :: shift
1078  real(dp) :: rhov_abi(:),rhov_big(:)
1079 
1080 !Local variables-------------------------------
1081 #if defined HAVE_BIGDFT
1082  integer :: nfft_abi,nfft_big,shift_
1083 #endif
1084 
1085 ! *************************************************************************
1086 
1087 #if defined HAVE_BIGDFT
1088  nfft_abi=size(rhov_abi) ; nfft_big=size(rhov_big)
1089  shift_=0;if (present(shift)) shift_=shift
1090  call wvl_rhov_abi2big_gen(nfft_abi,nfft_big,1,opt,rhov_abi,rhov_big,shift_)
1091 #else
1092  BIGDFT_NOTENABLED_ERROR()
1093  if (.false. .and. present(shift)) write(std_out,*)  opt,rhov_abi(1),rhov_big(1)
1094 #endif
1095 
1096 end subroutine wvl_rhov_abi2big_1D_1D

m_abi2big/wvl_rhov_abi2big_1D_2D [ Functions ]

[ Top ] [ m_abi2big ] [ Functions ]

NAME

  wvl_rhov_abi2big_1D_2D

FUNCTION

  Copies a density/potential from ABINIT to BigDFT or viceversa.
  Target : ABINIT 1D arrays (without spin), BigDFT 2D arrays (with spin)

INPUTS

  opt= 1: copy from ABINIT to BigDFT, 2: copy from BigDFT to ABINIT
  rhov_abi(:)   = density/potential array in ABINIT
  rhov_big(:,:) = density/potential array in BigDFT
  [shift] = shift to be applied in rhov_abi array (parallelism)

OUTPUT

SIDE EFFECTS

  At output rhov_abi is copied to rhov_abinit, or viceversa

PARENTS

CHILDREN

SOURCE

938 subroutine wvl_rhov_abi2big_1D_2D(opt,rhov_abi,rhov_big,shift)
939 
940 
941 !This section has been created automatically by the script Abilint (TD).
942 !Do not modify the following lines by hand.
943 #undef ABI_FUNC
944 #define ABI_FUNC 'wvl_rhov_abi2big_1D_2D'
945 !End of the abilint section
946 
947  implicit none
948 
949 !Arguments ------------------------------------
950  integer,intent(in) :: opt
951  integer,intent(in),optional :: shift
952  real(dp) :: rhov_abi(:),rhov_big(:,:)
953 
954 !Local variables-------------------------------
955 #if defined HAVE_BIGDFT
956  integer :: nfft_abi,nfft_big,shift_
957  character(len=100) :: message
958 #endif
959 
960 ! *************************************************************************
961 
962 #if defined HAVE_BIGDFT
963  if (size(rhov_big,2)/=1) then
964    message='wvl_rhov_abi2big: ABINIT and BigDFT objects do not have the same nspden'
965    MSG_BUG(message)
966  end if
967  nfft_abi=size(rhov_abi) ; nfft_big=size(rhov_big)
968  shift_=0;if (present(shift)) shift_=shift
969  call wvl_rhov_abi2big_gen(nfft_abi,nfft_big,1,opt,rhov_abi,rhov_big,shift_)
970 #else
971  BIGDFT_NOTENABLED_ERROR()
972  if (.false. .and. present(shift)) write(std_out,*)  opt,rhov_abi(1),rhov_big(1,1)
973 #endif
974 
975 end subroutine wvl_rhov_abi2big_1D_2D

m_abi2big/wvl_rhov_abi2big_1D_4D [ Functions ]

[ Top ] [ m_abi2big ] [ Functions ]

NAME

  wvl_rhov_abi2big_1D_4D

FUNCTION

  Copies a density/potential from ABINIT to BigDFT or viceversa.
  Target : ABINIT 1D arrays (without spin), BigDFT 4D arrays (with spin)

INPUTS

  opt= 1: copy from ABINIT to BigDFT, 2: copy from BigDFT to ABINIT
  rhov_abi(:)       = density/potential array in ABINIT
  rhov_big(:,:,:,:) = density/potential array in BigDFT
  [shift] = shift to be applied in rhov_abi array (parallelism)

OUTPUT

SIDE EFFECTS

  At output rhov_abi is copied to rhov_abinit, or viceversa

PARENTS

CHILDREN

SOURCE

807 subroutine wvl_rhov_abi2big_1D_4D(opt,rhov_abi,rhov_big,shift)
808 
809 
810 !This section has been created automatically by the script Abilint (TD).
811 !Do not modify the following lines by hand.
812 #undef ABI_FUNC
813 #define ABI_FUNC 'wvl_rhov_abi2big_1D_4D'
814 !End of the abilint section
815 
816  implicit none
817 
818 !Arguments ------------------------------------
819  integer,intent(in) :: opt
820  integer,intent(in),optional :: shift
821  real(dp) :: rhov_abi(:),rhov_big(:,:,:,:)
822 
823 !Local variables-------------------------------
824 #if defined HAVE_BIGDFT
825  integer :: nfft_abi,nfft_big,shift_
826  character(len=100) :: message
827 #endif
828 
829 ! *************************************************************************
830 
831 #if defined HAVE_BIGDFT
832  if (size(rhov_big,4)/=1) then
833    message='wvl_rhov_abi2big: ABINIT and BigDFT objects do not have the same nspden'
834    MSG_BUG(message)
835  end if
836  nfft_abi=size(rhov_abi) ; nfft_big=size(rhov_big)
837  shift_=0;if (present(shift)) shift_=shift
838  call wvl_rhov_abi2big_gen(nfft_abi,nfft_big,1,opt,rhov_abi,rhov_big,shift_)
839 #else
840  BIGDFT_NOTENABLED_ERROR()
841  if (.false. .and. present(shift)) write(std_out,*)  opt,rhov_abi(1),rhov_big(1,1,1,1)
842 #endif
843 
844 end subroutine wvl_rhov_abi2big_1D_4D

m_abi2big/wvl_rhov_abi2big_2D_1D [ Functions ]

[ Top ] [ m_abi2big ] [ Functions ]

NAME

  wvl_rhov_abi2big_2D_1D

FUNCTION

  Copies a density/potential from ABINIT to BigDFT or viceversa.
  Target : ABINIT 2D arrays (with spin), BigDFT 1D arrays (with or without spin)

INPUTS

  opt= 1: copy from ABINIT to BigDFT, 2: copy from BigDFT to ABINIT
  rhov_abi(:,:) = density/potential array in ABINIT
  rhov_big(:)   = density/potential array in BigDFT
  [shift] = shift to be applied in rhov_abi array (parallelism)

OUTPUT

SIDE EFFECTS

  At output rhov_abi is copied to rhov_abinit, or viceversa

PARENTS

CHILDREN

SOURCE

1003 subroutine wvl_rhov_abi2big_2D_1D(opt,rhov_abi,rhov_big,shift)
1004 
1005 
1006 !This section has been created automatically by the script Abilint (TD).
1007 !Do not modify the following lines by hand.
1008 #undef ABI_FUNC
1009 #define ABI_FUNC 'wvl_rhov_abi2big_2D_1D'
1010 !End of the abilint section
1011 
1012  implicit none
1013 
1014 !Arguments ------------------------------------
1015  integer,intent(in) :: opt
1016  integer,intent(in),optional :: shift
1017  real(dp) :: rhov_abi(:,:),rhov_big(:)
1018 
1019 !Local variables-------------------------------
1020 #if defined HAVE_BIGDFT
1021  integer :: nfft_abi,nfft_big,nspden,shift_
1022 #endif
1023 
1024 ! *************************************************************************
1025 
1026 #if defined HAVE_BIGDFT
1027  nspden=size(rhov_abi,2)
1028  nfft_abi=size(rhov_abi)/nspden ; nfft_big=size(rhov_big)/nspden
1029  shift_=0;if (present(shift)) shift_=shift
1030  call wvl_rhov_abi2big_gen(nfft_abi,nfft_big,nspden,opt,rhov_abi,rhov_big,shift_)
1031 #else
1032  BIGDFT_NOTENABLED_ERROR()
1033  if (.false. .and. present(shift)) write(std_out,*)  opt,rhov_abi(1,1),rhov_big(1)
1034 #endif
1035 
1036 end subroutine wvl_rhov_abi2big_2D_1D

m_abi2big/wvl_rhov_abi2big_2D_2D [ Functions ]

[ Top ] [ m_abi2big ] [ Functions ]

NAME

  wvl_rhov_abi2big_2D_2D

FUNCTION

  Copies a density/potential from ABINIT to BigDFT or viceversa.
  Target : ABINIT 2D arrays (with spin), BigDFT 2D arrays (with spin)

INPUTS

  opt= 1: copy from ABINIT to BigDFT, 2: copy from BigDFT to ABINIT
  rhov_abi(:,:) = density/potential array in ABINIT
  rhov_big(:,:) = density/potential array in BigDFT
  [shift] = shift to be applied in rhov_abi array (parallelism)

OUTPUT

SIDE EFFECTS

  At output rhov_abi is copied to rhov_abinit, or viceversa

PARENTS

CHILDREN

SOURCE

872 subroutine wvl_rhov_abi2big_2D_2D(opt,rhov_abi,rhov_big,shift)
873 
874 
875 !This section has been created automatically by the script Abilint (TD).
876 !Do not modify the following lines by hand.
877 #undef ABI_FUNC
878 #define ABI_FUNC 'wvl_rhov_abi2big_2D_2D'
879 !End of the abilint section
880 
881  implicit none
882 
883 !Arguments ------------------------------------
884  integer,intent(in) :: opt
885  integer,intent(in),optional :: shift
886  real(dp) :: rhov_abi(:,:),rhov_big(:,:)
887 
888 !Local variables-------------------------------
889 #if defined HAVE_BIGDFT
890  integer :: nfft_abi,nfft_big,nspden,shift_
891  character(len=100) :: message
892 #endif
893 
894 ! *************************************************************************
895 
896 #if defined HAVE_BIGDFT
897  nspden=size(rhov_abi,2)
898  if (size(rhov_big,2)/=nspden) then
899    message='wvl_rhov_abi2big: ABINIT and BigDFT objects do not have the same nspden'
900    MSG_BUG(message)
901  end if
902  nfft_abi=size(rhov_abi)/nspden ; nfft_big=size(rhov_big)/nspden
903  shift_=0;if (present(shift)) shift_=shift
904  call wvl_rhov_abi2big_gen(nfft_abi,nfft_big,nspden,opt,rhov_abi,rhov_big,shift_)
905 #else
906  BIGDFT_NOTENABLED_ERROR()
907  if (.false. .and. present(shift)) write(std_out,*)  opt,rhov_abi(1,1),rhov_big(1,1)
908 #endif
909 
910 end subroutine wvl_rhov_abi2big_2D_2D

m_abi2big/wvl_rhov_abi2big_2D_4D [ Functions ]

[ Top ] [ m_abi2big ] [ Functions ]

NAME

  wvl_rhov_abi2big_2D_4D

FUNCTION

  Copies a density/potential from ABINIT to BigDFT or viceversa.
  Target : ABINIT 2D arrays (with spin), BigDFT 4D arrays (with spin)

INPUTS

  opt= 1: copy from ABINIT to BigDFT, 2: copy from BigDFT to ABINIT
  rhov_abi(:,:)     = density/potential array in ABINIT
  rhov_big(:,:,:,:) = density/potential array in BigDFT
  [shift] = shift to be applied in rhov_abi array (parallelism)

OUTPUT

SIDE EFFECTS

  At output rhov_abi is copied to rhov_abinit, or viceversa

PARENTS

CHILDREN

SOURCE

741 subroutine wvl_rhov_abi2big_2D_4D(opt,rhov_abi,rhov_big,shift)
742 
743 
744 !This section has been created automatically by the script Abilint (TD).
745 !Do not modify the following lines by hand.
746 #undef ABI_FUNC
747 #define ABI_FUNC 'wvl_rhov_abi2big_2D_4D'
748 !End of the abilint section
749 
750  implicit none
751 
752 !Arguments ------------------------------------
753  integer,intent(in) :: opt
754  integer,intent(in),optional :: shift
755  real(dp) :: rhov_abi(:,:),rhov_big(:,:,:,:)
756 
757 !Local variables-------------------------------
758 #if defined HAVE_BIGDFT
759  integer :: nfft_abi,nfft_big,nspden,shift_
760  character(len=100) :: message
761 #endif
762 
763 ! *************************************************************************
764 
765 #if defined HAVE_BIGDFT
766  nspden=size(rhov_abi,2)
767  if (size(rhov_big,4)/=nspden) then
768    message='wvl_rhov_abi2big: ABINIT and BigDFT objects do not have the same nspden'
769    MSG_BUG(message)
770  end if
771  nfft_abi=size(rhov_abi)/nspden ; nfft_big=size(rhov_big)/nspden
772  shift_=0;if (present(shift)) shift_=shift
773  call wvl_rhov_abi2big_gen(nfft_abi,nfft_big,nspden,opt,rhov_abi,rhov_big,shift_)
774 #else
775  BIGDFT_NOTENABLED_ERROR()
776  if (.false. .and. present(shift)) write(std_out,*)  opt,rhov_abi(1,1),rhov_big(1,1,1,1)
777 #endif
778 
779 end subroutine wvl_rhov_abi2big_2D_4D

m_abi2big/wvl_rhov_abi2big_gen [ Functions ]

[ Top ] [ m_abi2big ] [ Functions ]

NAME

  wvl_rhov_abi2big_gen

FUNCTION

  Copies a density/potential from ABINIT to BigDFT or viceversa.
  This is a generic routine to copy objects.

INPUTS

  nfft_abi = size of rhov_abi
  nfft_big = size of rhov_big
  nspden =  number of spin components
  opt= 1) copy from ABINIT to BigDFT
       2) copy from BigDFT to ABINIT
  rhov_abi(nfft_abi,nspden) = density/potential array in ABINIT
  rhov_big(nfft_big,nspden) = density/potential array in BigDFT
  shift = shift to be applied in rhov_abi array (parallelism)

OUTPUT

SIDE EFFECTS

 At output rhov_abi is copied to rhov_abinit, or viceversa

NOTES

 This routine is duplicated:
 This option is faster but it requires more memory.
 Notice that we cannot point the variables since the spin convention is not
 the same in BigDFT and ABINIT.
 In ABINIT: index 1 is for the total spin (spin up + spin down) and index 2 is for spin up.
 In BigDFT: indices 1 and 2 are for spin up and down, respectively.

PARENTS

      m_abi2big

CHILDREN

SOURCE

1137 subroutine wvl_rhov_abi2big_gen(nfft_abi,nfft_big,nspden,opt,rhov_abi,rhov_big,shift)
1138 
1139 
1140 !This section has been created automatically by the script Abilint (TD).
1141 !Do not modify the following lines by hand.
1142 #undef ABI_FUNC
1143 #define ABI_FUNC 'wvl_rhov_abi2big_gen'
1144 !End of the abilint section
1145 
1146  implicit none
1147 
1148 !Arguments ------------------------------------
1149  integer,intent(in) :: nfft_abi,nfft_big,nspden,opt,shift
1150  real(dp) :: rhov_abi(nfft_abi,nspden),rhov_big(nfft_big,nspden)
1151 
1152 !Local variables-------------------------------
1153 #if defined HAVE_BIGDFT
1154  integer :: ifft,jfft,nfft
1155  real(dp) :: tmpUp,tmpDown,tmpTot
1156  character(len=100) :: message
1157  real(dp),allocatable :: rhoup(:),rhodn(:),rhotot(:)
1158 #endif
1159 
1160 ! *************************************************************************
1161 
1162  DBG_ENTER("COLL")
1163 
1164 #if defined HAVE_BIGDFT
1165 !No objects to copy; in BigDFT by default they have size of 1!
1166  if(size(rhov_big)==1.and.size(rhov_abi)==0) return
1167 
1168  nfft=nfft_big;if (nfft_big+shift>nfft_abi) nfft=nfft-shift
1169 
1170  if (nfft_abi<nfft+shift) then
1171    message='wvl_rhov_abi2big: cannot handle nfft(abi)<nfft(big)+shift case'
1172    MSG_BUG(message)
1173  end if
1174  if(nspden==4) then
1175    message='wvl_rhov_abi2big: nspden=4 not yet supported'
1176    MSG_ERROR(message)
1177  end if
1178 
1179  if (hmem.and.nspden==2) then
1180    if (opt==1) then
1181      ABI_ALLOCATE(rhoup,(nfft))
1182      ABI_ALLOCATE(rhodn,(nfft))
1183    else if (opt==2)  then
1184      ABI_ALLOCATE(rhotot,(nfft))
1185    end if
1186  end if
1187 
1188  if (opt==1) then !ABINIT -> BIGDFT
1189    if (nspden==2) then
1190      if (hmem) then
1191        rhoup(1:nfft)=rhov_abi(shift+1:shift+nfft,2)
1192        rhodn(1:nfft)=rhov_abi(shift+1:shift+nfft,1)-rhoup(1:nfft)
1193        rhov_big(:,1)=rhoup(:)
1194        rhov_big(:,2)=rhodn(:)
1195      else
1196        do ifft=1,nfft
1197          jfft=shift+ifft
1198 !        We change convention for BigDFT
1199          tmpDown=rhov_abi(jfft,1)-rhov_abi(jfft,2)
1200          tmpUp  =rhov_abi(jfft,2)
1201          rhov_big(ifft,1)=tmpUp
1202          rhov_big(ifft,2)=tmpDown
1203        end do
1204      end if !hmem
1205    else !nspden==1
1206      if (hmem) then
1207        rhov_big(1:nfft,1)=rhov_abi(shift+1:shift+nfft,1)
1208      else
1209        do ifft=1,nfft
1210          rhov_big(ifft,1)=rhov_abi(shift+ifft,1)
1211        end do
1212      end if!hmem
1213    end if !nspden
1214 
1215  else if (opt==2) then !BigDFT -> ABINIT
1216    if (nspden==2) then
1217      if (hmem) then
1218        rhotot(:)=rhov_big(:,1)+rhov_big(:,2)
1219        rhov_abi(shift+1:shift+nfft,1)=rhotot(1:nfft)
1220        rhov_abi(shift+1:shift+nfft,2)=rhov_big(1:nfft,1)
1221      else
1222        do ifft=1,nfft
1223          jfft=shift+ifft
1224 !        We change convention for BigDFT
1225          tmpTot=rhov_big(ifft,1)+rhov_big(ifft,2)
1226          rhov_abi(jfft,1)=tmpTot
1227          rhov_abi(jfft,2)=rhov_big(ifft,1) !Spin Up
1228        end do
1229      end if !hmem
1230    else if (nspden==1) then
1231      if (hmem) then
1232        rhov_abi(shift+1:shift+nfft,1)=rhov_big(1:nfft,1)
1233      else
1234        do ifft=1,nfft
1235          rhov_abi(shift+ifft,1)=rhov_big(ifft,1)
1236        end do
1237      end if !hmem
1238    end if !nspden
1239 
1240  else
1241    message='wvl_rhov_abi2big_gen: wrong option'
1242    MSG_BUG(message)
1243  end if
1244 
1245  if (hmem.and.nspden==2) then
1246    if (opt==1) then
1247      ABI_DEALLOCATE(rhoup)
1248      ABI_DEALLOCATE(rhodn)
1249    else if (opt==2) then
1250      ABI_DEALLOCATE(rhotot)
1251    end if !opt
1252  end if
1253 
1254 #else
1255  BIGDFT_NOTENABLED_ERROR()
1256  if (.false.) write(std_out,*)  nfft_abi,nfft_big,nspden,opt,shift,rhov_big(1,1),rhov_abi(1,1)
1257 #endif
1258 
1259  DBG_EXIT("COLL")
1260 
1261 end subroutine wvl_rhov_abi2big_gen

m_abi2big/wvl_vhartr_abi2big [ Functions ]

[ Top ] [ m_abi2big ] [ Functions ]

NAME

  wvl_vhartr_abi2big

FUNCTION

  Copies vhartree in ABINIT to BigDFT objects and viceversa.

INPUTS

  opt= 1) copy from ABINIT to BigDFT
       2) copy from BigDFT to ABINIT
  vhartr(nfft)= Hartree potential (ABINIT array)
  wvl_den= density-potential BigDFT object

OUTPUT

SIDE EFFECTS

 vhartr is copied to wvl_den, or viceversa, depending on "opt" (see above).
 It verifies that (or sets) wvl_den%rhov_is = HARTREE_POTENTIAL

NOTES

 It uses the generic routine wvl_rhov_abi2big.

PARENTS

CHILDREN

SOURCE

295 subroutine wvl_vhartr_abi2big(opt,vhartr,wvl_den)
296 
297 #if defined HAVE_BIGDFT
298   use BigDFT_API, only : HARTREE_POTENTIAL
299 #endif
300 
301 !This section has been created automatically by the script Abilint (TD).
302 !Do not modify the following lines by hand.
303 #undef ABI_FUNC
304 #define ABI_FUNC 'wvl_vhartr_abi2big'
305 !End of the abilint section
306 
307  implicit none
308 
309 !Arguments ------------------------------------
310  integer , intent(in)  :: opt
311  real(dp) , intent(inout)  :: vhartr(:)
312  type(wvl_denspot_type), intent(inout) :: wvl_den
313 
314 !Local variables-------------------------------
315 #if defined HAVE_BIGDFT
316  integer :: shiftV
317  character(len=100) :: message
318 #endif
319 
320 ! *************************************************************************
321 
322  DBG_ENTER("COLL")
323 
324 #if defined HAVE_BIGDFT
325 ! write(message,'(2a)') ch10, 'wvl_vhartr_abi2big: but why are you copying me :..o('
326 ! call wrtout(std_out,message,'COLL')
327 
328  shiftV=wvl_den%denspot%dpbox%ndims(1)*wvl_den%denspot%dpbox%ndims(2) &
329 &      *wvl_den%denspot%dpbox%i3xcsh
330 
331  if(opt==1) then !ABINIT -> BIGDFT
332 
333    call wvl_rhov_abi2big(opt,vhartr,wvl_den%denspot%rhov,shift=shiftV)
334    wvl_den%denspot%rhov_is = HARTREE_POTENTIAL
335 
336  elseif(opt==2) then !BigDFT -> ABINIT
337 
338    if(wvl_den%denspot%rhov_is .ne. HARTREE_POTENTIAL) then
339      message='wvl_vhartr_abi2big: rhov should contain the HARTREE_POTENTIAL'
340      MSG_BUG(message)
341    end if
342    call wvl_rhov_abi2big(opt,vhartr,wvl_den%denspot%rhov,shift=shiftV)
343 
344  else
345    message='wvl_vhartr_abi2big: wrong option'
346    MSG_BUG(message)
347  end if
348 
349 #else
350  BIGDFT_NOTENABLED_ERROR()
351  if (.false.) write(std_out,*)  opt,vhartr(1),wvl_den%symObj
352 #endif
353 
354  DBG_EXIT("COLL")
355 
356 end subroutine wvl_vhartr_abi2big

m_abi2big/wvl_vtrial_abi2big [ Functions ]

[ Top ] [ m_abi2big ] [ Functions ]

NAME

  wvl_vtrial_abi2big

FUNCTION

  Copies vtrial in ABINIT to BigDFT objects and viceversa.

INPUTS

  opt= 1) copy from ABINIT to BigDFT
       2) copy from BigDFT to ABINIT
  vtrial(nfft,nspden)= trial potential (ABINIT array)
  wvl_den= density-potential BigDFT object

OUTPUT

SIDE EFFECTS

 vtrial is copied to wvl_den, or viceversa, depending on "opt" (see above).
 It verifies that (or sets) wvl_den%rhov_is = KS_POTENTIAL.

NOTES

 It uses the generic routine wvl_rhov_abi2big.

PARENTS

      afterscfloop,newvtr,rhotov,setvtr,wvl_psitohpsi

CHILDREN

SOURCE

113 subroutine wvl_vtrial_abi2big(opt,vtrial,wvl_den)
114 
115 #if defined HAVE_BIGDFT
116   use BigDFT_API, only : KS_POTENTIAL
117 #endif
118 
119 !This section has been created automatically by the script Abilint (TD).
120 !Do not modify the following lines by hand.
121 #undef ABI_FUNC
122 #define ABI_FUNC 'wvl_vtrial_abi2big'
123 !End of the abilint section
124 
125  implicit none
126 
127 !Arguments ------------------------------------
128  integer,intent(in) :: opt
129  real(dp),intent(inout) :: vtrial(:,:)
130  type(wvl_denspot_type),intent(inout) :: wvl_den
131 
132 !Local variables-------------------------------
133 #if defined HAVE_BIGDFT
134  integer :: shiftV
135  character(len=100) :: message
136 #endif
137 
138 ! *************************************************************************
139 
140  DBG_ENTER("COLL")
141 
142 #if defined HAVE_BIGDFT
143 ! write(message,'(2a)') ch10,' wvl_vtrial_abi2big: but why are you copying me :..o('
144 ! call wrtout(std_out,message,'COLL')
145 
146  shiftV=wvl_den%denspot%dpbox%ndims(1)*wvl_den%denspot%dpbox%ndims(2) &
147 &      *wvl_den%denspot%dpbox%i3xcsh
148 
149  if(opt==1) then !ABINIT -> BIGDFT
150 
151    call wvl_rhov_abi2big(opt,vtrial,wvl_den%denspot%rhov,shift=shiftV)
152    wvl_den%denspot%rhov_is = KS_POTENTIAL
153 
154  elseif(opt==2) then !BigDFT -> ABINIT
155 
156    if(wvl_den%denspot%rhov_is .ne. KS_POTENTIAL) then
157      message='wvl_vtrial_abi2big: rhov should contain the KS_POTENTIAL'
158      MSG_BUG(message)
159    end if
160 
161    call wvl_rhov_abi2big(opt,vtrial,wvl_den%denspot%rhov,shift=shiftV)
162 
163  else
164    message='wvl_vtrial_abi2big: wrong option'
165    MSG_BUG(message)
166  end if
167 
168 #else
169  BIGDFT_NOTENABLED_ERROR()
170  if (.false.) write(std_out,*)  opt,vtrial(1,1),wvl_den%symObj
171 #endif
172 
173  DBG_EXIT("COLL")
174 
175 end subroutine wvl_vtrial_abi2big

m_abi2big/wvl_vxc_abi2big [ Functions ]

[ Top ] [ m_abi2big ] [ Functions ]

NAME

  wvl_vxc_abi2big

FUNCTION

  It copies the Vxc potential from ABINIT to BigDFT or viceversa.

INPUTS

  opt= 1) copy from ABINIT to BigDFT
       2) copy from BigDFT to ABINIT
  vxc(nfft,nspden)= trial potential (ABINIT array)
  wvl_den= density-potential BigDFT object

OUTPUT

SIDE EFFECTS

 vxc is copied to wvl_den, or viceversa, depending on "opt" (see above).

NOTES

 It uses the generic routine wvl_rhov_abi2big.

PARENTS

      wvl_psitohpsi

CHILDREN

SOURCE

387 subroutine wvl_vxc_abi2big(opt,vxc,wvl_den)
388 
389 
390 !This section has been created automatically by the script Abilint (TD).
391 !Do not modify the following lines by hand.
392 #undef ABI_FUNC
393 #define ABI_FUNC 'wvl_vxc_abi2big'
394 !End of the abilint section
395 
396  implicit none
397 
398 !Arguments ------------------------------------
399  integer,intent(in) :: opt
400  real(dp),intent(inout) :: vxc(:,:)
401  type(wvl_denspot_type), intent(inout) :: wvl_den
402 
403 !Local variables-------------------------------
404 #if defined HAVE_BIGDFT
405  integer :: shiftV
406 #endif
407 
408 ! *************************************************************************
409 
410  DBG_ENTER("COLL")
411 
412 #if defined HAVE_BIGDFT
413 ! write(message,'(2a)') ch10, 'wvl_vxc_abi2big: but why are you copying me :..o('
414 ! call wrtout(std_out,message,'COLL')
415 
416  shiftV=wvl_den%denspot%dpbox%ndims(1)*wvl_den%denspot%dpbox%ndims(2) &
417 &      *wvl_den%denspot%dpbox%i3xcsh
418 
419  call wvl_rhov_abi2big(opt,vxc,wvl_den%denspot%v_xc,shift=shiftV)
420 
421 #else
422  BIGDFT_NOTENABLED_ERROR()
423  if (.false.) write(std_out,*) opt,vxc(1,1),wvl_den%symObj
424 #endif
425 
426  DBG_EXIT("COLL")
427 
428 end subroutine wvl_vxc_abi2big