Skip to content
This repository was archived by the owner on Dec 6, 2023. It is now read-only.

Remove naked mallocs and use cython.cvarray +refcounting #5

Closed
wants to merge 2 commits into from
Closed

Remove naked mallocs and use cython.cvarray +refcounting #5

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
2150149
Remove naked mallocs and use numpy and refcounting
vene Aug 15, 2016
c5caa99
Use cvarray instead of np.empty
vene Aug 15, 2016
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
900 changes: 470 additions & 430 deletions polylearn/cd_direct_fast.cpp
View file
Open in desktop

Large diffs are not rendered by default.

62 changes: 27 additions & 35 deletions polylearn/cd_direct_fast.pyx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -7,8 +7,8 @@
# Author: Vlad Niculae
# License: BSD

from libc.stdlib cimport malloc, free
from libc.math cimport fabs
from cython.view cimport array

from lightning.impl.dataset_fast cimport ColumnDataset

Expand All @@ -19,7 +19,7 @@ from .cd_linear_fast cimport _cd_linear_epoch
cdef void _precompute(ColumnDataset X,
double[:, :, ::1] P,
Py_ssize_t order,
double* out,
double[:, ::1] out,
Py_ssize_t s,
unsigned int degree):

Expand All @@ -32,15 +32,17 @@ cdef void _precompute(ColumnDataset X,
cdef int n_nz

cdef Py_ssize_t i, j, ii
cdef unsigned int d
cdef double tmp

for i in range(n_samples):
out[i] = 0
out[degree - 1, i] = 0

for j in range(n_features):
X.get_column_ptr(j, &indices, &data, &n_nz)
for ii in range(n_nz):
i = indices[ii]
out[i] += (data[ii] * P[order, s, j]) ** degree
out[degree - 1, i] += (data[ii] * P[order, s, j]) ** degree


cdef inline double _update(int* indices,
Expand All @@ -50,12 +52,11 @@ cdef inline double _update(int* indices,
double[:] y,
double[:] y_pred,
LossFunction loss,
double* d1,
double* d2,
unsigned int degree,
double lam,
double beta,
double* cache_kp):
double[:, ::1] D,
double[:] cache_kp):

cdef double l1_reg = 2 * beta * fabs(lam)

Expand All @@ -70,10 +71,10 @@ cdef inline double _update(int* indices,
i = indices[ii]

if degree == 2:
kp = d1[i] - p_js * data[ii]
elif degree == 3:
kp = 0.5 * (d1[i] ** 2 - d2[i])
kp -= p_js * data[ii] * d1[i]
kp = D[0, i] - p_js * data[ii]
else: # degree == 3:
kp = 0.5 * (D[0, i] ** 2 - D[1, i])
kp -= p_js * data[ii] * D[0, i]
kp += p_js ** 2 * data[ii] ** 2

kp *= lam * data[ii]
Expand All @@ -97,12 +98,11 @@ cdef inline double _cd_direct_epoch(double[:, :, ::1] P,
double[:] y,
double[:] y_pred,
double[:] lams,
double* d1,
double* d2,
unsigned int degree,
double beta,
LossFunction loss,
double* cache_kp):
double[:, ::1] D,
double[:] cache_kp):

cdef Py_ssize_t s, j
cdef double p_old, update, offset
Expand All @@ -118,9 +118,9 @@ cdef inline double _cd_direct_epoch(double[:, :, ::1] P,
for s in range(n_components):

# initialize the cached ds for this s
_precompute(X, P, order, d1, s, 1)
_precompute(X, P, order, D, s, 1)
if degree == 3:
_precompute(X, P, order, d2, s, 2)
_precompute(X, P, order, D, s, 2)

for j in range(n_features):

Expand All @@ -129,7 +129,7 @@ cdef inline double _cd_direct_epoch(double[:, :, ::1] P,
# compute coordinate update
p_old = P[order, s, j]
update = _update(indices, data, n_nz, p_old, y, y_pred,
loss, d1, d2, degree, lams[s], beta, cache_kp)
loss, degree, lams[s], beta, D, cache_kp)
P[order, s, j] -= update
sum_viol += fabs(update)

Expand All @@ -138,9 +138,10 @@ cdef inline double _cd_direct_epoch(double[:, :, ::1] P,
i = indices[ii]

if degree == 3:
d2[i] -= (p_old ** 2 - P[order, s, j] ** 2) * data[ii] ** 2
D[1, i] -= ((p_old ** 2 - P[order, s, j] ** 2) *
data[ii] ** 2)

d1[i] -= update * data[ii]
D[0, i] -= update * data[ii]
y_pred[i] -= update * cache_kp[ii]
return sum_viol

Expand Down Expand Up @@ -169,11 +170,8 @@ def _cd_direct_ho(double[:, :, ::1] P not None,
cdef bint converged = False

# precomputed values
cdef double *d1 = <double *> malloc(n_samples * sizeof(double))
cdef double *d2
if degree == 3:
d2 = <double *> malloc(n_samples * sizeof(double))
cdef double *cache_kp = <double *> malloc(n_samples * sizeof(double))
cdef double[:, ::1] D = array((degree - 1, n_samples), sizeof(double), 'd')
cdef double[:] cache_kp = array((n_samples,), sizeof(double), 'd')

for it in range(max_iter):
viol = 0
Expand All @@ -182,11 +180,11 @@ def _cd_direct_ho(double[:, :, ::1] P not None,
viol += _cd_linear_epoch(w, X, y, y_pred, col_norm_sq, alpha, loss)

if fit_lower and degree == 3: # fit degree 2. Will be looped later.
viol += _cd_direct_epoch(P, 1, X, y, y_pred, lams, d1, d2,
2, beta, loss, cache_kp)
viol += _cd_direct_epoch(P, 1, X, y, y_pred, lams, 2, beta, loss,
D, cache_kp)

viol += _cd_direct_epoch(P, 0, X, y, y_pred, lams, d1, d2,
degree, beta, loss, cache_kp)
viol += _cd_direct_epoch(P, 0, X, y, y_pred, lams, degree, beta, loss,
D, cache_kp)

if verbose:
print("Iteration", it + 1, "violation sum", viol)
Expand All @@ -197,10 +195,4 @@ def _cd_direct_ho(double[:, :, ::1] P not None,
converged = True
break

# Free up cache
free(d1)
free(cache_kp)
if degree == 3:
free(d2)

return converged
Loading