A strange behavior I am observing is: Problem: I have to solve A*x=rhs, and currently I am currently trying to solve for a system where I know the exact solution. I have initialized the exact solution in the Vec x_exact.
MatMult(A, x_exact, dummy);// Storing the value of A*x_exact in dummy VecAXPY(dummy, -1.0, rhs); // dummy = dummy -rhs VecNorm(dummy, NORM_INFINITY, &norm_val); // norm_val = ||dummy||, which gives us the residual norm PetscPrintf(PETSC_COMM_SELF, "Norm = %f\n", norm_val); // Printing the norm. // Starting with the linear solver KSPCreate(PETSC_COMM_SELF, &ksp); KSPSetOperators(ksp, A, A); KSPSetFromOptions(ksp); KSPSolve(ksp,rhs,x_exact); // Solving the system A*x= rhs, with the given initial input x_exact. So the result will also be stored in x_exact On running with -pc_type lu -pc_factor_mat_solver_package superlu -ksp_monitor I get the following output: Norm = 0.000000 0 KSP Residual norm 4.371606462669e+04 1 KSP Residual norm 5.850058113796e+02 2 KSP Residual norm 5.832677911508e+02 3 KSP Residual norm 1.987386549571e+02 4 KSP Residual norm 1.220006530614e+02 . . . Since the initial guess is the exact solution should'nt the first residual itself be zero and converge in one iteration. Thanks, Kaushik On Tue, Apr 11, 2017 at 10:08 AM, Kaushik Kulkarni <[email protected]> wrote: > Thank you for the inputs. > I tried Barry' s suggestion to use SuperLU, but the solution does not > converge and on doing -ksp_monitor -ksp_converged_reason. I get the > following error:- > 240 KSP Residual norm 1.722571678777e+07 > Linear solve did not converge due to DIVERGED_DTOL iterations 240 > For some reason it is diverging, although I am sure that for the given > system a unique solution exists. > > Thanks, > Kaushik > > On Tue, Apr 11, 2017 at 1:04 AM, Xiaoye S. Li <[email protected]> wrote: > >> If you need to use SuperLU_DIST, the pivoting is done statically, using >> maximum weighted matching, so the small diagonals are usually taken care as >> well. It is not as good as partial pivoting, but works most of the time. >> >> Sherry >> >> On Mon, Apr 10, 2017 at 12:07 PM, Barry Smith <[email protected]> wrote: >> >>> >>> I would suggest using ./configure --download-superlu and then when >>> running the program -pc_type lu -pc_factor_mat_solver_package superlu >>> >>> Note that this is SuperLU, it is not SuperLU_DIST. Superlu uses >>> partial pivoting for numerical stability so should be able to handle the >>> small or zero diagonal entries. >>> >>> Barry >>> >>> > On Apr 10, 2017, at 1:17 PM, Kaushik Kulkarni <[email protected]> >>> wrote: >>> > >>> > Hello, >>> > I am trying to solve a 2500x2500 sparse matrix. To get an idea about >>> the matrix structure I have added a file matrix.log which contains the >>> output of MatView() and also the output of Matview_draw in the image file. >>> > >>> > From the matrix structure it can be seen that Jacobi iteration won't >>> work and some of the diagonal entries being very low(of the order of 1E-16) >>> LU factorization would also fail. >>> > >>> > Can someone please suggest what all could I try next, in order to >>> make the solution converge? >>> > >>> > Thanks, >>> > Kaushik >>> > >>> > -- >>> > Kaushik Kulkarni >>> > Fourth Year Undergraduate >>> > Department of Mechanical Engineering >>> > Indian Institute of Technology Bombay >>> > Mumbai, India >>> > https://kaushikcfd.github.io/About/ >>> > +91-9967687150 >>> > <matrix.log><matrix_pattern.png> >>> >>> >> > > > -- > Kaushik Kulkarni > Fourth Year Undergraduate > Department of Mechanical Engineering > Indian Institute of Technology Bombay > Mumbai, India > https://kaushikcfd.github.io/About/ > +91-9967687150 > -- Kaushik Kulkarni Fourth Year Undergraduate Department of Mechanical Engineering Indian Institute of Technology Bombay Mumbai, India https://kaushikcfd.github.io/About/ +91-9967687150
