Calculation of a reciprocal quantity form its real counterpart

Problems running VASP: crashes, internal errors, "wrong" results.


Moderators: Global Moderator, Moderator

Post Reply
Message
Author
admin
Administrator
Administrator
Posts: 2921
Joined: Tue Aug 03, 2004 8:18 am
License Nr.: 458

Calculation of a reciprocal quantity form its real counterpart

#1 Post by admin » Tue Aug 05, 2014 3:22 pm

Hello,

I have a question related to the calculation with VASP (in the source code) of a quantity in the reciprocal grid from its homologue in the real grid. In my case, the quantity is defined on a 3-D real space grid (with GRIDC%NGX,GRIDC%NGY,GRIDC%NGZ gridpoints along each of the three space directions). Imagine the quantity is denoted as 'AIN' and defined as a matrix AIN(I,J,K) where each index is defined such as 'I=1,GRIDC%NGX, J =1,GRIDC%NGY, K=1,GRIDC%NGZ. We want a matrix 'AOUT' which is the reciprocal of the quantity AOUT.

A possible way is calling the following subroutine RLTORC(GRID,AIN,AOUT):

********************************************************++

SUBROUTINE RLTORC(GRID,AIN,AOUT)
USE prec
USE mpimy
USE mgrid

TYPE (grid_3d) GRID
INTEGER NALLOC,NZ,N,NODE_ME,IONODE
INTEGER ISTAT,L
REAL(q), ALLOCATABLE :: AV(:)
COMPLEX(q) AOUT(GRID%RC%NP)
REAL(q) AIN(GRID%NGX,GRID%NGY,GRID%NGZ)

NALLOC=GRID%NGX*GRID%NGY
NODE_ME=0
IONODE =0
#ifdef MPI
NODE_ME=GRID%COMM%NODE_ME
IONODE =GRID%COMM%IONODE
#endif

ALLOCATE(AL(NALLOC),STAT=ISTAT)

DO NZ=1,GRID%NGZ
L=0
io_begin
DO NY=1,GRID%NGY
DO NX=1,GRID%NGX
L=L+1
V(L)=AIN(NX,NY,NZ)
END DO
ENDDO
io_end
CALL DIS_GRID_RL_PLANE(GRID,V, AOUT, .TRUE., NZ )
ENDDO

CALL FFT_RC_SCALE(AOUT,AOUT,GRID)
CALL SETUNB_COMPAT(AOUT,GRID)
DEALLOCATE(AV)
RETURN

END SUBROUTINE

**********************************************************************

This is similar than what is done in the subroutine 'INCHG' in the fileio.F file. However in the INCHG subroutine the quantity to convert is just a charge density. In my case I would like to convert any quantity from one grid to the other (a potential for instance). Then, it is possible that instead of using the subroutine 'FFT_RC_SCALE' in the subroutine 'RLTORC' that I propose, I would need just to call FFT3D(AOUT,GRID,-1) and the code would be like this:


********************************************************++

SUBROUTINE RLTORC(GRID,AIN,AOUT)
USE prec
USE mpimy
USE mgrid

TYPE (grid_3d) GRID
INTEGER NALLOC,NZ,N,NODE_ME,IONODE
INTEGER ISTAT,L
REAL(q), ALLOCATABLE :: AV(:)
COMPLEX(q) AOUT(GRID%RC%NP)
REAL(q) AIN(GRID%NGX,GRID%NGY,GRID%NGZ)

NALLOC=GRID%NGX*GRID%NGY
NODE_ME=0
IONODE =0
#ifdef MPI
NODE_ME=GRID%COMM%NODE_ME
IONODE =GRID%COMM%IONODE
#endif

ALLOCATE(AL(NALLOC),STAT=ISTAT)

DO NZ=1,GRID%NGZ
L=0
io_begin
DO NY=1,GRID%NGY
DO NX=1,GRID%NGX
L=L+1
V(L)=AIN(NX,NY,NZ)
END DO
ENDDO
io_end
CALL DIS_GRID_RL_PLANE(GRID,V, AOUT, .TRUE., NZ )
ENDDO

! (optional) CALL RL_ADD(AOUT,1.0_q,AOUT,0.0_q,AOUT,GRID)

CALL FFT3D(AOUT,GRID,-1)
CALL SETUNB_COMPAT(AOUT,GRID)
DEALLOCATE(AV)
RETURN

END SUBROUTINE

**********************************************************************

As the question is not complicated, and I expose the code, I would appreciate very much an answer.

Many thanks in advance
Last edited by admin on Tue Aug 05, 2014 3:22 pm, edited 1 time in total.

Post Reply