c_adaptive_mix Subroutine

subroutine c_adaptive_mix(z, Fz, alpha, iter)

Arguments

Type IntentOptional Attributes Name
complex(kind=8), intent(inout), dimension(:) :: z
complex(kind=8), intent(in), dimension(size(z)) :: Fz
real(kind=8), intent(in) :: alpha
integer, intent(in) :: iter

Calls

proc~~c_adaptive_mix~~CallsGraph proc~c_adaptive_mix c_adaptive_mix dimag dimag proc~c_adaptive_mix->dimag

Source Code

subroutine c_adaptive_mix(z,Fz,alpha,iter)
   complex(8),intent(inout),dimension(:)    :: z       !z_in
   complex(8),intent(in),dimension(size(z)) :: Fz      !z_out-z_in
   real(8),intent(in)                       :: alpha   !mixing α
   integer,intent(in)                       :: iter
   integer                                  :: N
   real(8),parameter                        :: alpha_max = 1d0 
   complex(8),allocatable,dimension(:),save :: Fz_prev !prev iter
   real(8),allocatable,dimension(:),save    :: beta    !mixing β
   integer                                  :: j
   logical                                  :: real_oscillates, imag_oscillates
   !
   N = size(z)
   !
   if(iter==1)then !reset persistent variables
      if(allocated(Fz_prev))deallocate(Fz_prev)
      allocate(Fz_prev(N))
      Fz_prev = zero
      !
      if(allocated(beta))deallocate(beta)
      allocate(beta(N))
      beta = alpha
   endif
   !
   Z = Z + beta*Fz ! = (1-β)z_in + βz_out (β: effective speed)            
   !
   if (iter > 1) then !update effective elemental mixing speed
      do j=1,N
         real_oscillates = dreal(Fz_prev(j))*dreal(Fz(j)) < 0
         imag_oscillates = dimag(Fz_prev(j))*dimag(Fz(j)) < 0
         if(real_oscillates.OR.imag_oscillates)then 
            beta(j) = alpha !reset to input mixing
         else !increase mixing speed by adding alpha
            beta(j) = beta(j) + alpha !but with a cutoff
            if (beta(j) > alpha_max) beta(j) = alpha_max
         end if
      end do
   end if
   !
   Fz_prev = Fz !store Fz_prev for future calls
   !
end subroutine c_adaptive_mix