Nothing
if(!require("GNE"))stop("this test requires package GNE.")
itermax <- 10
#-------------------------------------------------------------------------------
# (1) Example 5 of von Facchinei et al. (2007)
#-------------------------------------------------------------------------------
dimx <- c(1, 1)
#Gr_x_j O_i(x)
grobj <- function(x, i, j)
{
if(i == 1)
res <- c(2*(x[1]-1), 0)
if(i == 2)
res <- c(0, 2*(x[2]-1/2))
res[j]
}
#Gr_x_k Gr_x_j O_i(x)
heobj <- function(x, i, j, k)
2 * (i == j && j == k)
dimlam <- c(1, 1)
#constraint function g_i(x)
g <- function(x, i)
sum(x[1:2]) - 1
#Gr_x_j g_i(x)
grg <- function(x, i, j)
1
#Gr_x_k Gr_x_j g_i(x)
heg <- function(x, i, j, k)
0
x0 <- rep(0, sum(dimx))
z0 <- c(x0, 2, 2, max(10, 5-g(x0, 1) ), max(10, 5-g(x0, 2) ) )
#true value is (3/4, 1/4, 1/2, 1/2)
GNE.ceq(z0, dimx, dimlam, grobj=grobj, heobj=heobj,
constr=g, grconstr=grg, heconstr=heg, method="PR",
control=list(trace=0, maxit=itermax))
GNE.ceq(z0, dimx, dimlam, grobj=grobj, heobj=heobj,
constr=g, grconstr=grg, heconstr=heg, method="AS", global="pwldog",
xscalm="auto", control=list(trace=0, maxit=itermax))
#-------------------------------------------------------------------------------
# (2) Duopoly game of Krawczyk and Stanislav Uryasev (2000)
#-------------------------------------------------------------------------------
#constants
myarg <- list(d= 20, lambda= 4, rho= 1)
dimx <- c(1, 1)
#Gr_x_j O_i(x)
grobj <- function(x, i, j, arg)
{
res <- -arg$rho * x[i]
if(i == j)
res <- res + arg$d - arg$lambda - arg$rho*(x[1]+x[2])
-res
}
#Gr_x_k Gr_x_j O_i(x)
heobj <- function(x, i, j, k, arg)
arg$rho * (i == j) + arg$rho * (j == k)
dimlam <- c(1, 1)
#constraint function g_i(x)
g <- function(x, i)
-x[i]
#Gr_x_j g_i(x)
grg <- function(x, i, j)
-1*(i == j)
#Gr_x_k Gr_x_j g_i(x)
heg <- function(x, i, j, k)
0
#true value is (16/3, 16/3, 0, 0)
x0 <- rep(0, sum(dimx))
z0 <- c(x0, 2, 2, max(10, 5-g(x0, 1) ), max(10, 5-g(x0, 2) ) )
GNE.ceq(z0, dimx, dimlam, grobj=grobj, heobj=heobj, arggrobj=myarg,
argheobj=myarg, constr=g, grconstr=grg, heconstr=heg,
method="PR", control=list(trace=0, maxit=itermax))
GNE.ceq(z0, dimx, dimlam, grobj=grobj, heobj=heobj, arggrobj=myarg,
argheobj=myarg, constr=g, grconstr=grg, heconstr=heg,
method="AS", global="pwldog", xscalm="auto", control=list(trace=1, maxit=itermax))
#-------------------------------------------------------------------------------
# (3) River basin pollution game of Krawczyk and Stanislav Uryasev (2000)
#-------------------------------------------------------------------------------
myarg <- list(
C = cbind(c(.1, .12, .15), c(.01, .05, .01)),
U = cbind(c(6.5, 5, 5.5), c(4.583, 6.25, 3.75)),
K = c(100, 100),
E = c(.5, .25, .75),
D = c(3, .01)
)
dimx <- c(1, 1, 1)
#Gr_x_j O_i(x)
grobj <- function(x, i, j, arg)
{
dij <- 1*(i == j)
res <- -(-arg$D[2] - arg$C[i, 2]*dij) * x[i]
res - (arg$D[1] - arg$D[2]*sum(x[1:3]) - arg$C[i, 1] - arg$C[i, 2]*x[i]) * dij
}
#Gr_x_k Gr_x_j O_i(x)
heobj <- function(x, i, j, k, arg)
{
dij <- 1*(i == j)
dik <- 1*(i == k)
arg$D[2] * dik + arg$D[2] * dij + 2 * arg$C[i, 2] * dij * dik
}
dimlam <- c(2, 2, 2)
#g_i(x)
g <- function(x, i, arg)
c(sum(arg$U[, 1] * arg$E * x[1:3]) - arg$K[1],
sum(arg$U[, 2] * arg$E * x[1:3]) - arg$K[2])
#Gr_x_j g_i(x)
grg <- function(x, i, j, arg)
c(arg$U[j, 1] * arg$E[j], arg$U[j, 2] * arg$E[j])
#Gr_x_k Gr_x_j g_i(x)
heg <- function(x, i, j, k, arg)
c(0, 0)
#true value around (21.146, 16.027, 2.724, 0.574, 0.000)
x0 <- rep(0, sum(dimx))
z0 <- c(x0, rep(2, sum(dimlam)), pmax(10, 5-g(x0, 1, myarg) ),
pmax(10, 5-g(x0, 2, myarg) ), pmax(10, 5-g(x0, 3, myarg) ) )
# funCER(z0, dimx, dimlam, grobj=grobj,
# constr=g, grconstr=grg,
# arggrobj=myarg, argconstr=myarg,
# arggrconstr=myarg)
# jacCER(z0, dimx, dimlam, heobj=heobj,
# constr=g, grconstr=grg, heconstr=heg,
# argheobj=myarg, argconstr=myarg,
# arggrconstr=myarg, argheconstr=myarg)
GNE.ceq(z0, dimx, dimlam, grobj=grobj, heobj=heobj,
constr=g, grconstr=grg, heconstr=heg, method="PR", global="gline",
arggrobj=myarg, argheobj=myarg, argconstr=myarg,
arggrconstr=myarg, argheconstr=myarg, control=list(trace=1, maxit=itermax),
silent=TRUE)
GNE.ceq(z0, dimx, dimlam, grobj=grobj, heobj=heobj,
constr=g, grconstr=grg, heconstr=heg,
arggrobj=myarg, argheobj=myarg, argconstr=myarg,
arggrconstr=myarg, argheconstr=myarg, method="AS", global="pwldog",
xscalm="auto", control=list(trace=1, maxit=itermax))
#-------------------------------------------------------------------------------
# (4) Example of GNE with 4 solutions(!)
#-------------------------------------------------------------------------------
myarg <- list(C=c(2, 3), D=c(4,0))
dimx <- c(1, 1)
#Gr_x_j O_i(x)
grobj <- function(x, i, j, arg)
{
dij <- 1*(i == j)
other <- ifelse(i == 1, 2, 1)
2*(x[i] - arg$C[i])*(x[other] - arg$D[i])^4*dij + 4*(x[i] - arg$C[i])^2*(x[other] - arg$D[i])^3*(1-dij)
}
#Gr_x_k Gr_x_j O_i(x)
heobj <- function(x, i, j, k, arg)
{
dij <- 1*(i == j)
dik <- 1*(i == k)
other <- ifelse(i == 1, 2, 1)
res <- 2*(x[other] - arg$D[i])^4*dij*dik + 8*(x[i] - arg$C[i])*(x[other] - arg$D[i])^3*dij*(1-dik)
res <- res + 8*(x[i] - arg$C[i])*(x[other] - arg$D[i])^3*(1-dij)*dik
res + 12*(x[i] - arg$C[i])^2*(x[other] - arg$D[i])^2*(1-dij)*(1-dik)
}
dimlam <- c(1, 1)
#g_i(x)
g <- function(x, i)
ifelse(i == 1, sum(x[1:2]) - 1, 2*x[1]+x[2]-2)
#Gr_x_j g_i(x)
grg <- function(x, i, j)
ifelse(i == 1, 1, 1 + 1*(i == j))
#Gr_x_k Gr_x_j g_i(x)
heg <- function(x, i, j, k)
0
x0 <- c(2, 2)
z0 <- c(x0, 10, 10, max(10, 5-g(x0, 1) ), max(10, 5-g(x0, 2) ) )
GNE.ceq(z0, dimx, dimlam, grobj=grobj, heobj=heobj,
constr=g, grconstr=grg, heconstr=heg, arggrobj=myarg, argheobj=myarg,
method="PR", control=list(trace=0, maxit=itermax))
GNE.ceq(z0, dimx, dimlam, grobj=grobj, heobj=heobj,
constr=g, grconstr=grg, heconstr=heg, arggrobj=myarg, argheobj=myarg,
method="AS", global="pwldog",
xscalm="auto", control=list(trace=0, maxit=itermax))
#-------------------------------------------------------------------------------
# (5) another Example
#-------------------------------------------------------------------------------
# associated objective functions
dimx <- c(2, 2, 3)
#Gr_x_j O_i(x)
grfullob <- function(x, i, j)
{
x <- x[1:7]
if(i == 1)
grad <- 3*(x - 1:7)^2
if(i == 2)
grad <- 1:7*(x - 1:7)^(0:6)
if(i == 3)
{
s <- x[5]^2 + x[6]^2 + x[7]^2 - 5
grad <- c(1, 0, 1, 0, 4*x[5]*s, 4*x[6]*s, 4*x[7]*s)
}
grad[j]
}
#Gr_x_k Gr_x_j O_i(x)
hefullob <- function(x, i, j, k)
{
x <- x[1:7]
if(i == 1)
he <- diag( 6*(x - 1:7) )
if(i == 2)
he <- diag( c(0, 2, 6, 12, 20, 30, 42)*(x - 1:7)^c(0, 0:5) )
if(i == 3)
{
s <- x[5]^2 + x[6]^2 + x[7]^2
he <- rbind(rep(0, 7), rep(0, 7), rep(0, 7), rep(0, 7),
c(0, 0, 0, 0, 4*s+8*x[5]^2, 8*x[5]*x[6], 8*x[5]*x[7]),
c(0, 0, 0, 0, 8*x[5]*x[6], 4*s+8*x[6]^2, 8*x[6]*x[7]),
c(0, 0, 0, 0, 8*x[5]*x[7], 8*x[6]*x[7], 4*s+8*x[7]^2) )
}
he[j,k]
}
# constraint linked functions
dimlam <- c(1, 2, 2)
#constraint function g_i(x)
g <- function(x, i)
{
x <- x[1:7]
if(i == 1)
res <- sum( x^(1:7) ) -7
if(i == 2)
res <- c(sum(x) + prod(x) - 14, 20 - sum(x))
if(i == 3)
res <- c(sum(x^2) + 1, 100 - sum(x))
res
}
#Gr_x_j g_i(x)
grfullg <- function(x, i, j)
{
x <- x[1:7]
if(i == 1)
grad <- (1:7) * x ^ (0:6)
if(i == 2)
{
grad <- 1 + sapply(1:7, function(i) prod(x[-i]))
grad <- cbind(grad, -1)
}
if(i == 3)
grad <- cbind(2*x, -1)
if(i == 1)
res <- grad[j]
if(i != 1)
res <- grad[j,]
as.numeric(res)
}
#Gr_x_k Gr_x_j g_i(x)
hefullg <- function(x, i, j, k)
{
x <- x[1:7]
if(i == 1)
he1 <- diag( c(0, 2, 6, 12, 20, 30, 42) * x ^ c(0, 0, 1:5) )
if(i == 2)
{
he1 <- matrix(0, 7, 7)
he1[1, -1] <- sapply(2:7, function(i) prod(x[-c(1, i)]))
he1[2, -2] <- sapply(c(1, 3:7), function(i) prod(x[-c(2, i)]))
he1[3, -3] <- sapply(c(1:2, 4:7), function(i) prod(x[-c(3, i)]))
he1[4, -4] <- sapply(c(1:3, 5:7), function(i) prod(x[-c(4, i)]))
he1[5, -5] <- sapply(c(1:4, 6:7), function(i) prod(x[-c(5, i)]))
he1[6, -6] <- sapply(c(1:5, 7:7), function(i) prod(x[-c(6, i)]))
he1[7, -7] <- sapply(1:6, function(i) prod(x[-c(7, i)]))
he2 <- matrix(0, 7, 7)
}
if(i == 3)
{
he1 <- diag(rep(2, 7))
he2 <- matrix(0, 7, 7)
}
if(i != 1)
return( c(he1[j, k], he2[j, k]) )
else
return( he1[j, k] )
}
# (3) compute Phi
#
x0 <- rep(0, sum(dimx))
z0 <- c(x0, rep(2, sum(dimlam)), max(10, 5-g(x0, 1) ),
pmax(10, 5-g(x0, 2) ), pmax(10, 5-g(x0, 3) ) )
GNE.ceq(z0, dimx, dimlam, grobj=grfullob, heobj=hefullob,
constr=g, grconstr=grfullg, heconstr=hefullg,
method="PR", control=list(trace=0, maxit=itermax))
GNE.ceq(z0, dimx, dimlam, grobj=grfullob, heobj=hefullob,
constr=g, grconstr=grfullg, heconstr=hefullg,
method="AS", global="pwldog",
xscalm="auto", control=list(trace=3, maxit=itermax))
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.