LCOV - code coverage report
Current view: top level - wannier - wann_dipole_ionic.f (source / functions) Hit Total Coverage
Test: combined.info Lines: 0 53 0.0 %
Date: 2019-09-08 04:53:50 Functions: 0 1 0.0 %

          Line data    Source code
       1             : !--------------------------------------------------------------------------------
       2             : ! Copyright (c) 2016 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
       3             : ! This file is part of FLEUR and available as free software under the conditions
       4             : ! of the MIT license as expressed in the LICENSE file in more detail.
       5             : !--------------------------------------------------------------------------------
       6             : 
       7             :       module m_wann_dipole_ionic
       8             :       contains
       9           0 :       subroutine wann_dipole_ionic(
      10           0 :      >               natd,pos,omtil,
      11           0 :      >               amat,taual,ntype,
      12           0 :      >               neq,zatom,l_absolute,
      13             :      >               l_invsubtract,
      14             :      <               ionic_moment)
      15             :       use m_wann_ioncharge_gen
      16             :       implicit none
      17             :       real,intent(in)              :: omtil
      18             :       integer,intent(in)           :: natd
      19             :       real,intent(in)              :: pos(3,natd)
      20             :       real,intent(in)              :: amat(3,3)
      21             :       real,intent(in)              :: taual(3,natd)
      22             :       integer,intent(in)           :: ntype
      23             :       integer,intent(in)           :: neq(ntype)
      24             :       real,intent(in)              :: zatom(ntype)
      25             :       logical,intent(in)           :: l_absolute
      26             :       logical,intent(in)           :: l_invsubtract
      27             :       real,intent(out)             :: ionic_moment(3)
      28             : 
      29             :       integer                      :: i
      30             :       integer                      :: num_atoms
      31             :       integer                      :: j,k,ind,jj,jspin
      32           0 :       real,allocatable             :: ioncharge(:)
      33             :       real                         :: polarization
      34             :       real                         :: ionic_polarization(3)
      35             :       real,parameter               :: elemchargmu=1.60217646e-13
      36             :       real,parameter               :: bohrtocm=0.529177e-8
      37             :       character(len=6)             :: filename
      38             :       logical                      :: l_file
      39           0 :       real                         :: pos_inv(3,natd)
      40           0 :       real                         :: taual_inv(3,natd)
      41             :       real                         :: coordinate
      42             :       real                         :: ionchargesum
      43             : 
      44           0 :       write(666,*)"*****************************"
      45           0 :       write(666,*)" Ionic terms                 "
      46           0 :       write(666,*)"*****************************"
      47           0 :       write(*,*)  "*****************************"
      48           0 :       write(*,*)  " Ionic terms                 "
      49           0 :       write(*,*)  "*****************************"
      50             : 
      51             : c-----count atoms
      52           0 :       num_atoms=0
      53           0 :       do j=1,ntype
      54           0 :          do k=1,neq(j)
      55           0 :             num_atoms=num_atoms+1
      56             :          enddo
      57             :       enddo
      58             : 
      59             : c-----get charges of ions
      60           0 :       allocate( ioncharge(num_atoms) )
      61             :       call wann_ioncharge_gen(
      62             :      >         num_atoms,ntype,natd,
      63             :      >         neq,zatom,pos,
      64           0 :      <         ioncharge)
      65           0 :       ionchargesum=sum(ioncharge(1:num_atoms))
      66           0 :       write(666,*)"ionchargesum=",ionchargesum
      67             : 
      68             : c-----Try to find the atomic positions of the inversion symmetric counterpart.
      69             : c-----The ferroelectric system is assumed to be connected to an inversion 
      70             : c-----symmetric system by an adiabatic path. The polarization of the inversion
      71             : c-----symmetric system is zero and consequently ambiguities in the evaluation of
      72             : c-----the polarization can be removed by subtracting the calculated polarizations
      73             : c-----of the two systems. We do not care here, whether the calculated inversion
      74             : c-----symmetric positions are physically reasonable: For example two atoms might
      75             : c-----end up at the same location. But this does not matter.
      76           0 :       pos_inv=0.0
      77           0 :       if(l_invsubtract)then
      78           0 :         do j=1,num_atoms
      79           0 :           do k=1,3
      80           0 :                coordinate=taual(k,j)
      81           0 :                if(abs(coordinate).lt.0.125)then
      82           0 :                   taual_inv(k,j)=0.0
      83           0 :                elseif(abs(coordinate).lt.0.375)then
      84           0 :                   taual_inv(k,j)=0.25*coordinate/abs(coordinate)
      85           0 :                elseif(abs(coordinate).lt.0.625)then
      86           0 :                   taual_inv(k,j)=0.5*coordinate/abs(coordinate)
      87           0 :                elseif(abs(coordinate).lt.0.875)then
      88           0 :                   taual_inv(k,j)=0.75*coordinate/abs(coordinate)
      89             :                else
      90           0 :                   taual_inv(k,j)=coordinate/abs(coordinate)
      91             :                endif   
      92             :           enddo
      93             :         enddo
      94           0 :         do j=1,num_atoms
      95           0 :            pos_inv(:,j)=matmul(amat,taual_inv(:,j))
      96             :         enddo
      97           0 :         do j=1,num_atoms
      98           0 :           print*,"old position:"
      99           0 :           print*,pos(:,j)
     100           0 :           print*,"new position:"
     101           0 :           print*,pos_inv(:,j)
     102           0 :           print*,"************************"
     103             :         enddo
     104             :       endif
     105             : 
     106             : c-----compute the ionic moment
     107           0 :       ionic_moment=0.0
     108           0 :       do j=1,num_atoms
     109             :          ionic_moment(:)=
     110           0 :      &        ionic_moment(:)+ioncharge(j)*(pos(:,j)-pos_inv(:,j))
     111             :       enddo
     112             : 
     113             :       ionic_polarization =
     114           0 :      =   ionic_moment/omtil*elemchargmu/((bohrtocm)**2)
     115             : 
     116           0 :       write(*,  fmt=555)ionic_moment(:)
     117           0 :       write(666,fmt=555)ionic_moment(:)
     118           0 :       if(.not.l_absolute)then
     119           0 :         write(*,  fmt=777)ionic_polarization(:)
     120           0 :         write(666,fmt=777)ionic_polarization(:)
     121             :       endif
     122             : 
     123             :  555  format("ionic moment = (",3f12.6,")a.u.")
     124             :  777  format("ionic polarization = (",3f12.6,")uC/cm**2")
     125             : 
     126           0 :       end subroutine wann_dipole_ionic
     127             :       end module m_wann_dipole_ionic
     128             : 
     129             : 
     130             : 
     131             : 
     132             : 
     133             : 
     134             : 
     135             : 
     136             : 
     137             : 

Generated by: LCOV version 1.13