```import numpy as np

from pymoo.model.problem import Problem
from pymoo.problems.many import generic_sphere
from pymoo.util.misc import powerset

class WFG(Problem):

def __init__(self, n_var, n_obj, k=None, l=None, **kwargs):
super().__init__(n_var=n_var,
n_obj=n_obj,
n_constr=0,
xl=0,
xu=2 * np.arange(1, n_var + 1),
type_var=np.double,
**kwargs)

self.S = np.arange(2, 2 * self.n_obj + 1, 2)
self.A = np.ones(self.n_obj - 1)

if k:
self.k = k
else:
if n_obj == 2:
self.k = 4
else:
self.k = 2 * (n_obj - 1)

if l:
self.l = l
else:
self.l = n_var - self.k

self.validate(self.l, self.k, self.n_obj)

def validate(self, l, k, n_obj):
if n_obj < 2:
raise ValueError('WFG problems must have two or more objectives.')
if not k % (n_obj - 1) == 0:
raise ValueError('Position parameter (k) must be divisible by number of objectives minus one.')
if k < 4:
raise ValueError('Position parameter (k) must be greater or equal than 4.')
if (k + l) < n_obj:
raise ValueError('Sum of distance and position parameters must be greater than num. of objs. (k + l >= M).')

def _post(self, t, a):
x = []
for i in range(t.shape[1] - 1):
x.append(np.maximum(t[:, -1], a[i]) * (t[:, i] - 0.5) + 0.5)
x.append(t[:, -1])
return np.column_stack(x)

def _calculate(self, x, s, h):
return x[:, -1][:, None] + s * np.column_stack(h)

def _rand_optimal_position(self, n):
return np.random.random((n, self.k))

def _positional_to_optimal(self, K):
suffix = np.full((len(K), self.l), 0.35)
X = np.column_stack([K, suffix])
return X * self.xu

def _calc_pareto_set(self, n_pareto_points=500, *args, **kwargs):
ps = np.ones((2 ** self.k, self.k))
for i, s in enumerate(powerset(np.arange(self.k))):
ps[i, s] = 0

rnd = self._rand_optimal_position(n_pareto_points - len(ps))
ps = np.row_stack([ps, rnd])
ps = self._positional_to_optimal(ps)
return ps

def _calc_pareto_front(self, *args, **kwargs):
# ps = self.pareto_set(n_pareto_points=n_pareto_points)
# return self.evaluate(ps, return_values_of=["F"])
return None

class WFG1(WFG):

@staticmethod
def t1(x, n, k):
x[:, k:n] = _transformation_shift_linear(x[:, k:n], 0.35)
return x

@staticmethod
def t2(x, n, k):
x[:, k:n] = _transformation_bias_flat(x[:, k:n], 0.8, 0.75, 0.85)
return x

@staticmethod
def t3(x, n):
x[:, :n] = _transformation_bias_poly(x[:, :n], 0.02)
return x

@staticmethod
def t4(x, m, n, k):
w = np.arange(2, 2 * n + 1, 2)
gap = k // (m - 1)
t = []
for m in range(1, m):
_y = x[:, (m - 1) * gap: (m * gap)]
_w = w[(m - 1) * gap: (m * gap)]
t.append(_reduction_weighted_sum(_y, _w))
t.append(_reduction_weighted_sum(x[:, k:n], w[k:n]))
return np.column_stack(t)

def _evaluate(self, x, out, *args, **kwargs):
y = x / self.xu
y = WFG1.t1(y, self.n_var, self.k)
y = WFG1.t2(y, self.n_var, self.k)
y = WFG1.t3(y, self.n_var)
y = WFG1.t4(y, self.n_obj, self.n_var, self.k)

y = self._post(y, self.A)

h = [_shape_convex(y[:, :-1], m + 1) for m in range(self.n_obj - 1)]
h.append(_shape_mixed(y[:, 0], alpha=1.0, A=5.0))

out["F"] = self._calculate(y, self.S, h)

def _rand_optimal_position(self, n):
return np.power(np.random.random((n, self.k)), 50.0)

class WFG2(WFG):

def validate(self, l, k, n_obj):
super().validate(l, k, n_obj)
validate_wfg2_wfg3(l)

@staticmethod
def t2(x, n, k):
y = [x[:, i] for i in range(k)]

l = n - k
ind_non_sep = k + l // 2

i = k + 1
while i <= ind_non_sep:
head = k + 2 * (i - k) - 2
tail = k + 2 * (i - k)
i += 1

return np.column_stack(y)

@staticmethod
def t3(x, m, n, k):
ind_r_sum = k + (n - k) // 2
gap = k // (m - 1)

t = [_reduction_weighted_sum_uniform(x[:, (m - 1) * gap: (m * gap)]) for m in range(1, m)]
t.append(_reduction_weighted_sum_uniform(x[:, k:ind_r_sum]))

return np.column_stack(t)

def _evaluate(self, x, out, *args, **kwargs):
y = x / self.xu
y = WFG1.t1(y, self.n_var, self.k)
y = WFG2.t2(y, self.n_var, self.k)
y = WFG2.t3(y, self.n_obj, self.n_var, self.k)
y = self._post(y, self.A)

h = [_shape_convex(y[:, :-1], m + 1) for m in range(self.n_obj - 1)]
h.append(_shape_disconnected(y[:, 0], alpha=1.0, beta=1.0, A=5.0))

out["F"] = self._calculate(y, self.S, h)

class WFG3(WFG):

def __init__(self, n_var, n_obj, k=None, **kwargs):
super().__init__(n_var, n_obj, k=k, **kwargs)
self.A[1:] = 0

def validate(self, l, k, n_obj):
super().validate(l, k, n_obj)
validate_wfg2_wfg3(l)

def _evaluate(self, x, out, *args, **kwargs):
y = x / self.xu
y = WFG1.t1(y, self.n_var, self.k)
y = WFG2.t2(y, self.n_var, self.k)
y = WFG2.t3(y, self.n_obj, self.n_var, self.k)
y = self._post(y, self.A)

h = [_shape_linear(y[:, :-1], m + 1) for m in range(self.n_obj)]

out["F"] = self._calculate(y, self.S, h)

# from optproblems.wfg import WFG3 as WFG3opt
# out["F"] =  np.array([WFG3opt(self.n_obj, self.n_var, self.k).objective_function(_x) for _x in x])

class WFG4(WFG):

@staticmethod
def t1(x):
return _transformation_shift_multi_modal(x, 30.0, 10.0, 0.35)

@staticmethod
def t2(x, m, k):
gap = k // (m - 1)
t = [_reduction_weighted_sum_uniform(x[:, (m - 1) * gap: (m * gap)]) for m in range(1, m)]
t.append(_reduction_weighted_sum_uniform(x[:, k:]))
return np.column_stack(t)

def _evaluate(self, x, out, *args, **kwargs):
y = x / self.xu
y = WFG4.t1(y)
y = WFG4.t2(y, self.n_obj, self.k)
y = self._post(y, self.A)

h = [_shape_concave(y[:, :-1], m + 1) for m in range(self.n_obj)]

out["F"] = self._calculate(y, self.S, h)

def _calc_pareto_front(self, ref_dirs):
return generic_sphere(ref_dirs) * self.S

class WFG5(WFG):

@staticmethod
def t1(x):
return _transformation_param_deceptive(x, A=0.35, B=0.001, C=0.05)

def _evaluate(self, x, out, *args, **kwargs):
y = x / self.xu
y = WFG5.t1(y)
y = WFG4.t2(y, self.n_obj, self.k)
y = self._post(y, self.A)

h = [_shape_concave(y[:, :-1], m + 1) for m in range(self.n_obj)]

out["F"] = self._calculate(y, self.S, h)

def _calc_pareto_front(self, ref_dirs):
return generic_sphere(ref_dirs) * self.S

class WFG6(WFG):

@staticmethod
def t2(x, m, n, k):
gap = k // (m - 1)
t = [_reduction_non_sep(x[:, (m - 1) * gap: (m * gap)], gap) for m in range(1, m)]
t.append(_reduction_non_sep(x[:, k:], n - k))
return np.column_stack(t)

def _evaluate(self, x, out, *args, **kwargs):
y = x / self.xu
y = WFG1.t1(y, self.n_var, self.k)
y = WFG6.t2(y, self.n_obj, self.n_var, self.k)
y = self._post(y, self.A)

h = [_shape_concave(y[:, :-1], m + 1) for m in range(self.n_obj)]

out["F"] = self._calculate(y, self.S, h)

def _calc_pareto_front(self, ref_dirs):
return generic_sphere(ref_dirs) * self.S

class WFG7(WFG):

@staticmethod
def t1(x, k):
for i in range(k):
aux = _reduction_weighted_sum_uniform(x[:, i + 1:])
x[:, i] = _transformation_param_dependent(x[:, i], aux)
return x

def _evaluate(self, x, out, *args, **kwargs):
y = x / self.xu
y = WFG7.t1(y, self.k)
y = WFG1.t1(y, self.n_var, self.k)
y = WFG4.t2(y, self.n_obj, self.k)
y = self._post(y, self.A)

h = [_shape_concave(y[:, :-1], m + 1) for m in range(self.n_obj)]

out["F"] = self._calculate(y, self.S, h)

def _calc_pareto_front(self, ref_dirs):
return generic_sphere(ref_dirs) * self.S

class WFG8(WFG):

@staticmethod
def t1(x, n, k):
ret = []
for i in range(k, n):
aux = _reduction_weighted_sum_uniform(x[:, :i])
ret.append(_transformation_param_dependent(x[:, i], aux, A=0.98 / 49.98, B=0.02, C=50.0))
return np.column_stack(ret)

def _evaluate(self, x, out, *args, **kwargs):
y = x / self.xu
y[:, self.k:self.n_var] = WFG8.t1(y, self.n_var, self.k)
y = WFG1.t1(y, self.n_var, self.k)
y = WFG4.t2(y, self.n_obj, self.k)
y = self._post(y, self.A)

h = [_shape_concave(y[:, :-1], m + 1) for m in range(self.n_obj)]

out["F"] = self._calculate(y, self.S, h)

def _calc_pareto_front(self, ref_dirs):
return generic_sphere(ref_dirs) * self.S

def _positional_to_optimal(self, K):
k, l = self.k, self.l

for i in range(k, k + l):
u = K.sum(axis=1) / K.shape[1]
tmp1 = np.abs(np.floor(0.5 - u) + 0.98 / 49.98)
tmp2 = 0.02 + 49.98 * (0.98 / 49.98 - (1.0 - 2.0 * u) * tmp1)
suffix = np.power(0.35, np.power(tmp2, -1.0))

K = np.column_stack([K, suffix[:, None]])

ret = K * (2 * (np.arange(self.n_var) + 1))
return ret

class WFG9(WFG):

@staticmethod
def t1(x, n):
ret = []
for i in range(0, n - 1):
aux = _reduction_weighted_sum_uniform(x[:, i + 1:])
ret.append(_transformation_param_dependent(x[:, i], aux))
return np.column_stack(ret)

@staticmethod
def t2(x, n, k):
a = [_transformation_shift_deceptive(x[:, i], 0.35, 0.001, 0.05) for i in range(k)]
b = [_transformation_shift_multi_modal(x[:, i], 30.0, 95.0, 0.35) for i in range(k, n)]
return np.column_stack(a + b)

@staticmethod
def t3(x, m, n, k):
gap = k // (m - 1)
t = [_reduction_non_sep(x[:, (m - 1) * gap: (m * gap)], gap) for m in range(1, m)]
t.append(_reduction_non_sep(x[:, k:], n - k))
return np.column_stack(t)

def _evaluate(self, x, out, *args, **kwargs):
y = x / self.xu
y[:, :self.n_var - 1] = WFG9.t1(y, self.n_var)
y = WFG9.t2(y, self.n_var, self.k)
y = WFG9.t3(y, self.n_obj, self.n_var, self.k)

h = [_shape_concave(y[:, :-1], m + 1) for m in range(self.n_obj)]

out["F"] = self._calculate(y, self.S, h)

def _calc_pareto_front(self, ref_dirs):
return generic_sphere(ref_dirs) * self.S

def _positional_to_optimal(self, K):
k, l = self.k, self.l

suffix = np.full((len(K), self.l), 0.0)
X = np.column_stack([K, suffix])
X[:, self.k + self.l - 1] = 0.35

for i in range(self.k + self.l - 2, self.k - 1, -1):
m = X[:, i + 1:k + l]
val = m.sum(axis=1) / m.shape[1]
X[:, i] = 0.35 ** ((0.02 + 1.96 * val) ** -1)

ret = X * (2 * (np.arange(self.n_var) + 1))
return ret

# ---------------------------------------------------------------------------------------------------------
# TRANSFORMATIONS
# ---------------------------------------------------------------------------------------------------------

def _transformation_shift_linear(value, shift=0.35):
return correct_to_01(np.fabs(value - shift) / np.fabs(np.floor(shift - value) + shift))

def _transformation_shift_deceptive(y, A=0.35, B=0.005, C=0.05):
tmp1 = np.floor(y - A + B) * (1.0 - C + (A - B) / B) / (A - B)
tmp2 = np.floor(A + B - y) * (1.0 - C + (1.0 - A - B) / B) / (1.0 - A - B)
ret = 1.0 + (np.fabs(y - A) - B) * (tmp1 + tmp2 + 1.0 / B)
return correct_to_01(ret)

def _transformation_shift_multi_modal(y, A, B, C):
tmp1 = np.fabs(y - C) / (2.0 * (np.floor(C - y) + C))
tmp2 = (4.0 * A + 2.0) * np.pi * (0.5 - tmp1)
ret = (1.0 + np.cos(tmp2) + 4.0 * B * np.power(tmp1, 2.0)) / (B + 2.0)
return correct_to_01(ret)

def _transformation_bias_flat(y, a, b, c):
ret = a + np.minimum(0, np.floor(y - b)) * (a * (b - y) / b) \
- np.minimum(0, np.floor(c - y)) * ((1.0 - a) * (y - c) / (1.0 - c))
return correct_to_01(ret)

def _transformation_bias_poly(y, alpha):
return correct_to_01(y ** alpha)

def _transformation_param_dependent(y, y_deg, A=0.98 / 49.98, B=0.02, C=50.0):
aux = A - (1.0 - 2.0 * y_deg) * np.fabs(np.floor(0.5 - y_deg) + A)
ret = np.power(y, B + (C - B) * aux)
return correct_to_01(ret)

def _transformation_param_deceptive(y, A=0.35, B=0.001, C=0.05):
tmp1 = np.floor(y - A + B) * (1.0 - C + (A - B) / B) / (A - B)
tmp2 = np.floor(A + B - y) * (1.0 - C + (1.0 - A - B) / B) / (1.0 - A - B)
ret = 1.0 + (np.fabs(y - A) - B) * (tmp1 + tmp2 + 1.0 / B)
return correct_to_01(ret)

# ---------------------------------------------------------------------------------------------------------
# REDUCTION
# ---------------------------------------------------------------------------------------------------------

def _reduction_weighted_sum(y, w):
return correct_to_01(np.dot(y, w) / w.sum())

def _reduction_weighted_sum_uniform(y):
return correct_to_01(y.mean(axis=1))

def _reduction_non_sep(y, A):
n, m = y.shape
val = np.ceil(A / 2.0)

num = np.zeros(n)
for j in range(m):
num += y[:, j]
for k in range(A - 1):
num += np.fabs(y[:, j] - y[:, (1 + j + k) % m])

denom = m * val * (1.0 + 2.0 * A - 2 * val) / A

return correct_to_01(num / denom)

# ---------------------------------------------------------------------------------------------------------
# SHAPE
# ---------------------------------------------------------------------------------------------------------

def _shape_concave(x, m):
M = x.shape[1]
if m == 1:
ret = np.prod(np.sin(0.5 * x[:, :M] * np.pi), axis=1)
elif 1 < m <= M:
ret = np.prod(np.sin(0.5 * x[:, :M - m + 1] * np.pi), axis=1)
ret *= np.cos(0.5 * x[:, M - m + 1] * np.pi)
else:
ret = np.cos(0.5 * x[:, 0] * np.pi)
return correct_to_01(ret)

def _shape_convex(x, m):
M = x.shape[1]
if m == 1:
ret = np.prod(1.0 - np.cos(0.5 * x[:, :M] * np.pi), axis=1)
elif 1 < m <= M:
ret = np.prod(1.0 - np.cos(0.5 * x[:, :M - m + 1] * np.pi), axis=1)
ret *= 1.0 - np.sin(0.5 * x[:, M - m + 1] * np.pi)
else:
ret = 1.0 - np.sin(0.5 * x[:, 0] * np.pi)
return correct_to_01(ret)

def _shape_linear(x, m):
M = x.shape[1]
if m == 1:
ret = np.prod(x, axis=1)
elif 1 < m <= M:
ret = np.prod(x[:, :M - m + 1], axis=1)
ret *= 1.0 - x[:, M - m + 1]
else:
ret = 1.0 - x[:, 0]
return correct_to_01(ret)

def _shape_mixed(x, A=5.0, alpha=1.0):
aux = 2.0 * A * np.pi
ret = np.power(1.0 - x - (np.cos(aux * x + 0.5 * np.pi) / aux), alpha)
return correct_to_01(ret)

def _shape_disconnected(x, alpha=1.0, beta=1.0, A=5.0):
aux = np.cos(A * np.pi * x ** beta)
return correct_to_01(1.0 - x ** alpha * aux ** 2)

# ---------------------------------------------------------------------------------------------------------
# UTIL
# ---------------------------------------------------------------------------------------------------------

def validate_wfg2_wfg3(l):
if not l % 2 == 0:
raise ValueError('In WFG2/WFG3 the distance-related parameter (l) must be divisible by 2.')

def correct_to_01(X, epsilon=1.0e-10):
X[np.logical_and(X < 0, X >= 0 - epsilon)] = 0
X[np.logical_and(X > 1, X <= 1 + epsilon)] = 1
return X
```