[petsc-users] Solution of normal equations

[Jose E. Roman]

I guess you could try LSQR on the matrix M=[ sqrt(D); sqrt(R)*J ], given that M'*M=D+J'*R*J. This will attempt to solve a least squares system with coefficient matrix M, without forming the ill-conditioned M'*M. But the question is how to update your right-hand-side vector so that it solves the same system of equations that you want to solve.

Jose

> El 19 jul 2022, a las 19:58, Mark Adams <mfadams at lbl.gov> escribió:
> 
> This matrix is SPD so you want to use -ksp_type cg
> 
> Without a matrix, preconditioning is hard. You can give it a shot. You need to make a MatShell and use -pc_type none.
> 
> You will probably want to make a matrix and you should research solver methods that have been used for problems like yours.
> 
> Mark
> 
> On Tue, Jul 19, 2022 at 12:07 PM Alejandro Castro <alejandro.castro at tri.global> wrote:
> Hello PETSc users,
> 
> I am looking to solve a (sparse) system of equations H⋅x=b where matrix H has the form H=D+Jᵀ⋅R⋅J, where both D and R are positive diagonal matrices. This is very similar to a set of normal equations, though not quite. J is square and most often with a large condition number κ(J) (and thus H's condition number will scale with the square of κ(J)).
> 
> What method would you recommend for solving H⋅x=b? In particular, I imagine I do not want to form H directly. Also, what iterative method would you recommend given the large condition number?
> 
> Thank you in advance,
> 
> A.
> 
> 
> 
> -- 
> 	
> Alejandro Castro (he/him)
> Manager, Dynamics & Simulation
> TOYOTA RESEARCH INSTITUTE
> Los Altos / Cambridge
> 
> 
> Confidential or protected information may be contained in this email and/or attachment. Unless otherwise marked, all TRI email communications are considered "PROTECTED" and should not be shared or distributed. Thank you.