← Back to branch summary

~albertog/siesta/4.1-wfs

« back to all changes in this revision

Viewing changes to Util/Grid/m_gridfunc.F90

  • Committer: Alberto Garcia
  • Date: 2019-04-30 17:54:47 UTC
  • Revision ID: albertog@icmab.es-20190430175447-fc3hii6czqp1b256
Update manual for ol-stm. Print V(z) in vacpot

Re-designed vacpot using the functionality in m_gridfunc (also
available through Util/Grid/v_info).

Updated the manual.

Still need a proper set of examples.

Show diffs side-by-side

added added

removed removed

Lines of Context:
Util/Grid/m_gridfunc.F90
6
6
! See Docs/Contributors.txt for a list of contributors.
7
7
! ---
8
8
!
 
9
!> Utilities to dea with real-space grid functions
9
10
module m_gridfunc
10
11
 
11
12
  implicit          none
15
16
  integer, parameter          :: grid_p = sp  ! NOTE this
16
17
 
17
18
  type, public :: gridfunc_t
18
 
     real(dp)              ::  cell(3,3)   ! by columns, in bohr
19
 
     real(dp)              ::  origin(3) = [0.0_dp, 0.0_dp, 0.0_dp]  ! in bohr
 
19
     !> Lattice vectors (by columns, in bohr)
 
20
     real(dp)              ::  cell(3,3)   
 
21
     !> Possible shift for the origin of the box, useful in some
 
22
     !> situations (sub-boxes, for example) (in bohr)
 
23
     real(dp)              ::  origin(3) = [0.0_dp, 0.0_dp, 0.0_dp]
 
24
     !> The 'grid' format can in principle hold functions that
 
25
     !> are not periodic (for example, values in a sub-box)
20
26
     logical               ::  is_periodic(3) = [ .true., .true., .true.]
 
27
     !> Number of points along each cell direction
21
28
     integer               ::  n(3) = [ 0,0,0 ]
 
29
     !> Extra dimension, typically the spin index
22
30
     integer               ::  nspin
 
31
     !> Array holding the values of the grid function
23
32
     real(grid_p), allocatable ::  val(:,:,:,:)
24
33
  end type gridfunc_t
25
34
 
28
37
      public :: read_gridfunc_netcdf, write_gridfunc_netcdf
29
38
#endif
30
39
      public :: get_planar_average, grid_p, monoclinic_z
 
40
      public :: monoclinic
 
41
    
31
42
      private
32
43
 
33
44
 CONTAINS
71
82
 
72
83
   end subroutine write_gridfunc
73
84
 
 
85
   !> Computes a simple-minded 'planar average' of a grid
 
86
   !> function. This is a mathematical average. For a more
 
87
   !> physical average, see the 'macroave' utility
74
88
   subroutine get_planar_average(gf, dim, val)
 
89
     !> The object representing the grid function
75
90
     type(gridfunc_t), intent(in) :: gf
 
91
     !> The dimension that indexes the average values.
 
92
     !> These are obtained by summing over the two other dimensions.
76
93
     integer, intent(in) :: dim
 
94
     !> Array holding the averages (second dimension is the 'spin'))
77
95
     real(grid_p), allocatable :: val(:,:)
78
96
 
79
97
     integer :: isp, ix, iy, iz, n(3), spin_dim
81
99
 
82
100
     n = gf%n
83
101
     spin_dim = size(gf%val,dim=4)
84
 
     
85
 
     if (dim == 3) then
86
 
        allocate(val(n(3),spin_dim))
87
 
        do isp=1,gf%nspin
 
102
 
 
103
     allocate(val(n(dim),spin_dim))
 
104
 
 
105
     do isp=1,gf%nspin
 
106
        select case (dim)
 
107
        case (3)
88
108
           do iz=1,n(3)
89
109
              sum = 0.0_grid_p
90
110
              do iy=1,n(2)
94
114
              enddo
95
115
              val(iz,isp) = sum/(n(1)*n(2))
96
116
           enddo
97
 
        enddo
98
 
     else
99
 
        call die("non-z ave not implemented yet")
100
 
     endif
 
117
        case (2)
 
118
           do iy=1,n(2)
 
119
              sum = 0.0_grid_p
 
120
              do iz=1,n(3)
 
121
                 do ix=1,n(1)
 
122
                    sum = sum + gf%val(ix,iy,iz,isp)
 
123
                 enddo
 
124
              enddo
 
125
              val(iy,isp) = sum/(n(1)*n(3))
 
126
           enddo
 
127
        case (1)
 
128
           do ix=1,n(1)
 
129
              sum = 0.0_grid_p
 
130
              do iz=1,n(3)
 
131
                 do iy=1,n(2)
 
132
                    sum = sum + gf%val(ix,iy,iz,isp)
 
133
                 enddo
 
134
              enddo
 
135
              val(ix,isp) = sum/(n(2)*n(3))
 
136
           enddo
 
137
        end select
 
138
        
 
139
     enddo
101
140
 
102
141
   end subroutine get_planar_average
103
142
 
333
372
   stop
334
373
 end subroutine die
335
374
 
 
375
 !> Physically, it should not matter whether the 'c' axis points along
 
376
 !> z, since rotations are irrelevant. However, some programs are
 
377
 !> hard-wired for the z-direction. Hence the name. See the
 
378
 !> 'monoclinic' function for the rotationally invariant
 
379
 !> version.
 
380
 
336
381
 function monoclinic_z(cell)
337
382
   real(dp), intent(in) :: cell(3,3)
338
383
   logical monoclinic_z
348
393
 
349
394
   end function monoclinic_z
350
395
 
 
396
   !> Checks that the lattice vector of index 'dim' is orthogonal to
 
397
   !> the other two.
 
398
 
 
399
   function monoclinic(cell,dim)
 
400
     real(dp), intent(in) :: cell(3,3)
 
401
     integer, intent(in)  :: dim
 
402
     logical monoclinic
 
403
 
 
404
     real(dp), parameter :: tol = 1.0e-8_dp
 
405
 
 
406
   ! Check that the 'dim' vector is orthogonal to the other two
 
407
     monoclinic = .true.
 
408
     if (abs(dot_product(cell(:,dim),cell(:,next(dim,1)))) > tol) then
 
409
        monoclinic = .false.
 
410
     endif
 
411
     if (abs(dot_product(cell(:,dim),cell(:,next(dim,2)))) > tol) then
 
412
        monoclinic = .false.
 
413
     endif
 
414
 
 
415
   CONTAINS
 
416
     !> generates cyclic neighbors of dimension 'i' with shift 'm'
 
417
     !> For example:
 
418
     !>   x --> y:   next(1,1) = 2
 
419
     !>   y --> x:   next(2,2) = 1
 
420
     function next(i,m) result(ind)
 
421
       integer, intent(in) :: i, m
 
422
       integer :: ind
 
423
       
 
424
       ind = i + m
 
425
       if (ind > 3) then
 
426
          ind = ind - 3
 
427
       endif
 
428
     end function next
 
429
   end function monoclinic
 
430
 
351
431
 end module m_gridfunc
352
432