view release on metacpan or  search on metacpan

lbfgs.c  view on Meta::CPAN

		 */
		vecdot(&ys, it->y, it->s, n);
		vecdot(&yy, it->y, it->y, n);
		it->ys = ys;

		/*
			Recursive formula to compute dir = -(H \cdot g).
				This is described in page 779 of:
				Jorge Nocedal.
				Updating Quasi-Newton Matrices with Limited Storage.
				Mathematics of Computation, Vol. 35, No. 151,
				pp. 773--782, 1980.
		 */
		bound = (m <= k) ? m : k;
		++k;
		end = (end + 1) % m;

		if (param->orthantwise_c == 0.) {
			/* Compute the negative of gradients. */
			vecncpy(d, g, n);
		} else {
			/* Compute the negative of psuedo-gradients. */
			for (i = 0;i < n;++i) {
				if (x[i] < 0.) {
					/* Differentiable. */
					d[i] = -g[i] + param->orthantwise_c;
				} else if (0. < x[i]) {
					/* Differentiable. */
					d[i] = -g[i] - param->orthantwise_c;
				} else {
					if (g[i] < -param->orthantwise_c) {
						/* Take the right partial derivative. */
						d[i] = -g[i] - param->orthantwise_c;
					} else if (param->orthantwise_c < g[i]) {
						/* Take the left partial derivative. */
						d[i] = -g[i] + param->orthantwise_c;
					} else {
						d[i] = 0.;
					}
				}
			}
			/* Store the steepest direction.*/
			veccpy(w, d, n);
		}

		j = end;
		for (i = 0;i < bound;++i) {
			j = (j + m - 1) % m;	/* if (--j == -1) j = m-1; */
			it = &lm[j];
			/* \alpha_{j} = \rho_{j} s^{t}_{j} \cdot q_{k+1}. */
			vecdot(&it->alpha, it->s, d, n);
			it->alpha /= it->ys;
			/* q_{i} = q_{i+1} - \alpha_{i} y_{i}. */
			vecadd(d, it->y, -it->alpha, n);
		}

		vecscale(d, ys / yy, n);

		for (i = 0;i < bound;++i) {
			it = &lm[j];
			/* \beta_{j} = \rho_{j} y^t_{j} \cdot \gamma_{i}. */
			vecdot(&beta, it->y, d, n);
			beta /= it->ys;
			/* \gamma_{i+1} = \gamma_{i} + (\alpha_{j} - \beta_{j}) s_{j}. */
			vecadd(d, it->s, it->alpha - beta, n);
			j = (j + 1) % m;		/* if (++j == m) j = 0; */
		}

		/*
			Constrain the search direction for orthant-wise updates.
		 */
		if (param->orthantwise_c != 0.) {
			for (i = 0;i < n;++i) {
				if (d[i] * w[i] <= 0) {
					d[i] = 0;
				}
			}
		}

		/*
			Now the search direction d is ready. We try step = 1 first.
		 */
		step = 1.0;
	}

lbfgs_exit:
	/* Return the final value of the objective function. */
	if (ptr_fx != NULL) {
		*ptr_fx = fx;
	}

	/* Free memory blocks used by this function. */
	if (lm != NULL) {
		for (i = 0;i < m;++i) {
			vecfree(lm[i].s);
			vecfree(lm[i].y);
		}
		vecfree(lm);
	}
	vecfree(w);
	vecfree(d);
	vecfree(gp);
	vecfree(g);
	vecfree(xp);

	return ret;
}



static int line_search_backtracking(
	int n,
	lbfgsfloatval_t *x,
	lbfgsfloatval_t *f,
	lbfgsfloatval_t *g,
	lbfgsfloatval_t *s,
	lbfgsfloatval_t *stp,
	lbfgsfloatval_t *xp,
	lbfgs_evaluate_t proc_evaluate,
	void *instance,
	const lbfgs_parameter_t *param
	)
{
	int i, ret = 0, count = 0;