Solution (source code)

= Solution
{c}

Taking the <logarithm> of the normalized <probability density function>, then summing over the <independent> observations, gives the <log-likelihood>
$$
\boxed{\ell(p,\beta;x)
=n\{\log\beta-\log(-\log p)+\log(1-p)\}
-\beta\sum_{i=1}^nx_i
-\sum_{i=1}^n\log\!\left(1-(1-p)e^{-\beta x_i}\right).}
$$
Both $p$ and $\beta$ are parameters being inferred, so their contributions to the <normalizing constant> must be retained. The following <R (programming language)> implementation validates the arguments and returns one number:

``
ll_f <- function(p, beta, x) {
  scalar <- function(v) {
    is.numeric(v) && !is.complex(v) && length(v) == 1L &&
      is.finite(v)
  }
  if (!scalar(p) || p <= 0 || p >= 1)
    stop("p must lie strictly between 0 and 1")
  if (!scalar(beta) || beta <= 0)
    stop("beta must be positive")
  if (!is.numeric(x) || is.complex(x) || !is.null(dim(x)) ||
      any(!is.finite(x)) || any(x <= 0))
    stop("x must be a vector of finite positive numbers")
  if (length(x) == 0L) return(0)
  t <- log1p(-p) - beta * x
  length(x) * (log(beta) - log(-log(p))) +
    sum(t - log(-expm1(t)))
}
``

Here $t_i=\log(1-p)-\beta x_i<0$ and `-expm1(t)` computes $1-e^t$ accurately near zero. The empty-data <likelihood function> is one, so its <log-likelihood> is zero. Invalid parameter values produce errors in this evaluator; the next algorithm rejects candidates outside the parameter <support> before calling it.