anchor.diff:    end if
   
     call initialize_dirichlet(counts, nc)
  -  call compute_multiplicities(counts, nc)
  +  call compute_multiplicities(counts)
   
     call compute_integration_range()
   
  
positive.code:     log_fpxa = log_gamma(n_data + one) + log_gamma(alpha * alphabet_size) & 
           - alphabet_size * log_gamma(alpha) - log_gamma(n_data + alpha * alphabet_size)
      do i_ = 1, n_multi 
         a(i_) = multi(i_) * (log_gamma(multi_z(i_) - one + alpha) - log_gamma(multi_z(i_) * one))
      end do
      log_fpxa = log_fpxa + sum(a)
      
  
negative.code: 
      hdir = 0.0_real64
      do i = 1,n_multi 
         a(i) = - multi(i) * (multi_z(i) + alpha) * digamma(multi_z(i) + one + alpha) 
      end do
      hdir = sum(a)
      hdir = hdir / (n_data + alpha * alphabet_size)
  