Skip to content

Commit

Permalink
apply black
Browse files Browse the repository at this point in the history
  • Loading branch information
luk036 committed Feb 12, 2024
1 parent 4cf885a commit ffff6a1
Show file tree
Hide file tree
Showing 44 changed files with 1,752 additions and 1,733 deletions.
4 changes: 2 additions & 2 deletions docs/conf.py
Original file line number Diff line number Diff line change
Expand Up @@ -169,7 +169,7 @@ def setup(app):
# The theme to use for HTML and HTML Help pages. See the documentation for
# a list of builtin themes.
# html_theme = 'alabaster'
html_theme = 'sphinx_rtd_theme'
html_theme = "sphinx_rtd_theme"

# Theme options are theme-specific and customize the look and feel of a theme
# further. For a list of options available for each theme, see the
Expand Down Expand Up @@ -301,4 +301,4 @@ def setup(app):
"pandas": ("https://pandas.pydata.org/pandas-docs/stable", None),
"scipy": ("https://docs.scipy.org/doc/scipy/reference", None),
"pyscaffold": ("https://pyscaffold.org/en/stable", None),
}
}
23 changes: 12 additions & 11 deletions experiments/corr_fn.py
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,7 @@
corr_bspline,
corr_poly,
create_2d_isotropic,
create_2d_sites
create_2d_sites,
)
from ellpy.oracles.lsq_corr_oracle import lsq_oracle
from ellpy.oracles.mle_corr_oracle import mle_oracle
Expand All @@ -29,15 +29,15 @@ def lsq_corr_core2(Y: Arr, n: int, P: Callable):
Returns:
[type]: [description]
"""
normY = np.linalg.norm(Y, 'fro')
normY = np.linalg.norm(Y, "fro")
normY2 = 32 * normY * normY
val = 256 * np.ones(n + 1)
val[-1] = normY2 * normY2
x = np.zeros(n + 1) # cannot all zeros
x[0] = 1.
x[0] = 1.0
x[-1] = normY2 / 2
E = ell(val, x)
xb, _, ell_info = cutting_plane_dc(P, E, float('inf'))
xb, _, ell_info = cutting_plane_dc(P, E, float("inf"))
return xb[:-1], ell_info.num_iters, ell_info.feasible


Expand Down Expand Up @@ -67,13 +67,13 @@ def mle_corr_core(Y: Arr, n: int, P: Callable):
[type]: [description]
"""
x = np.zeros(n)
x[0] = 1.
E = ell(50., x)
x[0] = 1.0
E = ell(50.0, x)
# E.use_parallel_cut = False
# options = Options()
# options.max_it = 2000
# options.tol = 1e-8
xb, _, ell_info = cutting_plane_dc(P, E, float('inf'))
xb, _, ell_info = cutting_plane_dc(P, E, float("inf"))
# print(num_iters, feasible, status)
return xb, ell_info.num_iters, ell_info.feasible

Expand Down Expand Up @@ -124,10 +124,11 @@ def lsq_corr_bspline2(Y, s, m):

if __name__ == "__main__":
import matplotlib.pyplot as plt

# import matplotlib.pylab as lab
s = create_2d_sites(10, 8)
Y = create_2d_isotropic(s, 1000)
print('start ell...')
print("start ell...")
spl, num_iters, _ = lsq_corr_bspline2(Y, s, 5)
pol, num_iters, _ = lsq_corr_poly2(Y, s, 5)
# pol, num_iters, _ = mle_corr_poly(Y, s, 4)
Expand All @@ -141,9 +142,9 @@ def lsq_corr_bspline2(Y, s, m):
# h = s[-1] - s[0]
d = np.sqrt(10**2 + 8**2)
xs = np.linspace(0, d, 100)
plt.plot(xs, spl(xs), 'g', label='BSpline')
plt.plot(xs, spl(xs), "g", label="BSpline")
# plt.plot(xs, splcvx(xs), 'b', label='BSpline CVX')
plt.plot(xs, np.polyval(pol, xs), 'r', label='Polynomial')
plt.plot(xs, np.polyval(pol, xs), "r", label="Polynomial")
# plt.plot(xs, np.polyval(polcvx, xs), 'r', label='Polynomial CVX')
plt.legend(loc='best')
plt.legend(loc="best")
plt.show()
21 changes: 11 additions & 10 deletions experiments/corr_fn_cvx.py
Original file line number Diff line number Diff line change
Expand Up @@ -8,7 +8,7 @@
from ellpy.oracles.corr_oracle import (
construct_distance_matrix,
create_2d_isotropic,
create_2d_sites
create_2d_sites,
)

Arr = Union[np.ndarray]
Expand Down Expand Up @@ -41,11 +41,11 @@ def lsq_corr_poly(Y: Arr, s: Arr, n: int):
D = np.multiply(D, D1)
Sig += D * a[n - 2 - i]
constraints = [Sig >> 0]
prob = cvx.Problem(cvx.Minimize(cvx.norm(Sig - Y, 'fro')), constraints)
prob = cvx.Problem(cvx.Minimize(cvx.norm(Sig - Y, "fro")), constraints)
prob.solve(solver=cvx.CVXOPT)
# prob.solve()
if prob.status != cvx.OPTIMAL:
raise Exception('CVXPY Error')
raise Exception("CVXPY Error")
return np.poly1d(np.array(a.value).flatten())


Expand Down Expand Up @@ -90,18 +90,19 @@ def lsq_corr_bspline(Y: Arr, s: Arr, n: int):
constraints = [Sig >> 0]
for i in range(n - 1):
constraints += [c[i] >= c[i + 1]]
constraints += [c[-1] >= 0.]
constraints += [c[-1] >= 0.0]

prob = cvx.Problem(cvx.Minimize(cvx.norm(Sig - Y, 'fro')), constraints)
prob = cvx.Problem(cvx.Minimize(cvx.norm(Sig - Y, "fro")), constraints)
prob.solve(solver=cvx.CVXOPT)
# prob.solve()
if prob.status != cvx.OPTIMAL:
raise Exception('CVXPY Error')
raise Exception("CVXPY Error")
return BSpline(t, np.array(c.value).flatten(), k)


if __name__ == "__main__":
import matplotlib.pyplot as plt

# import matplotlib.pylab as lab
s = create_2d_sites(10, 8)
Y = create_2d_isotropic(s, 1000)
Expand All @@ -111,16 +112,16 @@ def lsq_corr_bspline(Y: Arr, s: Arr, n: int):
# # pol, num_iters, _ = mle_corr_poly(Y, s, 4)
# print(pol)
# print(num_iters)
print('start cvx...')
print("start cvx...")
splcvx = lsq_corr_bspline(Y, s, 5)
polcvx = lsq_corr_poly(Y, s, 5)

# h = s[-1] - s[0]
d = np.sqrt(10**2 + 8**2)
xs = np.linspace(0, d, 100)
# plt.plot(xs, spl(xs), 'g', label='BSpline')
plt.plot(xs, splcvx(xs), 'b', label='BSpline CVX')
plt.plot(xs, splcvx(xs), "b", label="BSpline CVX")
# plt.plot(xs, np.polyval(pol, xs), 'r', label='Polynomial')
plt.plot(xs, np.polyval(polcvx, xs), 'r', label='Polynomial CVX')
plt.legend(loc='best')
plt.plot(xs, np.polyval(polcvx, xs), "r", label="Polynomial CVX")
plt.legend(loc="best")
plt.show()
21 changes: 11 additions & 10 deletions experiments/exp2d.py
Original file line number Diff line number Diff line change
Expand Up @@ -2,6 +2,7 @@
from __future__ import print_function

import matplotlib.pylab as lab

# from pprint import pprint
import matplotlib.pyplot as plt
import numpy as np
Expand All @@ -11,13 +12,13 @@
# from corr_fn import lsq_corr_poly, lsq_corr_bspline

# a fake dataset to make the bumps with
nx = 80 # number of points
nx = 80 # number of points
ny = 64
n = nx*ny
s_end = [10., 8.]
n = nx * ny
s_end = [10.0, 8.0]
sdkern = 0.25 # width of kernel
var = 2. # standard derivation
tau = 0.00001 # standard derivation of white noise
var = 2.0 # standard derivation
tau = 0.00001 # standard derivation of white noise
N = 4 # number of samples

# create sites s
Expand All @@ -39,19 +40,19 @@
# ym = np.random.randn(n)
for k in range(N):
x = var * np.random.randn(n)
y = A @ x + tau*np.random.randn(n)
y = A @ x + tau * np.random.randn(n)
Ys[:, k] = y

Y = np.cov(Ys, bias=True)

plt.subplot(2, 2, 1)
lab.contourf(xx, yy, np.reshape(Ys[:, 0], (ny, nx)), cmap='Greens')
lab.contourf(xx, yy, np.reshape(Ys[:, 0], (ny, nx)), cmap="Greens")
plt.subplot(2, 2, 2)
lab.contourf(xx, yy, np.reshape(Ys[:, 1], (ny, nx)), cmap='Greens')
lab.contourf(xx, yy, np.reshape(Ys[:, 1], (ny, nx)), cmap="Greens")
plt.subplot(2, 2, 3)
lab.contourf(xx, yy, np.reshape(Ys[:, 2], (ny, nx)), cmap='Greens')
lab.contourf(xx, yy, np.reshape(Ys[:, 2], (ny, nx)), cmap="Greens")
plt.subplot(2, 2, 4)
lab.contourf(xx, yy, np.reshape(Ys[:, 3], (ny, nx)), cmap='Greens')
lab.contourf(xx, yy, np.reshape(Ys[:, 3], (ny, nx)), cmap="Greens")
plt.show()


Expand Down
15 changes: 8 additions & 7 deletions src/ellpy/cutting_plane.py
Original file line number Diff line number Diff line change
Expand Up @@ -29,8 +29,7 @@ def __init__(self, feasible: bool, num_iters: int, status: CUTStatus):
self.status: CUTStatus = status


def cutting_plane_feas(Omega: Callable[[Any], Any], S,
options=Options()) -> CInfo:
def cutting_plane_feas(Omega: Callable[[Any], Any], S, options=Options()) -> CInfo:
"""Find a point in a convex set (defined through a cutting-plane oracle).
Description:
Expand Down Expand Up @@ -75,8 +74,9 @@ def cutting_plane_feas(Omega: Callable[[Any], Any], S,
return CInfo(feasible, niter + 1, status)


def cutting_plane_dc(Omega: Callable[[Any, Any], Any], S, t,
options=Options()) -> Tuple[Any, Any, CInfo]:
def cutting_plane_dc(
Omega: Callable[[Any, Any], Any], S, t, options=Options()
) -> Tuple[Any, Any, CInfo]:
"""Cutting-plane method for solving convex optimization problem
Arguments:
Expand Down Expand Up @@ -134,7 +134,7 @@ def cutting_plane_q(Omega, S, t, options=Options()):
status = CUTStatus.nosoln

for niter in range(options.max_it):
retry = (status == CUTStatus.noeffect)
retry = status == CUTStatus.noeffect
cut, x0, t1, more_alt = Omega(S.xc, t, retry)
if t1 is not None: # better t obtained
t = t1
Expand All @@ -153,8 +153,9 @@ def cutting_plane_q(Omega, S, t, options=Options()):
return x_best, t, ret


def bsearch(Omega: Callable[[Any], bool], Interval: Tuple,
options=Options()) -> Tuple[Any, CInfo]:
def bsearch(
Omega: Callable[[Any], bool], Interval: Tuple, options=Options()
) -> Tuple[Any, CInfo]:
"""[summary]
Arguments:
Expand Down
23 changes: 12 additions & 11 deletions src/ellpy/oracles/chol_ext.py
Original file line number Diff line number Diff line change
Expand Up @@ -10,15 +10,16 @@
class chol_ext:
"""LDLT factorization (mainly for LMI oracles)
- LDL^T square-root-free version
- Option allow semidefinite
- Cholesky–Banachiewicz style, row-based
- Lazy evaluation
- A matrix A in R^{m x m} is positive definite
iff v' A v > 0 for all v in R^n.
- O(p^3) per iteration, independent of N
- LDL^T square-root-free version
- Option allow semidefinite
- Cholesky–Banachiewicz style, row-based
- Lazy evaluation
- A matrix A in R^{m x m} is positive definite
iff v' A v > 0 for all v in R^n.
- O(p^3) per iteration, independent of N
"""
__slots__ = ('p', 'v', '_n', '_T', 'allow_semidefinite')

__slots__ = ("p", "v", "_n", "_T", "allow_semidefinite")

def __init__(self, N: int):
"""Construct a new chol ext object
Expand Down Expand Up @@ -65,9 +66,9 @@ def factor(self, getA: Callable[[int, int], float]) -> bool:
s = j + 1
d = getA(i, s) - (self._T[i, start:s] @ self._T[start:s, s])
self._T[i, i] = d
if d > 0.:
if d > 0.0:
continue
if d < 0. or not self.allow_semidefinite:
if d < 0.0 or not self.allow_semidefinite:
self.p = start, i + 1
break
start = i + 1 # T[i, i] == 0 (very unlikely), restart at i+1
Expand Down Expand Up @@ -98,7 +99,7 @@ def witness(self):
raise AssertionError()
start, n = self.p
m = n - 1
self.v[m] = 1.
self.v[m] = 1.0
for i in range(m, start, -1):
self.v[i - 1] = -(self._T[i:n, i - 1] @ self.v[i:n])
return -self._T[m, m]
Expand Down
7 changes: 4 additions & 3 deletions src/ellpy/oracles/corr_bspline_oracle.py
Original file line number Diff line number Diff line change
Expand Up @@ -25,8 +25,8 @@ def mono_oracle(x):
for i in range(n - 1):
fj = x[i + 1] - x[i]
if fj > 0:
g[i] = -1.
g[i + 1] = 1.
g[i] = -1.0
g[i + 1] = 1.0
return g, fj


Expand All @@ -36,6 +36,7 @@ class mono_decreasing_oracle2:
Returns:
[type]: [description]
"""

def __init__(self, basis):
"""[summary]
Expand All @@ -61,7 +62,7 @@ def __call__(self, x: Arr, t: float) -> Tuple[Cut, Optional[float]]:
if cut:
g1, fj = cut
g[:-1] = g1
g[-1] = 0.
g[-1] = 0.0
return (g, fj), None
return self.basis(x, t)

Expand Down
4 changes: 2 additions & 2 deletions src/ellpy/oracles/corr_oracle.py
Original file line number Diff line number Diff line change
Expand Up @@ -19,7 +19,7 @@ def create_2d_sites(nx=10, ny=8) -> Arr:
Arr: location of sites
"""
n = nx * ny
s_end = np.array([10., 8.])
s_end = np.array([10.0, 8.0])
hgen = halton([2, 3])
s = s_end * np.array([hgen() for _ in range(n)])
return s
Expand All @@ -39,7 +39,7 @@ def create_2d_isotropic(s: Arr, N=3000) -> Arr:
"""
n = s.shape[0]
sdkern = 0.12 # width of kernel
var = 2. # standard derivation
var = 2.0 # standard derivation
tau = 0.00001 # standard derivation of white noise
np.random.seed(5)

Expand Down
1 change: 1 addition & 0 deletions src/ellpy/oracles/csdlowpass_oracle.py
Original file line number Diff line number Diff line change
Expand Up @@ -18,6 +18,7 @@ class csdlowpass_oracle:
Returns:
[type]: [description]
"""

rcsd = None

def __init__(self, nnz, lowpass):
Expand Down
21 changes: 11 additions & 10 deletions src/ellpy/oracles/cycle_ratio_oracle.py
Original file line number Diff line number Diff line change
Expand Up @@ -12,16 +12,17 @@
class cycle_ratio_oracle:
"""Oracle for minimum ratio cycle problem.
This example solves the following convex problem:
This example solves the following convex problem:
max t
s.t. u[j] - u[i] ≤ mij - sij * x,
t ≤ x
max t
s.t. u[j] - u[i] ≤ mij - sij * x,
t ≤ x
where sij is not necessarily positive.
The problem could be unbounded???
where sij is not necessarily positive.
The problem could be unbounded???
"""

class ratio:
def __init__(self, G):
"""[summary]
Expand All @@ -42,8 +43,8 @@ def eval(self, e, x):
Any: function evaluation
"""
u, v = e
cost = self.G[u][v]['cost']
time = self.G[u][v]['time']
cost = self.G[u][v]["cost"]
time = self.G[u][v]["time"]
return cost - time * x

def grad(self, e, x):
Expand All @@ -57,7 +58,7 @@ def grad(self, e, x):
[type]: [description]
"""
u, v = e
time = self.G[u][v]['time']
time = self.G[u][v]["time"]
return -time

def __init__(self, G, u):
Expand Down Expand Up @@ -89,4 +90,4 @@ def __call__(self, x, t) -> Tuple[Cut, Any]:
if cut:
return cut, None

return (-1, 0.), x
return (-1, 0.0), x
Loading

0 comments on commit ffff6a1

Please sign in to comment.