LCOV - code coverage report
Current view: top level - mix - a_pulay.F90 (source / functions) Hit Total Coverage
Test: combined.info Lines: 0 38 0.0 %
Date: 2019-09-08 04:53:50 Functions: 0 2 0.0 %

          Line data    Source code
       1             : MODULE m_a_pulay
       2             :   USE m_juDFT
       3             :   PRIVATE
       4             :   REAL :: distance(3)
       5             :   PUBLIC a_pulay
       6             : CONTAINS
       7           0 :   SUBROUTINE a_pulay(alpha,fm,sm)
       8             :     USE m_pulay
       9             :     USE m_types_mat
      10             :     USE m_types_mixvector
      11             :     USE m_mixing_history
      12             :     IMPLICIT NONE
      13             :     REAL,INTENT(IN)                 :: alpha
      14             :     TYPE(t_mixvector),INTENT(IN)    :: fm(:)
      15             :     TYPE(t_mixvector),INTENT(INOUT) :: sm(:)
      16             : 
      17             : 
      18             :     REAL    :: fac(2)=(/1.0,0.5/)
      19             :     INTEGER,parameter :: parallel_images=2
      20             :     INTEGER,PARAMETER :: maxiter=4,mode=2
      21             :     INTEGER :: hlen,image,i,ii,local_length
      22             :     INTEGER :: local_hist(maxiter*parallel_images)
      23             : 
      24           0 :     hlen=SIZE(sm) !Number of results in history
      25           0 :     image=MOD(hlen-1,parallel_images)+1
      26           0 :     PRINT *,"HI:",hlen,image
      27           0 :     IF (hlen<parallel_images) THEN
      28             :        !create inital parallel images by simple mixing of first result
      29           0 :        sm(hlen)=sm(1)+(alpha*fac(image))*fm(1)
      30           0 :        RETURN
      31             :     ENDIF
      32             : 
      33             :     SELECT CASE(mode)
      34             :     CASE(1)
      35             :        local_length=0
      36             :        DO i=image,hlen,parallel_images
      37             :           local_length=local_length+1
      38             :           local_hist(local_length)=i
      39             :        ENDDO
      40             :        
      41             :        IF (local_length>maxiter) THEN !recreate to enable cross-image mixing
      42             :           local_length=hlen+1-parallel_images
      43             :           local_hist(1)=image
      44             :           local_hist(2:local_length)=(/(i,i=parallel_images+1,hlen)/)
      45             :        END IF
      46             :        sm(hlen)=simple_pulay(alpha*fac(image),fm(local_hist(:local_length)),sm(local_hist(:local_length)))
      47           0 :        IF (hlen==parallel_images*(maxiter+1)) CALL mixing_history_limit(parallel_images)
      48             :     case(2)
      49           0 :        local_length=hlen-image
      50           0 :        local_hist(1:local_length)=(/(i,i=1,local_length)/)
      51           0 :        local_length=local_length+1
      52           0 :        local_hist(local_length)=hlen
      53           0 :        sm(hlen)=simple_pulay(alpha*fac(image),fm(local_hist(:local_length)),sm(local_hist(:local_length)))
      54           0 :        IF (hlen==parallel_images*(maxiter)) CALL mixing_history_limit(parallel_images)
      55             :     END SELECT
      56             : 
      57           0 :     PRINT *,"P:",local_hist(:local_length)
      58             :     
      59             :     
      60             :   END SUBROUTINE a_pulay
      61             : 
      62           0 :   FUNCTION simple_pulay(alpha,fm,sm)RESULT(sm_out)
      63             :     USE m_types_mixvector
      64             :     IMPLICIT NONE
      65             :     REAL,INTENT(IN)                 :: alpha
      66             :     TYPE(t_mixvector),INTENT(IN)    :: fm(:)
      67             :     TYPE(t_mixvector),INTENT(IN)    :: sm(:)
      68             :     TYPE(t_mixvector)               :: sm_out
      69             :     
      70             :     
      71           0 :     TYPE(t_mixvector)   :: mdf
      72           0 :     REAL,ALLOCATABLE    :: a(:,:),b(:),work(:)
      73           0 :     INTEGER,ALLOCATABLE :: ipiv(:)
      74             :     INTEGER             :: n,nn,info,lwork,hlen
      75             :     
      76             :     
      77             : 
      78           0 :     hlen=SIZE(sm)
      79           0 :     IF (hlen==1) THEN
      80           0 :        sm_out=sm(1)+alpha*fm(1)
      81           0 :        RETURN
      82             :     ENDIF
      83             : 
      84           0 :     ALLOCATE(a(hlen+1,hlen+1),b(hlen+1),ipiv(hlen+1),work((hlen+1)**2))
      85           0 :     DO n=1,hlen
      86           0 :        mdf=fm(n)%apply_metric()
      87           0 :        DO nn=1,n
      88           0 :           a(n,nn)=mdf.dot.fm(nn)
      89           0 :           a(nn,n)=a(n,nn)
      90             :        ENDDO
      91             :     ENDDO
      92             : 
      93           0 :     a(:,hlen+1)=1.0
      94           0 :     a(hlen+1,:)=1.0
      95           0 :     a(hlen+1,hlen+1)=0.0
      96           0 :     b=0.0;b(hlen+1)=1.
      97           0 :     CALL  DSYSV( 'Upper', hlen+1, 1, a, SIZE(a,1), ipiv, b, SIZE(b,1), work, SIZE(work), INFO )
      98           0 :     sm_out=0.0*mdf
      99           0 :     DO n=1,hlen
     100           0 :        sm_out=sm_out+b(n)*(sm(n)+alpha*fm(n))
     101             :     ENDDO
     102             :   
     103           0 :   END FUNCTION SIMPLE_PULAY
     104             : 
     105             : END MODULE m_a_pulay

Generated by: LCOV version 1.13