anchor.diff:    
           C_center(1:3) = nucl_coord(k,1:3)
     
  -        v_k_dump = v_k(k,1:klocmax)
  -        n_k_dump = n_k(k,1:klocmax)
  -        dz_k_dump =  dz_k(k,1:klocmax)
  +        v_k_dump = pseudo_v_k(k,1:pseudo_klocmax)
  +        n_k_dump = pseudo_n_k(k,1:pseudo_klocmax)
  +        dz_k_dump =  pseudo_dz_k(k,1:pseudo_klocmax)
   
  -        c = c + Vloc(klocmax, v_k_dump,n_k_dump, dz_k_dump, &
  +        c = c + Vloc(pseudo_klocmax, v_k_dump,n_k_dump, dz_k_dump, &
                        A_center,power_A,alpha,B_center,power_B,beta,C_center)
     
  +
  +        n_kl_dump = pseudo_n_kl(k,1:pseudo_kmax,0:pseudo_lmax)
  +        v_kl_dump = pseudo_v_kl(k,1:pseudo_kmax,0:pseudo_lmax)
  +        dz_kl_dump = pseudo_dz_kl(k,1:pseudo_kmax,0:pseudo_lmax)
  +
  +        c = c + Vpseudo(pseudo_lmax,pseudo_kmax,v_kl_dump,n_kl_dump,dz_kl_dump,A_center,power_A,alpha,B_center,power_B,beta,C_center)
  +  
         enddo
  -      ao_nucl_elec_integral_pseudo_local(i,j) = ao_nucl_elec_integral_pseudo_local(i,j) + &
  +      ao_nucl_elec_integral_pseudo(i,j) = ao_nucl_elec_integral_pseudo(i,j) + &
                                      ao_coef_normalized_ordered_transp(l,j)*ao_coef_normalized_ordered_transp(m,i)*c
        enddo
        enddo
  
positive.code: BEGIN_PROVIDER [ double precision, ao_nucl_elec_integral_pseudo, (ao_num_align,ao_num) ]
   implicit none
  ! Nucleus-pseudopotential interaction
   END_DOC
   write(output_monoints,*) 'Providing the pseudopotential integrals '
   ao_nucl_elec_integral_pseudo = ao_nucl_elec_integral_pseudo_local + ao_nucl_elec_integral_pseudo_non_local
  END_PROVIDER
  
  
  BEGIN_PROVIDER [ double precision, ao_nucl_elec_integral_pseudo_local, (ao_num_align,ao_num)]
   implicit none
   BEGIN_DOC
  ! Local component of the pseudopotential
   END_DOC
   double precision  :: alpha, beta, gama, delta
   integer           :: num_A,num_B
   double precision  :: A_center(3),B_center(3),C_center(3)
  
negative.code:  double precision  :: cpu_1, cpu_2, wall_1, wall_2, wall_0
   integer           :: thread_num
  
    ao_nucl_elec_integral_pseudo = 0.d0
  
    !                 
    ! |   _   _  _. | 
    ! |_ (_) (_ (_| | 
    !                 
    !! Parameters of the local part of pseudo:
  
  !  integer klocmax
  !  integer, allocatable ::  n_k(:,:)
  !  double precision, allocatable ::  v_k(:,:), dz_k(:,:)
  !
  !  call ezfio_get_pseudo_klocmax(klocmax)
  !
  !  allocate(n_k(nucl_num,klocmax),v_k(nucl_num,klocmax), dz_k(nucl_num,klocmax))
  !
  !  call ezfio_get_pseudo_v_k(v_k)
  !  call ezfio_get_pseudo_n_k(n_k)
  !  call ezfio_get_pseudo_dz_k(dz_k)
  !
    !! Dump array 
    integer, allocatable ::  n_k_dump(:)
    double precision, allocatable ::  v_k_dump(:), dz_k_dump(:) 
   
    allocate(n_k_dump(1:pseudo_klocmax), v_k_dump(1:pseudo_klocmax), dz_k_dump(1:pseudo_klocmax)) 
  
  
    !                               
    ! |\ |  _  ._    |  _   _  _. | 
    ! | \| (_) | |   | (_) (_ (_| | 
    !                              
    !! Parameters of non local part of pseudo:
  
  !  integer :: kmax,lmax
  !  integer, allocatable ::  n_kl(:,:,:)
  !  double precision, allocatable ::  v_kl(:,:,:), dz_kl(:,:,:) 
  !
  !  call ezfio_get_pseudo_lmaxpo(lmax)
  !  call ezfio_get_pseudo_kmax(kmax)
  !  !lmax plus one -> lmax
  !  lmax = lmax - 1 
  ! 
  !  allocate(n_kl(nucl_num,kmax,0:lmax), v_kl(nucl_num,kmax,0:lmax), dz_kl(nucl_num,kmax,0:lmax)) 
  !
  !  call ezfio_get_pseudo_n_kl(n_kl)
  !  call ezfio_get_pseudo_v_kl(v_kl)
  !  call ezfio_get_pseudo_dz_kl(dz_kl)
  !
  !
    !! Dump array 
    integer, allocatable ::  n_kl_dump(:,:)
    double precision, allocatable ::  v_kl_dump(:,:), dz_kl_dump(:,:) 
   
    allocate(n_kl_dump(pseudo_kmax,0:pseudo_lmax), v_kl_dump(pseudo_kmax,0:pseudo_lmax), dz_kl_dump(pseudo_kmax,0:pseudo_lmax)) 
  
    !  _                
    ! /   _. |  _     | 
    ! \_ (_| | (_ |_| | 
    !                   
  
    write(output_monoints,*) 'Providing the nuclear electron pseudo integrals '
  
    call wall_time(wall_1)
    call cpu_time(cpu_1)
  