1!
2! Copyright (C) 1996-2016	The SIESTA group
3!  This file is distributed under the terms of the
4!  GNU General Public License: see COPYING in the top directory
5!  or http://www.gnu.org/copyleft/gpl.txt.
6! See Docs/Contributors.txt for a list of contributors.
7!
8      subroutine savepsi(psiprev, psi, nuo, nuotot, nocc)
9C *********************************************************************
10C Copies psi into psiprev allowing for the fact that the parallel form
11C requires a re-distribution of the data.
12C Written by J.D. Gale, March 2000
13C **************************** INPUT **********************************
14C real*8  psi(2,nuotot,nuo)     : Wavefunctions in current k point
15C integer nuo                   : Number of (local) orbitals in the unit cell
16C integer nuotot                : Number of orbitals in the unit cell
17C integer nocc                  : number of occupied states
18C *************************** OUTPUT **********************************
19C real*8  psiprev(2,nuo,nuotot) : Wavefunctions in previous k point
20C *************************** UNITS ***********************************
21C Lengths in atomic units (Bohr).
22C k vectors in reciprocal atomic units.
23C Energies in Rydbergs.
24C *********************************************************************
25      use precision
26      use parallel,     only : Node, Nodes
27      use parallelsubs, only : GetNodeOrbs, LocalToGlobalOrb
28      use alloc,        only : re_alloc, de_alloc
29#ifdef MPI
30      use mpi_siesta
31#endif
32      implicit none
33
34      integer nuo, nuotot, nocc
35
36      real(dp)
37     .  psiprev(2,nuo,nuotot),
38     .  psi(2,nuotot,nuo)
39
40C**** Internal variables ***********************************************
41
42      integer
43     .  iuo, juo
44#ifdef MPI
45      integer
46     .  MPIerror, iuog, juog, n, noccloc
47      real(dp), dimension(:,:,:), pointer ::  psitmp
48#endif
49
50#ifdef MPI
51! AG
52! Allocate as in detover, using the number of orbitals on the first node,
53! as some of the nodes might have zero orbitals!
54
55      call GetNodeOrbs(nocc,0,Nodes,noccloc)
56      nullify( psitmp )
57      call re_alloc( psitmp, 1, 2, 1, nuotot, 1, noccloc,
58     &                 name='psitmp', routine='savepsi' )
59
60      do n = 0,Nodes-1
61
62C Broadcast copy of psi on node n to all other nodes
63        call GetNodeOrbs(nocc,n,Nodes,noccloc)
64
65        if (Node .eq. n) then
66           psitmp(1:2,1:nuotot,1:noccloc) = psi(1:2,1:nuotot,1:noccloc)
67        endif
68        call MPI_Bcast(psitmp(1,1,1),2*nuotot*noccloc,
69     .    MPI_double_precision,n,MPI_Comm_World,MPIerror)
70
71C Save local part of psiprev
72        do iuo = 1,noccloc
73          call LocalToGlobalOrb(iuo,n,Nodes,iuog)
74          do juo = 1,nuo
75            call LocalToGlobalOrb(juo,Node,Nodes,juog)
76            psiprev(1,juo,iuog) = psitmp(1,juog,iuo)
77            psiprev(2,juo,iuog) = psitmp(2,juog,iuo)
78          enddo
79        enddo
80
81
82      enddo
83      call de_alloc( psitmp,  name='psitmp' , routine='savepsi')
84
85#else
86C Straight serial copy
87
88      psiprev(1:2,1:nuotot,1:nocc) = psi(1:2,1:nuotot,1:nocc)
89
90#endif
91
92      end
93