User federico poloni - MathOverflow most recent 30 from http://mathoverflow.net 2013-05-18T07:51:26Z http://mathoverflow.net/feeds/user/1898 http://www.creativecommons.org/licenses/by-nc/2.5/rdf http://mathoverflow.net/questions/129741/spectral-radius-monotonicity/129789#129789 Answer by Federico Poloni for spectral radius monotonicity Federico Poloni 2013-05-05T23:10:55Z 2013-05-15T07:17:51Z <p>Not true in general, as noted by @SergeiIvanov, but true for (element-wise) nonnegative matrices.</p> <p>Note that if $\rho(S) &lt; b$, then $b(bI-S)^{-1}=(I-\frac{S}{b})^{-1}=\sum_{i=0}^\infty \frac{S^i}{b^i}$. In particular, thanks to this expansion, if $\rho(S)&lt; a &lt; b$, then $b(bI-S)^{-1}&lt; a(aI-S)^{-1}$ in the componentwise ordering, and thus also $b(bI-S)^{-1}T \leq a(aI-S)^{-1}T$ for any nonnegative $T$. Now, it is a part of the Perron-Frobenius theorem that for any $A,B$ with $0 \leq A \leq B$ then $\rho(A) \leq \rho(B)$, and that's all we need here.</p> http://mathoverflow.net/questions/128285/the-first-eigenvalue-of-a-branching-process-matrix/128320#128320 Answer by Federico Poloni for The first eigenvalue of a branching process matrix Federico Poloni 2013-04-22T07:31:06Z 2013-04-22T07:31:06Z <p>You can see this process as a dynamical system or a Markov chain without normalization. If the matrix is irreducible, starting from every initial distribution of number of individuals $w_0$, the process will "converge" (in some suitable sense) to $w_{k}=\alpha \lambda^k v$, for some $\alpha\in\mathbb{R}$, and $(\lambda,v)$ the Perron eigenpair.</p> <p>Thus, in the stationary limit, the ratios among the number of individuals of different types at each time step $k$ $(w_k)_i/(w_k)_j$ are the ratios of components of the Perron vector $v_i/v_j$, while the number of individuals is multiplied by $\lambda$ at each iteration. So $\lambda$ is a growth factor for the number of individuals at each iteration.</p> http://mathoverflow.net/questions/127204/efficient-computation-of-markov-chain-transition-probability-matrix/127626#127626 Answer by Federico Poloni for Efficient computation of Markov chain transition probability matrix Federico Poloni 2013-04-15T14:35:07Z 2013-04-15T14:35:07Z <p>What you need is called "computing the action of the matrix exponential" (that is, computing $\exp(A)b$ without forming $\exp(A)$ explicitly. There are techniques based on complex integrals and Krylov subspaces. See <a href="http://eprints.ma.man.ac.uk/1426/" rel="nofollow">http://eprints.ma.man.ac.uk/1426/</a> and the references included there.</p> http://mathoverflow.net/questions/117415/old-books-still-used/117737#117737 Answer by Federico Poloni for Old books still used Federico Poloni 2012-12-31T15:39:39Z 2013-04-14T21:21:02Z <p>In numerical linear algebra, Gantmacher's <em>The theory of matrices</em> is still a widely read and cited text (see <a href="http://www.ams.org/mathscinet/search/publications.html?refcit=1657129&amp;loc=refcit" rel="nofollow">MathSciNet</a> <a href="http://www.ams.org/mathscinet/search/publications.html?refcit=107647%2520107648%2520107649&amp;loc=refcit" rel="nofollow">citations</a>). The Russian original dates back to 1953 (thanks @Giuseppe), and the first English translation is from 1959.</p> http://mathoverflow.net/questions/126164/strassens-algorithm/126208#126208 Answer by Federico Poloni for Strassen's algorithm Federico Poloni 2013-04-01T20:51:07Z 2013-04-01T20:51:07Z <p>Quid's reply already gives you a complete answer for the $2\times 2$ case; however, there is an additional consideration that adds to the picture of the optimal exponent of matrix multiplication. </p> <p>One can prove that, for any $t$, if there exists a straight-line program that multiplies $n\times n$ matrices in $O(n^t)$ ops (without any bilinearity assumption or restriction on the operations made: divisions, additions, products of three more terms, everything is allowed), then there is a bilinear algorithm (i.e., one that can be "formulated with tensors") that achieves the same asymptotic complexity $O(n^t)$. So one does not lose generality by restricting to bilinear algorithms, when looking for the optimal exponent of matrix multiplication.</p> <p>I had this same doubt myself a few years ago, and I found a proof of this statement in Borodin-Munro, <em>The computational complexity of algebraic and numeric problems</em>.</p> http://mathoverflow.net/questions/122818/nth-root-of-a-matrix-as-an-analytic-function/122823#122823 Answer by Federico Poloni for Nth root of a matrix as an analytic function? Federico Poloni 2013-02-24T20:08:03Z 2013-02-24T20:08:03Z <p>If I am reading this correctly, you are fine with a power series whose (scalar) coefficients depend on the matrix $A$. In this case, it suffices to take a polynomial $p$ that interpolates $\sqrt[n]{x}$, such that for each eigenvalue $\lambda$ with multiplicity $k_\lambda$, the first $k_\lambda-1$ derivatives of $p$ coincide with those of $\sqrt[n]{x}$ (Hermite interpolant). A degree-$k$ polynomial will always do the job.</p> http://mathoverflow.net/questions/101644/fiction-books-about-mathematicians/101688#101688 Answer by Federico Poloni for Fiction books about mathematicians? Federico Poloni 2012-07-08T17:01:25Z 2013-02-19T09:33:08Z <p><a href="https://en.wikipedia.org/wiki/Greg_Egan" rel="nofollow">Greg Egan</a> is a science-fiction writer that holds a B.S. in Mathematics (and has co-authored a paper with John Baez). He often manages to insert some advanced maths, physics and computer science content in his novels: for instance, listing only mathematics, fiber bundles in <em>Diaspora</em>, Einstein's equation for general relativity in <em>Incandescence</em>, Cantor sets and commutative hypercubes in the short stories <em>The Infinite Assassin</em> and <a href="http://outofthiseos.typepad.com/blog/files/GregEganGlory.pdf" rel="nofollow"><em>Glory</em></a>. </p> <p>His story <a href="http://www.asimovs.com/_issue_0805/DarkINtegers.shtml" rel="nofollow"><em>Dark Integers</em></a> deserves special mention; it is a sequel to <em>Luminous</em>, best read in order.</p> <p>It is truly science fiction written for scientists and mathematicians in particular; they are the only readers that are able to grasp fully both the casual references to advanced mathematical content and the grand ideas underlying his stories. Even after a master in pure maths and a phd in numerical analysis, often I feel that I do not know enough geometry and theoretical physics to get all the facets and implications of what he writes. </p> <p>This feature sets him apart from most other writers in this list, who address maths from a popular-science point of view.</p> <p>On the top of my head, I find it difficult to name a novel of his that does <em>not</em> feature a scientist among the protagonists.</p> http://mathoverflow.net/questions/101644/fiction-books-about-mathematicians/122281#122281 Answer by Federico Poloni for Fiction books about mathematicians? Federico Poloni 2013-02-19T09:20:44Z 2013-02-19T09:20:44Z <p>I was a bit disappointed by Iain M. Banks' <a href="https://en.wikipedia.org/wiki/The_Algebraist" rel="nofollow"><em>The Algebraist</em></a>. Apart from the nice title, which is what made me buy the book, there is very little mathematical content in it.</p> <p>(Well, I guess I am spoiled by Greg Egan.)</p> http://mathoverflow.net/questions/122186/sum-of-elements-of-inverse-matrix/122188#122188 Answer by Federico Poloni for Sum of elements of inverse matrix Federico Poloni 2013-02-18T16:21:04Z 2013-02-18T16:21:04Z <p>The sum of the element of a matrix $M$ is $e^T M e$, where $e$ is the vector of all ones.</p> <p>So, instead of computing the inverse, you should solve the system $Ax=e$ and then compute $e^Tx$. This might look like a simple trick, but solving linear systems is faster than computing inverses in basically all settings.</p> <p>Of course you should then use a method to solve this linear system which is appropriate to the matrix that you are dealing with (but there is a large amount of literature on that).</p> <p>I don't think that you can get the quantity you want any faster than this, unless your matrix has very special properties.</p> http://mathoverflow.net/questions/120545/compute-roots-of-sum-i-c-i-a-i-b-i-xp/120547#120547 Answer by Federico Poloni for Compute roots of sum_i c_i/(a_i + b_i x)^p Federico Poloni 2013-02-01T20:56:48Z 2013-02-01T20:56:48Z <p>This recent <a href="http://git.phc.unipi.it/tesi-equazioni-secolari-e-polinomi/tesi-equazioni-secolari-e-polinomi" rel="nofollow">master thesis</a> by Leonardo Robol treats the case $p=1$ in a numerically sound way. I think they are going to release some code soon, so you might want to contact the author.</p> http://mathoverflow.net/questions/119765/kronecker-structured-matrix-kernel Kronecker-structured matrix kernel Federico Poloni 2013-01-24T15:20:38Z 2013-01-25T08:51:57Z <p>Let $A,B\in\mathbb{C}^{n\times 3n}$ be two matrices, and denote the Kronecker matrix product by $\otimes$. The matrix $$M= \begin{bmatrix} A \otimes I_n \\ I_n \otimes B\end{bmatrix}$$ has size $2n^2 \times 3n^2$, therefore in the generic case its kernel has dimension $n^2$.</p> <blockquote> <p>Is there a simple characterization of this kernel, or a way to compute it by exploiting this structure? </p> </blockquote> <p>This problem comes from a numerical linear algebra problem (two-parameter eigenvalues problems). Do matrices with this structure appear in other fields, or do they ring a bell to you?</p> <p>It looks appealing to look for vectors in this kernel in the form $v=a \otimes b \otimes c$, with $a,c\in\mathbb{C}^{n}$, $b\in\mathbb{C}^{3}$, since in this way they would be compatible with both Kronecker product structures. We normalize by restricting one entry of each of $a,b,c$ to be $1$, so we are left with $(n-1) + 2 + (n-1) = 2n$ unknowns. The condition $Mv=0$ is equivalent to $M(a\otimes b)=0$ and $N(b\otimes c)=0$, which are precisely $2n$ equations. So it looks like we should have several complex solutions in the generic case. Is there a viable way to compute them?</p> http://mathoverflow.net/questions/117963/linear-maps-preserving-positive-semidefiniteness Linear maps preserving positive semidefiniteness Federico Poloni 2013-01-03T14:35:34Z 2013-01-03T15:43:03Z <p>I know of Choi's theorem and some related problems, but not a solution to this exact problem:</p> <blockquote> <p>Characterize the linear maps from the space $S_n$ of symmetric $n \times n$ matrices to itself that preserve positive semidefiniteness.</p> </blockquote> <p>It looks a natural question; has a simple characterization been found? Where can I find it?</p> http://mathoverflow.net/questions/117635/strictly-diagonally-dominant-hermitian-matrices-eigenvalues-sign/117639#117639 Answer by Federico Poloni for Strictly diagonally dominant hermitian matrices eigenvalues sign Federico Poloni 2012-12-30T15:10:35Z 2012-12-30T15:10:35Z <p>You do not need to show that the discs are disjoints. In fact, this won't hold for most diagonally dominant matrices, unlike the main result that you wish to prove.</p> <p>What you need is a stronger form of the Gerschgorin disc thorem, which is due to O. Taussky-Todd and is today normally taught alongside the standard version:</p> <p><em><strong>Theorem</em></strong> If the union of $k$ Gerschgorin discs is disjoint from the union of the other $n-k$ discs then the former union contains exactly $k$ and the latter $n-k$ eigenvalues of $A$. (<a href="https://en.wikipedia.org/wiki/Gershgorin_circle_theorem#Strengthening_of_the_theorem" rel="nofollow">Wikipedia reference with proof</a>).</p> http://mathoverflow.net/questions/117302/canonical-form-of-a-general-bilinear-form/117316#117316 Answer by Federico Poloni for Canonical form of a general Bilinear Form Federico Poloni 2012-12-27T10:45:23Z 2012-12-27T10:45:23Z <p>People in matrix analysis would call this a "canonical form under congruence". Take a look at <a href="http://arxiv.org/abs/0709.2473" rel="nofollow">http://arxiv.org/abs/0709.2473</a>; the solution is stated there.</p> http://mathoverflow.net/questions/117125/a-novice-question-on-quantum-mechanics/117133#117133 Answer by Federico Poloni for A novice question on Quantum Mechanics Federico Poloni 2012-12-24T09:21:58Z 2012-12-24T14:52:16Z <p>As a linear algebraist, the following interpretation is the one that helped me the most: the set of (non-necessarily pure) quantum states is in a one-to-one correspondence with trace-1 symmetric positive-semidefinite matrices; pure states are the ones that correspond to rank-1 matrices $\mid \! u \rangle \langle u \! \mid$, while non-pure states correspond to their convex combinations.</p> <p>So instead of thinking about vectors and their "formal convex combinations", you just have to think about symmetric rank-1 matrices and their convex combinations in the usual sense. This interpretation will also come in useful later on in your study of QM.</p> <p>[EDIT: fixed two mistakes pointed out in the comments]</p> http://mathoverflow.net/questions/110872/condition-number-related-to-root-finding-problems/110895#110895 Answer by Federico Poloni for Condition Number related to Root finding problems Federico Poloni 2012-10-28T11:11:44Z 2012-10-28T11:11:44Z <p>In addition to the convergence speed (and radius) mentioned in Pietro Majer's answer, there is another factor: if the problem is ill-conditioned, the solution is sensitive to perturbations. </p> <p>If you make an error of magnitude $\varepsilon$ in computing the iteration or the parameters data, then the solution is perturbed by a quantity $\kappa \varepsilon$, where $\kappa$ is the condition number. In practice, you have an error of at least one part in $10^{-16}$ in most computer implementations, and often an even larger one if your function depends on measured data or approximations. So Newton's method may converge without trouble, but the solution that you get could be total garbage.</p> http://mathoverflow.net/questions/110816/nonlinear-matrix-equation-2/110818#110818 Answer by Federico Poloni for Nonlinear matrix equation 2 Federico Poloni 2012-10-27T09:37:03Z 2012-10-27T09:37:03Z <p>Only a partial answer:</p> <p>Suppose $\alpha:=v^TAw\neq 0$, $\beta:=v^T B w \neq 0$. Then the equations are $$(\alpha A + \beta B)v=\lambda_1v+\lambda_2w,$$ $$(\alpha A + \beta B)w=\lambda_2v+\lambda_1w.$$ The idea in the answer to your <a href="http://mathoverflow.net/questions/109792/nonlinear-matrix-equation" rel="nofollow">other question</a> still applies, you can tell that $v$ and $w$ are linear combinations of two eigenvectors $x_1,x_2$ of $\alpha A +\beta B$. I am not sure that this is enough in this case though, since getting the eigenpairs of $\alpha A + \beta B$ for yet-unknown $\alpha,\beta$ is not as easy as for a known matrix. </p> http://mathoverflow.net/questions/110706/fixed-points-of-system-of-quadratic-equations/110725#110725 Answer by Federico Poloni for fixed points of system of quadratic equations Federico Poloni 2012-10-26T06:16:25Z 2012-10-26T06:16:25Z <p>Are you sure that this holds in general? I have studied a very similar problem a few years ago (<a href="https://www.sciencedirect.com/science/article/pii/S0024379511004484" rel="nofollow">paper</a>, <a href="http://arxiv.org/abs/1004.1500" rel="nofollow">arxiv</a>), and there you can have up to $2^n$ solutions and you need a minimality property to characterize the physically meaningful solution. Maybe you can get something useful out of the paper.</p> http://mathoverflow.net/questions/110189/number-of-parameters-needed-to-specify-a-hermitian-matrix-of-rank-r/110194#110194 Answer by Federico Poloni for Number of parameters needed to specify a Hermitian matrix of rank r. Federico Poloni 2012-10-20T21:17:22Z 2012-10-20T21:17:22Z <p>Not sure if I am missing something here...</p> <p>1) Rank-$r$ Hermitian matrices are determined uniquely by their image $U$ and how they act when restricted to $U$. The image can be any dimension-$r$ subspace. Almost all subspaces have a basis in the form $\begin{bmatrix}I_r \\ X \end{bmatrix}$, with $X$ any $(n-r)\times r$ matrix, so you have $2(n-r)r$ degrees of freedom for the image. Possible actions on $U$ are isomorphic to $r\times r$ Hermitian matrices, so $r(r-1)+r=r^2$ real dof's. Overall this makes $2nr-r^2$ parameters.</p> <p>2) If you know that $H$ is spd, then you have to restrict the second part to positive-definite matrices, but they still have the same number of parameters, so you still get the same answer.</p> http://mathoverflow.net/questions/109949/how-to-efficiently-compute-the-generalized-cross-product/109995#109995 Answer by Federico Poloni for How to efficiently compute the generalized cross product? Federico Poloni 2012-10-18T09:29:01Z 2012-10-18T09:29:01Z <p>Regarding your question 2: the approach I'd take is computing the first determinant from a RQ factorization of the leading $(n-1)\times(n-1)$ matrix, and then each other by replacing in turn each row of $R$ with $(last row)Q$ and re-orthogonalizing manually with $O(n)$ Givens transformations on the left. In this way you pay $O(n^3)$ for the first determinant and then $O(n^2)$ (instead of the usual $O(n^3)$) for each subsequent one.</p> <p>Normwise stability should be ensured since we are only making orthogonal transformations.</p> http://mathoverflow.net/questions/109792/nonlinear-matrix-equation/109984#109984 Answer by Federico Poloni for Nonlinear matrix equation Federico Poloni 2012-10-18T07:50:41Z 2012-10-18T07:50:41Z <p>First of all, note that $w^TAv=v^TAw$ is a scalar.</p> <p>Here is an idea that should greatly simplify the equation:</p> <p>Your equations say that $Aw$ and $Av$ are both contained in $U=\operatorname{span}(v,w)$, therefore $U$ is an invariant subspace of $A$. You can get all two-dimensional invariant subspaces by taking $U=\operatorname{span}(x_1,x_2)$, where $x_1$ and $x_2$ are eigenvectors of $A$ (proof: consider $A$ restricted to the subspace $U$; it is a symmetric linear operator, so it has two eigenvalues which are also eigenvalues of $A$).</p> <p>So all solutions must be of the form $v=\alpha x_1 +\beta x_2$ and $w=\gamma x_1+\delta x_2$, where $x_1$ and $x_2$ are two eigenvectors of $A$. Making this ansatz the problem becomes a $2\times 2$ one in $\alpha,\beta,\gamma,\delta$ and should be easy to solve explicitly.</p> http://mathoverflow.net/questions/109845/what-dose-rank-two-modification-w-r-t-matrix-mean/109857#109857 Answer by Federico Poloni for What dose "Rank-two Modification w.r.t Matrix" Mean? Federico Poloni 2012-10-16T22:26:42Z 2012-10-16T22:26:42Z <p>A matrix $B$ is said to be a rank-$k$ modification of $A$ if $B-A$ has rank $k$ (or sometimes rank <em>at most</em> $k$).</p> http://mathoverflow.net/questions/109569/constructing-equivalent-algebraic-expressions-for-matrix-equations/109588#109588 Answer by Federico Poloni for Constructing equivalent algebraic expressions for matrix equations Federico Poloni 2012-10-14T07:59:50Z 2012-10-14T07:59:50Z <p>Another approach would be an improved version of what René Pannekoek suggests in a comment: compute the Schur form (not the Jordan form, which is hopelessly numerically unstable) of $A=QTQ^H$; get for free the Schur form of $A+kI=Q(T+kI)Q^H$; set $y=Q^H x$; rewrite as $f(x)=y^H(T+kI)^{-1}Ty$. This lowers your computational cost for evaluating the function in $K$ different points from $O(n^3K)$ to $O(n^3+n^2K)$ (precompute the Schur form, and do for each point one matrix product with $Q$, one with $T$, one triangular system to solve with $T+kI$).</p> <p>I did not check the details, but I assume that Dima Pasechnik's method is $O(n^4+Kn)$, since you need $O(n)$ matmuls to precompute the Padé approximation. If you combine the two ideas (first reduce to Schur form, then compute a Padé approximation of $g(y)=y^H(T+kI)^{-1}Ty$), I think that you can get $O(n^3+nK)$.</p> <p>I am afraid that the Padé part could be more numerically unstable though. Since $f(k)$ has poles in the eigenvalues of $A$, I guess that the denominator of the order-$n$ Padé approximant is the characteristic polynomial of $A$, and generally working with characteristic polynomials isn't a good idea numerically.</p> http://mathoverflow.net/questions/85251/non-computational-software-useful-to-mathematicians Non-computational software useful to mathematicians Federico Poloni 2012-01-09T12:23:17Z 2012-08-28T08:21:05Z <p>The MathOverflow question <a href="http://mathoverflow.net/questions/19046/open-source-mathematical-software" rel="nofollow">Open source mathematical software</a> contains a list of programs that are useful to perform various computational tasks, such as computer algebra systems.</p> <p>However, evaluating complicated formulas is not all that a professional mathematician needs to do. For instance, another important part of it is <em>communicating</em> results, producing papers and slides. There was a Mathoverflow question devoted specifically to <a href="http://mathoverflow.net/questions/3044/tools-for-collaborative-paper-writing" rel="nofollow">tools for collaborative paper writing.</a></p> <p>It would be useful to identify and gather on MO a list of <strong>areas of activities</strong> where research mathematicians can use software as part of their professional activity.</p> <p>I think it would also be useful to gather on MO a list of <strong>software</strong> that does not strictly do <em>computations</em>, but is nevertheless useful to those who research and teach mathematics.</p> <p>The first example that springs to mind is of course $\LaTeX$, but there is much more:</p> <ul> <li>citation and literature management software, such as Jabref, Zotero, Mendeley</li> <li>conference management software, such as <a href="http://pkp.sfu.ca/?q=ocs" rel="nofollow">Open Conference Systems</a> (never actually used it, but it seems interesting)</li> <li>reference tools such as the <a href="https://oeis.org/" rel="nofollow">Online Encyclopedia of Integer Sequences</a> and <a href="http://pi.lacim.uqam.ca/eng/" rel="nofollow">Plouffe's Inverter</a></li> <li>Euclidean geometry software such as Geogebra</li> <li>diff'ing and merging tools, such as <a href="http://www.ctan.org/tex-archive/support/latexdiff" rel="nofollow">latexdiff</a>, <a href="http://kdiff3.sourceforge.net/" rel="nofollow">kdiff3</a></li> <li>specific LaTeX packages such as Beamer for presentations or several drawing packages (TikZ, Eukleides).</li> </ul> <p>I find some of those real gems, and I'd like to find out more examples. So, the question is:</p> <blockquote> <p>Can you provide (other) examples of programs that are useful to professional mathematicians in their job, while not being strictly speaking "software that does complicated computations"?</p> </blockquote> <p>The question is a little broad and perhaps if there is much software relevant to a specific activity it will be wise to ask, based on input given to the question, a more specific question. Also please do not interpret this question too broadly (see discussion below).</p> <p>Also check these MO questions: <a href="http://mathoverflow.net/questions/3044/tools-for-collaborative-paper-writing" rel="nofollow">Tools for collaborative paper-writing</a> (mainly regarding revision control software), <a href="http://mathoverflow.net/questions/2147/most-helpful-math-resources-on-the-web" rel="nofollow">Most helpful math resources on the web</a> (mainly regarding online databases of something).</p> http://mathoverflow.net/questions/105384/name-for-a-matrix-operation/105386#105386 Answer by Federico Poloni for name for a matrix operation Federico Poloni 2012-08-24T13:23:50Z 2012-08-24T13:23:50Z <p>it is called "diagonal congruence" <a href="http://dx.doi.org/10.1080/03081080600872327" rel="nofollow">here</a>. This makes sense, at least when $D$ is real, since it is a <a href="http://en.wikipedia.org/wiki/Matrix_congruence" rel="nofollow">congruence</a>. "Conjugate" sounds more like $D^{-1}AD$ or $\overline{A}$ to me.</p> http://mathoverflow.net/questions/104803/optimizing-the-condition-number/104834#104834 Answer by Federico Poloni for Optimizing the condition number Federico Poloni 2012-08-16T13:52:07Z 2012-08-16T13:52:07Z <p>I have been working recently on a problem that appears to be related, namely, maximizing the absolute value of the determinant of a chosen $m\times m$ submatrix $S$ of a $m\times N$ submatrix $V$ (think to the columns of $V$ as your vectors). I can tell you the following.</p> <ul> <li>Finding the maximum volume submatrix is an <a href="http://dx.doi.org/10.1007/s00453-011-9582-6" rel="nofollow">NP-hard problem</a>, so I guess that your problem may suffer the same fate</li> <li>a useful relaxation, however, is finding a submatrix that has <em>locally</em> maximum volume, i.e., larger than all those that can be obtained by changing one vector only. There is a <a href="http://dx.doi.org/10.1080/03081088508817636" rel="nofollow">paper</a> by Knuth (yes, <em>that</em> Knuth) that studies the problem. The interesting feature is that in this way the matrix $S^{-1}V$ has a submatrix equal to the identity matrix (obvious) and <em>all other entries bounded in modulus by 1</em>. This is enough to ensure good conditioning, at least in the applications that we were <a href="http://opus4.kobv.de/opus4-matheon/frontdoor/index/index/docId/906" rel="nofollow">investigating</a>. One can prove that the rows of the matrix $V$ are a well-conditioned basis for its row space, for instance. Another reference is <a href="http://spring.inm.ras.ru/osel/?p=13" rel="nofollow">this paper</a>.</li> <li>a second relaxation is finding a submatrix $S$ such that $S^{-1}V$ has all entries bounded by some real $\tau>1$. This problem can be solved explicitly in $O(Nm^2\frac{\log m}{\log\tau})$. The good conditioning properties carry over to this relaxation, up to a factor $\tau$. Putting everything together, I think that one can prove using the techniques in our paper that the matrix chosen by the algorithm has conditioning within $O(\sqrt{Nm}\tau)$ from the conditioning of the rectangular matrix $V$ (defined as largest over smallest singular value), which is a theoretical maximum for the conditioning of $S$ (I think).</li> <li>if you have a practical computation and you want to check how things work with this solution, I have some <a href="https://bitbucket.org/fph/pgdoubling/overview" rel="nofollow">Matlab code</a> in a packaged and usable state that you might wish to try out. Feel free to ask for explanation if the documentation is too poor.</li> </ul> http://mathoverflow.net/questions/104810/stability-of-levinson-durbin-method-for-toeplitz-system-solutions/104818#104818 Answer by Federico Poloni for Stability of Levinson-Durbin method for Toeplitz system solutions ? Federico Poloni 2012-08-16T07:35:23Z 2012-08-16T07:35:23Z <p>It can be unstable, and the condition number is not an adequate measure to tell when it fails. You may want to check the numerical experiments in <a href="http://www.jstor.org/stable/2153371" rel="nofollow">http://www.jstor.org/stable/2153371</a> for specific examples and discussion.</p> <p>If I remember correctly, you can prove stability only for symmetric matrices whose principal minors are all positive (totally positive matrices), otherwise things can go wrong. More precisely, there are numbers called <em>reflection coefficients</em> that can be defined along the algorithm; you get stability only if they are all positive, as for instance in the totally positive case (sorry but I do not remember the details exactly, it's some years since I worked on that).</p> <p>The GKO algorithm, developed in the article I linked above, is generally considered to be a more stable alternative. However, you have to pay for that with larger computational times and memory overhead (you can save the memory overhead though).</p> http://mathoverflow.net/questions/103294/solving-non-linear-equations/103296#103296 Answer by Federico Poloni for solving non linear equations Federico Poloni 2012-07-27T11:44:04Z 2012-07-27T11:44:04Z <p>To answer completely, one should know what is the application that you have in mind and if there is any specific convergence result for the Newton method.</p> <p>But, in general, I think that you are right: solving an inexact model to full precision is not needed. If you are solving an inexact problem $P'$ whose solution $x'$ is at distance $d$ from the solution $x$ to the true problem $P$, then it is useless to determine $x'$ more precisely than within a radius $O(d)$. However, a good deal of error analysis is required to determine precisely how much imprecise you can be in solving the intermediate problems.</p> <p>Of course this is only needed if you wish to prove things in a rigorous setting (and if your problem allows it): if, as it often happens in engineering, the approach is "run Newton, cross fingers and it will converge to the solution", then all these strategies are just heuristics, too, and nothing more precise can be proved.</p> http://mathoverflow.net/questions/101651/inverse-m-matrix-times-mixed-sign-vector inverse M-matrix times mixed-sign vector Federico Poloni 2012-07-08T12:00:54Z 2012-07-08T12:00:54Z <p>Recently a colleague and I came across this unusual phenomenon.</p> <p>Take $M\in\mathbb{R}^{n\times n}$ a singular irreducible M-matrix, and $b\in\mathbb{R}^{n}$ such that the system $Mx=b$ is solvable (so, in general, $b$ is going to contain both positive and negative entries). For each $i=1,2,\dots,n$, let $M^{(i)}\in\mathbb{R}^{n-1\times n-1}$ and $b^{(i)}\in\mathbb{R}^{n-1}$ be the matrices and vectors obtained by eliminating the $i$th row and column from $M$ and $b$. Then, at least one of the (nonsingular) systems $M^{(i)}x^{(i)}=b^{(i)}$ has a solution $x^{(i)}\geq 0$.</p> <p>The proof in itself is quite easy: notice that the $x^{(i)}$ correspond to solutions of $Mx=b$ with one zero component; the solutions of this system can be written as $\alpha z +y$, with $z,y\geq 0$, and by minimality there is an $\alpha$ that makes one component zero and the rest positive.</p> <p>However, I find this unusual because the typical argument in this field is finding out that a vector is nonnegative if it can be written as $A^{-1}c$, with $A^{-1}$ and $c$ nonnegative. Here we have a $c$ with mixed signs instead.</p> <p>What I would be interested in is a better characterization of the entry(or entries) $i$ that has $x^{(i)}\geq 0$. I hope it satisfies some minimality property more meaningful than the obvious $i=\arg\min \min y_i/z_i$. The ideal for me would be finding a characterization that lets me operate on $M$ and $b$ only and be able to tell what is the "correct" entry to take.</p> <p>Is this topic already studied? Maybe the magic $i$ corresponds to some special property of the graph associated to $M$.</p> http://mathoverflow.net/questions/101527/singular-value-decomposition/101565#101565 Answer by Federico Poloni for singular value decomposition Federico Poloni 2012-07-07T09:40:13Z 2012-07-07T09:40:13Z <p>Let $A_2A_1^{-1}=U_A DV_A$ and $B_2B_1^{-1}=U_B DV_B$ be the two SVDs with the same $D$.</p> <p>Set $U_2'=U_AU_B'$, $U_1=V_B'V_A$, $V=B_1^{-1} U_1 A_1$.</p> <p>The first equality is clear. The second one is proved by $$U_2'B_2V=U_A U_B'B_2 B_1^{-1} V_B' V_A A_1=U_ADV_A A_1=A_2A_1^{-1}A_1=A_2.$$</p> http://mathoverflow.net/questions/130916/how-to-maximize-the-determinant-of-a-matrix-of-the-form-vdvh Comment by Federico Poloni Federico Poloni 2013-05-17T07:49:02Z 2013-05-17T07:49:02Z Not sure about that: $V^HV$ is a rank-$M$ $2M\times 2M$ matrix, so its determinant would be zero. http://mathoverflow.net/questions/129741/spectral-radius-monotonicity/129789#129789 Comment by Federico Poloni Federico Poloni 2013-05-15T07:18:29Z 2013-05-15T07:18:29Z You are right! One needs an additional hypothesis (irreducibility) to get strict inequalities in P-F. http://mathoverflow.net/questions/130574/why-do-i-get-estimated-error-1-ind-when-doing-bicgstab-linear-solver-using-il Comment by Federico Poloni Federico Poloni 2013-05-14T12:58:42Z 2013-05-14T12:58:42Z I think it's off-topic even there. This is for the program's mailing list or in its manual; the point is deciphering an obscure error message, there is no mathematics or science involved here. http://mathoverflow.net/questions/130520/how-to-use-a-numerical-method-similar-to-newtons-method-to-solve-multiple-roots Comment by Federico Poloni Federico Poloni 2013-05-14T06:39:59Z 2013-05-14T06:39:59Z You are probably misunderstanding the term &quot;multiple roots&quot;. By that, we usually mean &quot;roots with multiplicity higher than 1&quot;. If you wish to find several roots of the same equation, Newton's method has no guarantees to find them all; it might almost always converge to the same one. http://mathoverflow.net/questions/120397/eigenvalues-of-random-hamiltonian-matrices/130458#130458 Comment by Federico Poloni Federico Poloni 2013-05-13T13:31:49Z 2013-05-13T13:31:49Z I must be missing something then. Doesn't conjecture (1) (proved in a paper) say that they should be uniformly distributed in the limit? http://mathoverflow.net/questions/120397/eigenvalues-of-random-hamiltonian-matrices/130458#130458 Comment by Federico Poloni Federico Poloni 2013-05-13T08:47:38Z 2013-05-13T08:47:38Z It is curious that there seems to be an area with less eigenvalues around the imaginary axis (for Hamiltonian) or the real axis (for Ginibre). Is it an artifact of the numerical methods? I can imagine them having trouble in identifying properly purely imaginary/real eigenvalues. http://mathoverflow.net/questions/88725/spectral-properties-of-ai-a-1/88727#88727 Comment by Federico Poloni Federico Poloni 2013-05-06T21:12:13Z 2013-05-06T21:12:13Z T is triangular. <a href="https://en.wikipedia.org/wiki/Schur_form" rel="nofollow">en.wikipedia.org/wiki/Schur_form</a> http://mathoverflow.net/questions/129636/can-anyone-recommend-a-package-in-c-which-has-most-efficiency-in-doing-sparse-m Comment by Federico Poloni Federico Poloni 2013-05-04T12:43:17Z 2013-05-04T12:43:17Z I can recommend you to ask on <a href="http://scicomp.stackexchange.com/" rel="nofollow">scicomp.stackexchange.com</a>; this question is off-topic here. http://mathoverflow.net/questions/129277/problems-where-conjugate-gradient-works-much-better-than-gmres Comment by Federico Poloni Federico Poloni 2013-05-01T07:14:34Z 2013-05-01T07:14:34Z From the look of it, yes, but there is probably a completely theoretical question in matrix approximation theory hidden behind this one. http://mathoverflow.net/questions/128993/on-solution-of-a-discrete-time-equation/129019#129019 Comment by Federico Poloni Federico Poloni 2013-04-29T08:55:18Z 2013-04-29T08:55:18Z For a general linear matrix equation with more than two terms, I am not aware of any method faster than the $O(n^6)$ of the naive &quot;vectorize everything&quot; algorithm. (caveat: I know almost nothing of LMIs, so there might be a method in that area that I do not know about). Maybe you'll have more luck posting this as a separate question, either here or on <a href="http://scicomp.stackexchange.com/" rel="nofollow">scicomp.stackexchange.com</a>. http://mathoverflow.net/questions/128993/on-solution-of-a-discrete-time-equation/129019#129019 Comment by Federico Poloni Federico Poloni 2013-04-29T07:10:14Z 2013-04-29T07:10:14Z it more stable, and ultimately obtain a discrete-time version of the ADI method for Lyapunov equations. In short, with some algebraic manipulations everything that works in the continuous-time case rates to work in the discrete-time case as well. http://mathoverflow.net/questions/128993/on-solution-of-a-discrete-time-equation/129019#129019 Comment by Federico Poloni Federico Poloni 2013-04-29T07:09:07Z 2013-04-29T07:09:07Z Yep. <code>dlyap</code> will be fine for a small-scale problem; it should use a variant of the same Bartels-Stewart method that is used for continous-time lyapunov eqs; essentially, take a Schur form of $F$ and solve directly via back-substitution for each entry of $X$ &quot;in the right order&quot;. For large-scale problems, you can truncate the series $X=\sum_{i=0}^{\infty} F^{i}CF^{Ti}$, or obtain the partial sum truncated at the term $2^{k}$ directly from the one truncated at $2^{k-1}$ with some manipulations (Smith methods). You can apply some M&#246;bius transforms to $F$ without changing the solution to make... http://mathoverflow.net/questions/128993/on-solution-of-a-discrete-time-equation/129019#129019 Comment by Federico Poloni Federico Poloni 2013-04-28T18:37:08Z 2013-04-28T18:37:08Z There are numerical methods dealing directly with the equation $X-FXF^T=C$, known as discrete-time Lyapunov equation, or Stein equation. There is no need to invert $F$ and convert it to a Sylvester; sometimes algorithms even go the other way round and convert a continuous-time Lyapunov equation to a discrete-time one to solve it. http://mathoverflow.net/questions/33353/naming-for-the-map-t-x-mapsto-a-x-b Comment by Federico Poloni Federico Poloni 2013-04-25T06:51:45Z 2013-04-25T06:51:45Z Since you are back here: thought that the name &quot;discrete-time Sylvester operator&quot; makes a lot of sense for the map $S:X\to X-AXB$, since a discrete-time Sylvester equation can be written $S(x)=C$, analogously with other similarly-named operator. The name with commas returns a whopping two Google results, but I think that a person working in matrix equation would understand immediately what you're talking about if you use it. http://mathoverflow.net/questions/128285/the-first-eigenvalue-of-a-branching-process-matrix/128320#128320 Comment by Federico Poloni Federico Poloni 2013-04-22T09:02:31Z 2013-04-22T09:02:31Z You are right, sorry for forgetting this part!