glimix_core/cov/_kron2sum.py
def LdKL_dot(self, v, v1=None):
        """
        Implements L(∂K)Lᵀv.

        The array v can have one or two dimensions and the first dimension has to have
        size n⋅p.

        Let vec(V) = v. We have

            L(∂K)Lᵀ⋅v = ((Lₕ∂C₀Lₕᵀ) ⊗ (LₓGGᵀLₓᵀ))vec(V) = vec(LₓGGᵀLₓᵀVLₕ∂C₀Lₕᵀ),

        when the derivative is over the parameters of C₀. Similarly,

            L(∂K)Lᵀv = ((Lₕ∂C₁Lₕᵀ) ⊗ (LₓLₓᵀ))vec(V) = vec(LₓLₓᵀVLₕ∂C₁Lₕᵀ),

        over the parameters of C₁.
        """
        self._init_svd()

        def dot(a, b):
            r = tensordot(a, b, axes=([1], [0]))
            if a.ndim > b.ndim:
                return r.transpose([0, 2, 1])
            return r

        Lh = self.Lh
        V = unvec(v, (self.Lx.shape[0], -1) + v.shape[1:])
        LdKL_dot = {
            "C0.Lu": empty((v.shape[0],) + v.shape[1:] + (self._C0.Lu.shape[0],)),
            "C1.Lu": empty((v.shape[0],) + v.shape[1:] + (self._C1.Lu.shape[0],)),
        }

        dC0 = self._C0.gradient()["Lu"]
        for i in range(self._C0.Lu.shape[0]):
            t = dot(self._LxG, dot(self._LxG.T, dot(V, Lh @ dC0[..., i] @ Lh.T)))
            LdKL_dot["C0.Lu"][..., i] = t.reshape((-1,) + t.shape[2:], order="F")

        dC1 = self._C1.gradient()["Lu"]
        for i in range(self._C1.Lu.shape[0]):
            t = dot(V, Lh @ dC1[..., i] @ Lh.T)
            LdKL_dot["C1.Lu"][..., i] = t.reshape((-1,) + t.shape[2:], order="F")

        return LdKL_dot
def LdKL_dot(self, v, v1=None):
        """
        Implements L(∂K)Lᵀv.

        The array v can have one or two dimensions and the first dimension has to have
        size n⋅p.

        Let vec(V) = v. We have

            L(∂K)Lᵀ⋅v = ((Lₕ∂C₀Lₕᵀ) ⊗ (LₓGGᵀLₓᵀ))vec(V) = vec(LₓGGᵀLₓᵀVLₕ∂C₀Lₕᵀ),

        when the derivative is over the parameters of C₀. Similarly,

            L(∂K)Lᵀv = ((Lₕ∂C₁Lₕᵀ) ⊗ (LₓLₓᵀ))vec(V) = vec(LₓLₓᵀVLₕ∂C₁Lₕᵀ),

        over the parameters of C₁.
        """
        self._init_svd()

        def dot(a, b):
            r = tensordot(a, b, axes=([1], [0]))
            if a.ndim > b.ndim:
                return r.transpose([0, 2, 1])
            return r

        Lh = self.Lh
        V = unvec(v, (self.Lx.shape[0], -1) + v.shape[1:])
        LdKL_dot = {
            "C0.Lu": empty((v.shape[0],) + v.shape[1:] + (self._C0.Lu.shape[0],)),
            "C1.Lu": empty((v.shape[0],) + v.shape[1:] + (self._C1.Lu.shape[0],)),
        }

        dC0 = self._C0.gradient()["Lu"]
        for i in range(self._C0.Lu.shape[0]):
            t = dot(self._LxG, dot(self._LxG.T, dot(V, Lh @ dC0[..., i] @ Lh.T)))
            LdKL_dot["C0.Lu"][..., i] = t.reshape((-1,) + t.shape[2:], order="F")

        dC1 = self._C1.gradient()["Lu"]
        for i in range(self._C1.Lu.shape[0]):
            t = dot(V, Lh @ dC1[..., i] @ Lh.T)
            LdKL_dot["C1.Lu"][..., i] = t.reshape((-1,) + t.shape[2:], order="F")

        return LdKL_dot
def LdKL_dot ( self , v , v1 = None ) : self . _init_svd ( ) def dot ( a , b ) : r = tensordot ( a , b , axes = ( [ 1 ] , [ 0 ] ) ) if a . ndim > b . ndim : return r . transpose ( [ 0 , 2 , 1 ] ) return r Lh = self . Lh V = unvec ( v , ( self . Lx . shape [ 0 ] , - 1 ) + v . shape [ 1 : ] ) LdKL_dot = { "C0.Lu" : empty ( ( v . shape [ 0 ] , ) + v . shape [ 1 : ] + ( self . _C0 . Lu . shape [ 0 ] , ) ) , "C1.Lu" : empty ( ( v . shape [ 0 ] , ) + v . shape [ 1 : ] + ( self . _C1 . Lu . shape [ 0 ] , ) ) , } dC0 = self . _C0 . gradient ( ) [ "Lu" ] for i in range ( self . _C0 . Lu . shape [ 0 ] ) : t = dot ( self . _LxG , dot ( self . _LxG . T , dot ( V , Lh @ dC0 [ ... , i ] @ Lh . T ) ) ) LdKL_dot [ "C0.Lu" ] [ ... , i ] = t . reshape ( ( - 1 , ) + t . shape [ 2 : ] , order = "F" ) dC1 = self . _C1 . gradient ( ) [ "Lu" ] for i in range ( self . _C1 . Lu . shape [ 0 ] ) : t = dot ( V , Lh @ dC1 [ ... , i ] @ Lh . T ) LdKL_dot [ "C1.Lu" ] [ ... , i ] = t . reshape ( ( - 1 , ) + t . shape [ 2 : ] , order = "F" ) return LdKL_dot
Implements L(∂K)Lᵀv.

        The array v can have one or two dimensions and the first dimension has to have
        size n⋅p.

        Let vec(V) = v. We have

            L(∂K)Lᵀ⋅v = ((Lₕ∂C₀Lₕᵀ) ⊗ (LₓGGᵀLₓᵀ))vec(V) = vec(LₓGGᵀLₓᵀVLₕ∂C₀Lₕᵀ),

        when the derivative is over the parameters of C₀. Similarly,

            L(∂K)Lᵀv = ((Lₕ∂C₁Lₕᵀ) ⊗ (LₓLₓᵀ))vec(V) = vec(LₓLₓᵀVLₕ∂C₁Lₕᵀ),

        over the parameters of C₁.
Implements L ( ∂K ) Lᵀv .