Backward pass:
dG/dx_i += dG/dy_i dy_i/dx_i = dG/dy_i exp(x_i)
174667
@Override
    public void backward() {
        Tensor tmp = new Tensor(yAdj); // copy
        tmp.elemMultiply(y);
        modInX.getOutputAdj().elemAdd(tmp);
    }
1287249
static void backsubstitution0134(DMatrixRMaj P_plus, DMatrixRMaj P , DMatrixRMaj X ,
									 double H[] ) {
		final int N = P.numRows;
		DMatrixRMaj tmp = new DMatrixRMaj(N*2, 1);

		double H6 = H[6];
		double H7 = H[7];
		double H8 = H[8];

		for (int i = 0, index = 0; i < N; i++) {
			double x = -X.data[index],y = -X.data[index+1];
			double sum = P.data[index++] * H6 + P.data[index++] * H7 + H8;
			tmp.data[i] = x * sum;
			tmp.data[i+N] = y * sum;
		}

		double h0=0,h1=0,h3=0,h4=0;
		for (int i = 0; i < N; i++) {
			double P_pls_0 = P_plus.data[i];
			double P_pls_1 = P_plus.data[i+N];

			double tmp_i = tmp.data[i];
			double tmp_j = tmp.data[i+N];

			h0 += P_pls_0*tmp_i;
			h1 += P_pls_1*tmp_i;
			h3 += P_pls_0*tmp_j;
			h4 += P_pls_1*tmp_j;
		}

		H[0] = -h0;
		H[1] = -h1;
		H[3] = -h3;
		H[4] = -h4;
	}
645002
def adapt(self, d, x):
        
        y = np.dot(self.w, x)
        e = d - y
        self.eps = self.eps - self.ro * self.mu * e * self.last_e * \
            np.dot(x, self.last_x) / \
            (np.dot(self.last_x, self.last_x) + self.eps)**2
        nu = self.mu / (self.eps + np.dot(x, x))
        self.w += nu * e * x        
        self.last_e = e
989821
def policy_backward(eph, epx, epdlogp, model):
    """"""
    dW2 = np.dot(eph.T, epdlogp).ravel()
    dh = np.outer(epdlogp, model["W2"])
    # Backprop relu.
    dh[eph <= 0] = 0
    dW1 = np.dot(dh.T, epx)
    return {"W1": dW1, "W2": dW2}
1087893
def _Backward3_T_Ps(P, s):
    
    sc = 4.41202148223476
    if s <= sc:
        T = _Backward3a_T_Ps(P, s)
    else:
        T = _Backward3b_T_Ps(P, s)
    return T
195172
public static double getGradDotDirApprox(Function fn, IntDoubleVector x, IntDoubleVector d, double c) {
        double dot = 0;
        {
            // L(\theta + c * d)
            IntDoubleVector d1 = d.copy();
            d1.scale(c);
            IntDoubleVector x1 = x.copy();
            x1.add(d1);
            dot += fn.getValue(x1);
        }
        {
            // - L(\theta - c * d)
            IntDoubleVector d1 = d.copy();
            d1.scale(-c);
            IntDoubleVector x1 = x.copy();
            x1.add(d1);        
            dot -= fn.getValue(x1);
        }
        dot /= (2.0 * c);
        return dot;
    }
886371
def xstep(self):
        
        """

        YU = self.Y - self.U
        self.X = np.asarray(sl.cho_solve_ATAI(
            self.D, 1.0, self.block_sep0(YU) +
            self.D.T.dot(self.block_sep1(YU)), self.lu, self.piv),
            dtype=self.dtype)
333466
def _call(self, x):
        """"""
        # The algorithm takes finite differences in one dimension at a time
        # reusing the dx vector to improve efficiency.
        dfdx = self.domain.zero()
        dx = self.domain.zero()

        if self.method == 'backward':
            fx = self.functional(x)
            for i in range(self.domain.size):
                dx[i - 1] = 0  # reset step from last iteration
                dx[i] = self.step
                dfdx[i] = fx - self.functional(x - dx)
        elif self.method == 'forward':
            fx = self.functional(x)
            for i in range(self.domain.size):
                dx[i - 1] = 0  # reset step from last iteration
                dx[i] = self.step
                dfdx[i] = self.functional(x + dx) - fx
        elif self.method == 'central':
            for i in range(self.domain.size):
                dx[i - 1] = 0  # reset step from last iteration
                dx[i] = self.step / 2
                dfdx[i] = self.functional(x + dx) - self.functional(x - dx)
        else:
            raise RuntimeError('unknown method')

        dfdx /= self.step
        return dfdx
444449
def forward(inputs_i, Q, G, A, b, h, U_Q, U_S, R, verbose=False):
    
    nineq, nz, neq, _ = get_sizes(G, A)

    # find initial values
    d = torch.ones(nineq).type_as(Q)
    nb = -b if b is not None else None
    factor_kkt(U_S, R, d)
    x, s, z, y = solve_kkt(
        U_Q, d, G, A, U_S,
        inputs_i, torch.zeros(nineq).type_as(Q), -h, nb)
    # x1, s1, z1, y1 = factor_solve_kkt(Q, torch.eye(nineq).type_as(Q), G, A, inputs_i,
    # torch.zeros(nineq).type_as(Q), -h, nb)

    if torch.min(s) < 0:
        s -= torch.min(s) - 1
    if torch.min(z) < 0:
        z -= torch.min(z) - 1

    prev_resid = None
    for i in range(20):
        # affine scaling direction
        rx = (torch.mv(A.t(), y) if neq > 0 else 0.) + \
            torch.mv(G.t(), z) + torch.mv(Q, x) + inputs_i
        rs = z
        rz = torch.mv(G, x) + s - h
        ry = torch.mv(A, x) - b if neq > 0 else torch.Tensor([0.])
        mu = torch.dot(s, z) / nineq
        pri_resid = torch.norm(ry) + torch.norm(rz)
        dual_resid = torch.norm(rx)
        resid = pri_resid + dual_resid + nineq * mu
        d = z / s
        if verbose:
            print(("primal_res = {0:.5g}, dual_res = {1:.5g}, " +
                   "gap = {2:.5g}, kappa(d) = {3:.5g}").format(
                pri_resid, dual_resid, mu, min(d) / max(d)))
        # if (pri_resid < 5e-4 and dual_resid < 5e-4 and mu < 4e-4):
        improved = (prev_resid is None) or (resid < prev_resid + 1e-6)
        if not improved or resid < 1e-6:
            return x, y, z
        prev_resid = resid

        factor_kkt(U_S, R, d)
        dx_aff, ds_aff, dz_aff, dy_aff = solve_kkt(U_Q, d, G, A, U_S,
                                                   rx, rs, rz, ry)
        # D = torch.diag((z/s).cpu()).type_as(Q)
        # dx_aff1, ds_aff1, dz_aff1, dy_aff1 = factor_solve_kkt(Q, D, G, A, rx, rs, rz, ry)

        # compute centering directions
        alpha = min(min(get_step(z, dz_aff), get_step(s, ds_aff)), 1.0)
        sig = (torch.dot(s + alpha * ds_aff, z +
                         alpha * dz_aff) / (torch.dot(s, z)))**3
        dx_cor, ds_cor, dz_cor, dy_cor = solve_kkt(
            U_Q, d, G, A, U_S, torch.zeros(nz).type_as(Q),
            (-mu * sig * torch.ones(nineq).type_as(Q) + ds_aff * dz_aff) / s,
            torch.zeros(nineq).type_as(Q), torch.zeros(neq).type_as(Q))
        # dx_cor, ds_cor, dz_cor, dy_cor = factor_solve_kkt(Q, D, G, A,
        #     torch.zeros(nz).type_as(Q),
        #     (-mu*sig*torch.ones(nineq).type_as(Q) + ds_aff*dz_aff)/s,
        #     torch.zeros(nineq).type_as(Q), torch.zeros(neq).type_as(Q))

        dx = dx_aff + dx_cor
        ds = ds_aff + ds_cor
        dz = dz_aff + dz_cor
        dy = dy_aff + dy_cor if neq > 0 else None
        alpha = min(1.0, 0.999 * min(get_step(s, ds), get_step(z, dz)))
        dx_norm = torch.norm(dx)
        dz_norm = torch.norm(dz)
        if np.isnan(dx_norm) or dx_norm > 1e5 or dz_norm > 1e5:
            # Overflow, return early
            return x, y, z

        x += alpha * dx
        s += alpha * ds
        z += alpha * dz
        y = y + alpha * dy if neq > 0 else None

    return x, y, z
1153747
def Itaax(mt, x, t):
    
    return (Sx(mt, x) - Sx(mt, x + t)) / Dx(mt, x)
165240
@Override
    public void backward() {
        Tensor tmp1 = new Tensor(yAdj); // copy
        tmp1.multiply(weightX);
        modInX.getOutputAdj().elemAdd(tmp1);

        Tensor tmp2 = new Tensor(yAdj); // copy
        tmp2.multiply(weightW);
        modInW.getOutputAdj().elemAdd(tmp2);
    }
528003
def forward(self, j=1):
        
        if j < 0:
            return self.backward(-j)
        self.__i += j
        return self[self.__i-j:self.__i]
468597
def _call(self, dx):
        """"""
        x = self.point

        dx_norm = dx.norm()
        if dx_norm == 0:
            return 0

        scaled_dx = dx * (self.step / dx_norm)

        if self.method == 'backward':
            dAdx = self.operator(x) - self.operator(x - scaled_dx)
        elif self.method == 'forward':
            dAdx = self.operator(x + scaled_dx) - self.operator(x)
        elif self.method == 'central':
            dAdx = (self.operator(x + scaled_dx / 2) -
                    self.operator(x - scaled_dx / 2))
        else:
            raise RuntimeError('unknown method')

        return dAdx * (dx_norm / self.step)
948225
def DebyeChain(q, Rg):
    
    a = (q * Rg) ** 2
    return 2 * (np.exp(-a) - 1 + a) / a ** 2
1166329
def chain_2(d2f_dg2, dg_dx, df_dg, d2g_dx2):
    
    if np.all(dg_dx==1.) and np.all(d2g_dx2 == 0):
        return d2f_dg2
    dg_dx_2 = np.clip(dg_dx, -np.inf, _lim_val_square)**2
    #dg_dx_2 = dg_dx**2
    return d2f_dg2*(dg_dx_2) + df_dg*d2g_dx2
366856
def _deriv_logaddexp2(x1, x2):
    """"""
    y1 = np.exp2(x1)
    y2 = np.exp2(x2)
    df_dx1 = y1 / (y1 + y2)
    df_dx2 = y2 / (y1 + y2)
    return np.vstack([df_dx1, df_dx2]).T
1247906
def xstep(self):
        
        """

        ngsit = 0
        gsrrs = np.inf
        while gsrrs > self.opt['GSTol'] and ngsit < self.opt['MaxGSIter']:
            self.X = self.GaussSeidelStep(self.S, self.X,
                                          self.cnst_AT(self.Y-self.U),
                                          self.rho, self.lcw, self.Wdf2)
            gsrrs = sl.rrs(
                self.rho*self.cnst_AT(self.cnst_A(self.X)) +
                self.Wdf2*self.X, self.Wdf2*self.S +
                self.rho*self.cnst_AT(self.Y - self.U))
            ngsit += 1

        self.xs = (ngsit, gsrrs)