Running examples in ‘MOSAlloc-Ex.R’ failed
The error most likely occurred in:

> ### Name: mosalloc
> ### Title: Multiobjective sample allocation for constraint multivariate and
> ###   multidomain optimal allocation in survey sampling
> ### Aliases: mosalloc
> 
> ### ** Examples
> 
> # Artificial population of 50 568 business establishments and 5 business
> # sectors (data from Valliant, R., Dever, J. A., & Kreuter, F. (2013).
> # Practical tools for designing and weighting survey samples. Springer.
> # https://doi.org/10.1007/978-1-4614-6449-5, Example 5.2 pages 133-9)
> 
> # See also <https://umd.app.box.com/s/9yvvibu4nz4q6rlw98ac/file/297813512360>
> # file: Code 5.3 constrOptim.example.R
> 
> Nh <- c(6221, 11738, 4333, 22809, 5467) # stratum sizes
> ch <- c(120, 80, 80, 90, 150) # stratum-specific cost of surveying
> 
> # Revenues
> mh.rev <- c(85, 11, 23, 17, 126) # mean revenue
> Sh.rev <- c(170.0, 8.8, 23.0, 25.5, 315.0) # standard deviation revenue
> 
> # Employees
> mh.emp <- c(511, 21, 70, 32, 157) # mean number of employees
> Sh.emp <- c(255.50, 5.25, 35.00, 32.00, 471.00) # std. dev. employees
> 
> # Proportion of estabs claiming research credit
> ph.rsch <- c(0.8, 0.2, 0.5, 0.3, 0.9)
> 
> # Proportion of estabs with offshore affiliates
> ph.offsh <- c(0.06, 0.03, 0.03, 0.21, 0.77)
> 
> budget <- 300000 # overall available budget
> n.min  <- 100 # minimum stratum-specific sample size
> 
> # Examples
> #----------------------------------------------------------------------------
> # Example 1: Minimization of the variation of estimates for revenue subject
> # to cost restrictions and precision restrictions to the coefficient of
> # variation of estimates for the proportion of businesses with offshore
> # affiliates.
> 
> l <- rep(n.min, 5) # minimum sample size per stratum
> u <- Nh            # maximum sample size per stratum
> C <- rbind(ch,
+            ch * c(-1, -1, -1, 0, 0))
> c <- c(budget,                        # Maximum overall survey budget
+        - 0.5 * budget)                # Minimum overall budget for strata 1-3
> 
> # We require at maximum 5 % relative standard error for estimates of
> # proportion of businesses with offshore affiliates
> A <- matrix(ph.offsh * (1 - ph.offsh) * Nh**3/(Nh - 1)/sum(Nh * ph.offsh)**2,
+  nrow = 1)
> a <- sum(ph.offsh * (1 - ph.offsh) * Nh**2/(Nh - 1)
+ )/sum(Nh * ph.offsh)**2 + 0.05**2
> 
> D <- matrix(Sh.rev**2 * Nh**2, nrow = 1) # objective variance components
> d <- sum(Sh.rev**2 * Nh) # finite population correction
> 
> opts = list(sense = "max_precision",
+             f = NULL, df = NULL, Hf = NULL,
+             init_w = 1,
+             mc_cores = 1L, pm_tol = 1e-05,
+             max_iters = 100L, print_pm = FALSE)
> 
> sol <- mosalloc(D = D, d = d, A = A, a = a, C = C, c = c, l = l, u = u,
+                 opts = opts)
> 
> # Check solution statement of the internal solver to verify feasibility
> sol$Ecosolver$Ecoinfostring # [1] "Optimal solution found"
[1] "Optimal solution found"
> 
> # Check constraints
> c(C[1, ] %*% sol$n) # [1] 3e+05
[1] 3e+05
> c(C[2, ] %*% sol$n) # [1] -150000
[1] -150000
> c(sqrt(A %*% (1 / sol$n) - A %*% (1 / Nh))) # 5 % rel. std. err.
[1] 0.05
> 
> #----------------------------------------------------------------------------
> # Example 2: Minimization of the maximum relative variation of estimates for
> # the total revenue, the number of employee, the number of businesses claimed
> # research credit, and the number of businesses with offshore affiliates
> # subject to cost restrictions
> 
> l <- rep(n.min, 5) # minimum sample size ber stratum
> u <- Nh            # maximum sample size per stratum
> C <- rbind(ch, ch * c(-1, -1, -1, 0, 0))
> c <- c(budget, - 0.5 * budget)
> A <- NULL # no precision constraint
> a <- NULL # no precision constraint
> 
> # Precision components (Variance / Totals^2) for multidimensional objective
> D <- rbind(Sh.rev**2 * Nh**2/sum(Nh * mh.rev)**2,
+            Sh.emp**2 * Nh**2/sum(Nh * mh.emp)**2,
+            ph.rsch * (1 - ph.rsch) * Nh**3/(Nh - 1)/sum(Nh * ph.rsch)**2,
+            ph.offsh * (1 - ph.offsh) * Nh**3/(Nh - 1)/sum(Nh * ph.offsh)**2)
> 
> d <- as.vector(D %*% (1 / Nh)) # finite population correction
> 
> opts = list(sense = "max_precision",
+             f = NULL, df = NULL, Hf = NULL,
+             init_w = 1,
+             mc_cores = 1L, pm_tol = 1e-05,
+             max_iters = 100L, print_pm = FALSE)
> 
> sol <- mosalloc(D = D, d = d, C = C, c = c, l = l, u = u, opts = opts)
> 
> # Obtain optimal objective value
> sol$J # [1] 0.0017058896 0.0004396972 0.0006428475 0.0017058896
[1] 0.0017058896 0.0004396972 0.0006428475 0.0017058896
> 
> # Obtain corresponding normal vector
> sol$Normal # [1] 6.983113e-01 1.337310e-11 1.596167e-11 3.016887e-01
[1] 6.983113e-01 1.337467e-11 1.596354e-11 3.016887e-01
> 
> # => Revenue and offshore affiliates are dominating the solution with a
> #    ratio of approximately 2:1 (sol$Normal[1] / sol$Normal[4])
> 
> #----------------------------------------------------------------------------
> # Example 3: Example 2 with preference weighting
> 
> w <- c(1, 3.85, 3.8, 1.3) # preference weighting
> l <- rep(n.min, 5) # minimum sample size ber stratum
> u <- Nh            # maximum sample size per stratum
> C <- rbind(ch, ch * c(-1, -1, -1, 0, 0))
> c <- c(budget, - 0.5 * budget)
> A <- NULL # no precision constraint
> a <- NULL # no precision constraint
> 
> D <- rbind(Sh.rev**2 * Nh**2/sum(Nh * mh.rev)**2,
+            Sh.emp**2 * Nh**2/sum(Nh * mh.emp)**2,
+            ph.rsch * (1 - ph.rsch) * Nh**3/(Nh - 1)/sum(Nh * ph.rsch)**2,
+            ph.offsh * (1 - ph.offsh) * Nh**3/(Nh - 1)/sum(Nh * ph.offsh)**2)
> 
> d <- as.vector(D %*% (1 / Nh))
> 
> opts = list(sense = "max_precision",
+             f = NULL, df = NULL, Hf = NULL,
+             init_w = w,
+             mc_cores = 1L, pm_tol = 1e-05,
+             max_iters = 100L, print_pm = FALSE)
> 
> mosalloc(D = D, d = d, C = C, c = c, l = l, u = u, opts = opts)
$w
[1] 1.00 3.85 3.80 1.30

$n
[1] 755.8117 499.7609 241.5215 691.6864 584.9882

$J
[1] 0.0018907289 0.0004917598 0.0004982304 0.0014563656

$Objective
NULL

$Utiopian
NULL

$Normal
[1] 2.363801e-08 2.237152e+00 1.186254e+00 1.387742e-01

$dfJ
NULL

$Sensitivity
$Sensitivity$D
[1] 2.363801e-08 5.810784e-01 3.121722e-01 1.067494e-01

$Sensitivity$A
numeric(0)

$Sensitivity$C
[1]  1.050676e-08 -7.307475e-09

$Sensitivity$lbox
[1] 9.093240e-17 1.547014e-16 4.498037e-16 1.019954e-16 1.253083e-16

$Sensitivity$ubox
[1] 1.129080e-17 5.617199e-18 1.503004e-17 2.807207e-18 1.251986e-17


$Qbounds
[1] 1 1 1 1 1

$Dbounds
[1] 0.001890729 0.001893275 0.001893275 0.001893275

$Scalingfactor
[1] 0.0016654338 0.0020770702 0.0022015450 0.0002711488 0.0005464200

$Ecosolver
$Ecosolver$Ecoinfostring
[1] "Optimal solution found"

$Ecosolver$Ecoredcodes
exitFlag     iter  mi_iter   numerr 
       0       12       -1        0 

$Ecosolver$Ecosummary
       pcost        dcost         pres         dres         pinf         dinf 
1.069071e+00 1.069071e+00 1.191108e-12 4.899284e-12 0.000000e+00 0.000000e+00 
     pinfres      dinfres          gap       relgap           r0 
7.107497e-02          NaN 2.309435e-09 2.160227e-09 1.000000e-10 


$Timing
     TotalTime InnerTime ECOS_runtime ECOS_tsetup ECOS_tsolve
[1,]     0.001         0   8.9041e-05   8.625e-06  8.0416e-05

$Iteration
NULL

> 
> #----------------------------------------------------------------------------
> # Example 4: Example 2 with multiple preference weightings for simultaneous
> # evaluation
> 
> w <- matrix(c(1.0, 1.0, 1.0, 1.0,       # matrix of preference weightings
+               1.0, 3.9, 3.9, 1.3,
+               0.8, 4.2, 4.8, 1.5,
+               1.2, 3.5, 4.8, 2.0,
+               2.0, 1.0, 1.0, 2.0), 5, 4, byrow = TRUE)
> w <- w / w[,1]     # rescale w (ensure the first weighting to be one)
> l <- rep(n.min, 5) # minimum sample size ber stratum
> u <- Nh            # maximum sample size per stratum
> C <- rbind(ch, ch * c(-1, -1, -1, 0, 0))
> c <- c(budget, - 0.5 * budget)
> A <- NULL # no precision constraint
> a <- NULL # no precision constraint
> 
> D <- rbind(Sh.rev**2 * Nh**2/sum(Nh * mh.rev)**2,
+            Sh.emp**2 * Nh**2/sum(Nh * mh.emp)**2,
+            ph.rsch * (1 - ph.rsch) * Nh**3/(Nh - 1)/sum(Nh * ph.rsch)**2,
+            ph.offsh * (1 - ph.offsh) * Nh**3/(Nh - 1)/sum(Nh * ph.offsh)**2)
> 
> d <- as.vector(D %*% (1 / Nh))
> 
> opts = list(sense = "max_precision",
+             f = NULL, df = NULL, Hf = NULL,
+             init_w = w,
+             mc_cores = 1L, pm_tol = 1e-05,
+             max_iters = 100L, print_pm = FALSE)
> 
> sols <- mosalloc(D = D, d = d, C = C, c = c, l = l, u = u, opts = opts)
> lapply(sols, function(sol){sol$Qbounds})
[[1]]
[1] 1 1 1 1 1

[[2]]
[1] 1 1 1 1 1

[[3]]
[1] 1 1 1 1 1

[[4]]
[1] 1 1 1 1 1

[[5]]
[1] 1 1 1 1 1

> 
> #----------------------------------------------------------------------------
> # Example 5: Example 2 where a weighted sum scalarization of the objective
> # components is minimized
> 
> l <- rep(n.min, 5) # minimum sample size ber stratum
> u <- Nh            # maximum sample size per stratum
> C <- matrix(ch, nrow = 1)
> c <- budget
> A <- NULL # no precision constraint
> a <- NULL # no precision constraint
> 
> # Objective variance components
> D <- rbind(Sh.rev**2 * Nh**2/sum(Nh * mh.rev)**2,
+            Sh.emp**2 * Nh**2/sum(Nh * mh.emp)**2,
+            ph.rsch * (1 - ph.rsch) * Nh**3/(Nh - 1)/sum(Nh * ph.rsch)**2,
+            ph.offsh * (1 - ph.offsh) * Nh**3/(Nh - 1)/sum(Nh * ph.offsh)**2)
> 
> d <- as.vector(D %*% (1 / Nh)) # finite population correction
> 
> # Simple weighted sum as decision functional
> wss <- c(1, 1, 0.5, 0.5) # preference weighting (weighted sum scalarization)
> 
> Dw <- wss %*% D
> dw <- as.vector(wss %*% d)
> 
> opts = list(sense = "max_precision",
+             f = NULL, df = NULL, Hf = NULL,
+             init_w = 1,
+             mc_cores = 1L, pm_tol = 1e-05,
+             max_iters = 1000L, print_pm = FALSE)
> 
> # Solve weighted sum scalarization (WSS) via mosalloc
> sol_wss <- mosalloc(D = Dw, d = dw, C = C, c = c, l = l, u = u, opts = opts)
Error in mosalloc(D = Dw, d = dw, C = C, c = c, l = l, u = u, opts = opts) : 
  d is not an utopian vector! d is too large.
Execution halted