Transpose instead of PseudoInverse

[unrealwill]

In the blog https://ocramz.github.io/posts/2023-12-21-assignment-riemann-opt.html for the definition of the projection you wrote

α = (I − xx⊤)⊤(z − xz⊤)1

The middle T should be a dagger meaning "left pseudo-inverse" according to the referenced paper [4] https://arxiv.org/abs/1802.02628

In the corresponding code

assignment-riemann-opt/mctorch/manifolds/doublystochastic.py

Line 178 in aed67ed

alpha, beta = self._lsolve(x, b)

you seem to be doing some sort of inversion I don't understand but is probably the intended left pseudo inverse.

assignment-riemann-opt/mctorch/manifolds/doublystochastic.py

Lines 162 to 172 in aed67ed

	def _lsolve(self, x, b):
	# TODO: A better way to solve it is implemented in `_optimLSolve`
	# Once Pytorch gains support for `LinearOperator`/scipy like `cg`
	# function, it can be used.
	alpha = torch.empty((self._k, self._n))
	beta = torch.empty((self._k, self._m))
	for i in range(self._k):
	A = torch.cat((torch.cat((torch.diag(self._p[i]), x[i]), dim=-1), torch.cat((x[i].T, torch.diag(self._q[i])), dim=-1)), dim=0)
	zeta = torch.linalg.solve(A, b[i])
	alpha[i], beta[i] = zeta[:self._n], zeta[self._n:]
	return alpha, beta

Doesn't this inversion computational O(n^3) cost per iteration, render the whole exercise pointless (aka complexity >> Munkres ) ?

[ocramz]

Hi, thank you for looking into this! The code you are pointing at is actually an old version, the current results use the doublystochastic.py module in the root directory. In there you'll see we do a simple transpose, which works because the projector I - x x^T is orthogonal. I should probably clean up this repo a bit! Thank you again, glad to hear the topic is interesting.

…

On Wed, 8 Jan 2025 at 13:30, unrealwill ***@***.***> wrote: In the blog https://ocramz.github.io/posts/2023-12-21-assignment-riemann-opt.html for the definition of the projection you wrote α = (I − xx⊤)⊤(z − xz⊤)1 The middle T should be a dagger meaning "left pseudo-inverse" according to the referenced paper [4] https://arxiv.org/abs/1802.02628 In the corresponding code https://github.com/ocramz/assignment-riemann-opt/blob/aed67ed8be296f4e5850797c841cbcf19338192e/mctorch/manifolds/doublystochastic.py#L178 you seem to be doing some sort of inversion I don't understand but is probably the intended left pseudo inverse. https://github.com/ocramz/assignment-riemann-opt/blob/aed67ed8be296f4e5850797c841cbcf19338192e/mctorch/manifolds/doublystochastic.py#L162-L172 Doesn't this inversion computational O(n^3) cost per iteration, render the whole exercise pointless (aka complexity >> Munkres ) ? — Reply to this email directly, view it on GitHub <#2>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ABNBDKGEKB5EQ22GEY76UZD2JUK6TAVCNFSM6AAAAABUZ2APKGVHI2DSMVQWIX3LMV43ASLTON2WKOZSG43TKMRXHA4TGNI> . You are receiving this because you are subscribed to this thread.Message ID: ***@***.***>

[unrealwill]

Thanks for your quick response, HN front page brought me to your your blog post.

Have you empirically checked that the projection on the tangent is indeed correct when you use transpose instead of dagger (I would find it very surprising even when the matrix is symmetric which I-XXT is, the left pseudo inverse is something different) ?

The tangent space is defined by equation (21), meaning all rows and all columns sum to 0, so it should be quite easy to check whether or not the point is on the tangent space or not.

Coming back to your lsolve of a previous version, I think you were trying to implement the "Remark 2" of the reference paper, but you were doing some splitting along "k" I don't understand.