> ------------------------------------------------------------------------ > [0]PETSC ERROR: MatILUFactorSymbolic() line 4573 in > src/mat/interface/e:\cygwin\ > home\PETSC-~1.1-P\src\mat\INTERF~1\matrix.c > [0]PETSC ERROR: No support for this operation for this object type! > [0]PETSC ERROR: Matrix type sbaij symbolic ILU!
Using sbaij format, you must use preconditioner icc, not ilu. Hong > [0]PETSC ERROR: PCSetUp_ILU() line 532 in > src/ksp/pc/impls/factor/ilu/e:\cygwin\ > home\PETSC-~1.1-P\src\ksp\pc\impls\factor\ilu\ilu.c > [0]PETSC ERROR: PCSetUp() line 798 in > src/ksp/pc/interface/e:\cygwin\home\PETSC- > ~1.1-P\src\ksp\pc\INTERF~1\precon.c > [0]PETSC ERROR: KSPSetUp() line 234 in > src/ksp/ksp/interface/e:\cygwin\home\PETS > C-~1.1-P\src\ksp\ksp\INTERF~1\itfunc.c > [0]PETSC ERROR: KSPSolve() line 334 in > src/ksp/ksp/interface/e:\cygwin\home\PETS > C-~1.1-P\src\ksp\ksp\INTERF~1\itfunc.c > > > -----Original Message----- > > From: owner-petsc-users at mcs.anl.gov > > [mailto:owner-petsc-users at mcs.anl.gov] On Behalf Of Hong Zhang > > Sent: Wednesday, October 11, 2006 2:14 PM > > To: petsc-users at mcs.anl.gov > > Subject: RE: General matrix interative solver > > > > > > Julian, > > > > If you do not perform matrix reordering, you can use petsc > > sbaij format, which only stores the upper triangular entries. > > aij format has more algorithmic options and is more efficient > > and supports matrix reordering. > > > > Hong > > > > On Wed, 11 Oct 2006, Julian wrote: > > > > > I'm using '-ksp_type cg' now to solve the symmetric > > matric...But then > > > do I need to store all the non-zero elements of a symmetric > > matrix ? > > > i.e. A[i,j] as well as A[j,i] ? > > > > > > Shouldn't I be able to store just the upper or lower > > triangle or the > > > matrix ? Is that possible with petsc ? > > > > > > Julian. > > > > > > > > > > -----Original Message----- > > > > From: owner-petsc-users at mcs.anl.gov > > > > [mailto:owner-petsc-users at mcs.anl.gov] On Behalf Of Julian > > > > Sent: Wednesday, October 11, 2006 12:37 PM > > > > To: petsc-users at mcs.anl.gov > > > > Subject: RE: General matrix interative solver > > > > > > > > Barry, > > > > > > > > I tried sending the commands from a file for now... And > > once I used > > > > '-ksp_type cg', I got pretty much the same performance as that of > > > > the inhouse code ! > > > > Now, I'll try the performance with larger cases... > > > > > > > > Thanks for your help! > > > > > > > > Julian. > > > > > > > > > > > > > -----Original Message----- > > > > > From: owner-petsc-users at mcs.anl.gov > > > > > [mailto:owner-petsc-users at mcs.anl.gov] On Behalf Of Barry Smith > > > > > Sent: Wednesday, October 11, 2006 11:25 AM > > > > > To: petsc-users at mcs.anl.gov > > > > > Subject: RE: General matrix interative solver > > > > > > > > > > > > > > > PetscOptionsSetValue() but recommend putting them on > > the command > > > > > line or in a file (pass the name of the file into > > > > PetscInitialize()). > > > > > Having to recompile for every option is too painful. > > > > > > > > > > Barry > > > > > > > > > > > > > > > On Wed, 11 Oct 2006, Julian wrote: > > > > > > > > > > > Hong, > > > > > > > > > > > > The times are for the solution process alone and the > > > > > initial guess is > > > > > > the same, i.e. zero. > > > > > > > > > > > > As for the algorithm, they are probably different. > > > > > > Can you tell me how to pass 'runtime options' like -help > > > > > and -ksp_view > > > > > > from within the code... I mean, at compile time. > > > > > > > > > > > > thanks > > > > > > > > > > > > > > > > > > > > Julian, > > > > > > > > > > > > > > When comparing the inhouse solver and petsc solver, you > > > > need make > > > > > > > sure > > > > > > > > > > > > > > 1. Timings are collected for solution process. The matrix > > > > > and vector > > > > > > > assembly should be excluded. > > > > > > > 2. They should use same iterative algorithm. By default, > > > > > petsc uses > > > > > > > gmres with restart=30 and ilu(0) preconditioner. > > > > Petsc supports > > > > > > > symmetric matrices, e.g., runtime option '-ksp_type cg > > > > > -pc_type > > > > > > > icc' > > > > > > > might give better performance 3. They should start from > > > > > > > same intial guess. By default, petsc > > > > > > > initial guess is zero. > > > > > > > > > > > > > > You can use '-ksp_view' to see what algorithm and options > > > > > are used > > > > > > > in petsc. > > > > > > > > > > > > > > Hong > > > > > > > > > > > > > > On Wed, 11 Oct 2006, Julian wrote: > > > > > > > > > > > > > > > Hello, > > > > > > > > > > > > > > > > I implemented the iterative sparse matrix solver in PetSc > > > > > > > into my FEM > > > > > > > > code recently. I compared the results from a problem > > > > with 1317 > > > > > > > > unknowns. I used a direct solver to obtain the reference > > > > > > > solution. I > > > > > > > > have another in-house sparse iterative solver that I have > > > > > > > been using > > > > > > > > so far. It was written by someone else but I have > > access to > > > > > > > the source for that solver. > > > > > > > > I find the 'error norm' in the solution by taking the > > > > > > > square root of > > > > > > > > the sum of the squares of the absolute differences > > > > between the > > > > > > > > solution from the direct solver and the iterative > > > > solver. I am > > > > > > > > ignoring the numerical zeros in the solutions > > when doing this. > > > > > > > > I find that in order to get same order of the error > > > > > norm (1e-13) > > > > > > > > as the in-house iterative solver, the petsc solver > > > > takes a much > > > > > > > > longer time and larger number of iterations. While > > > > the inhouse > > > > > > > > solver took less than one second, the petsc > > solver took 13 > > > > > > > > seconds. The inhouse solver took 476 iterations whereas > > > > > the petsc > > > > > > > > solver took > > > > > > > 4738 iterations. > > > > > > > > I'm guessing this has to do with different setting of > > > > > the solver > > > > > > > > in petsc such as the preconditioner etc. > > > > > > > > Can you tell me what the different settings are? > > And how to > > > > > > > tweak them > > > > > > > > so that I can atleast get as good as a performance as the > > > > > > > inhouse code ? > > > > > > > > Given below is how I have implemented the petsc solver: > > > > > > > > > > > > > > > > /////initialization > > > > > > > > PetscInitializeNoArguments(); > > > > > > > > Assert( mat = (double*)malloc(sizeof(Mat)) ); > > > > > > > > MatCreateSeqAIJ(PETSC_COMM_SELF, L, L, > > > > > > > > PETSC_DEFAULT, PETSC_NULL, (Mat*)mat); > > > > > > > > > > > > > > > > ////// this is the function I use to populate the matrix > > > > > > > > MatSetValue(*(Mat*)mat, ii, jj, value, ADD_VALUES); > > > > > > > > > > > > > > > > ////// this is how I actaully solve the matrix > > > > > > > > MatAssemblyBegin(*(Mat*)mat, MAT_FINAL_ASSEMBLY); > > > > > > > > MatAssemblyEnd(*(Mat*)mat, MAT_FINAL_ASSEMBLY); > > > > > > > > > > > > > > > > double iter_error = 1e-10; > > > > > > > > int max_iter_num = 10000; > > > > > > > > int num_of_iter; > > > > > > > > > > > > > > > > Vec rhs, x; > > > > > > > > VecCreateSeqWithArray(PETSC_COMM_SELF, L, b, &rhs); > > > > > > > > VecDuplicate(rhs, &x); > > > > > > > > > > > > > > > > KSP ksp; > > > > > > > > KSPCreate(PETSC_COMM_SELF, &ksp); > > > > > > > > KSPSetTolerances(ksp, iter_error, PETSC_DEFAULT, > > > > > > > > PETSC_DEFAULT, max_iter_num); > > > > > > > > KSPSetFromOptions(ksp); > > > > > > > > KSPSetOperators(ksp, *(Mat*)mat, *(Mat*)mat, > > > > > > > > SAME_PRECONDITIONER); > > > > > > > > KSPSolve(ksp,rhs,x); > > > > > > > > > > > > > > > > PetscReal r_norm; > > > > > > > > KSPGetResidualNorm(ksp, &r_norm); > > > > > > > > KSPGetIterationNumber(ksp, &num_of_iter); > > > > > > > > > > > > > > > > cout << "max_iter_num\t" << max_iter_num << endl; > > > > > > > > cout << "iter_error\t" << iter_error << endl; > > > > > > > > > > > > > > > > cout << "Matrix solver step " << num_of_iter << ", > > > > > > > residual " > > > > > > > > << r_norm << ".\n"; > > > > > > > > > > > > > > > > PetscScalar *p; > > > > > > > > VecGetArray(x, &p); > > > > > > > > for(int i=0; i<L; i++) { > > > > > > > > b[i] = p[i]; > > > > > > > > } > > > > > > > > VecRestoreArray(x, &p); > > > > > > > > > > > > > > > > KSPDestroy(ksp); > > > > > > > > VecDestroy(rhs); > > > > > > > > VecDestroy(x); > > > > > > > > > > > > > > > > cout <<"Iterations for convergence="<< num_of_iter << " - > > > > > > > Residual Norm = " > > > > > > > > << r_norm << endl; > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > If this is not the typical method to be used to solve > > > > > this kind of > > > > > > > > problem, please let me know what functions I should use. > > > > > > > > I should mention that the inhouse code is for symmetric > > > > > > > matrices and > > > > > > > > from what I understand, the petsc solver works for general > > > > > > > unsymmetric matrices. > > > > > > > > But I think for iterative solvers, it should still give > > > > > around the > > > > > > > > same performance. > > > > > > > > I tested the solvers against some other problems as > > > > well, and I > > > > > > > > got the same performance.. In some cases, no > > matter how many > > > > > > > iterations it > > > > > > > > goes through, the petsc solver would not go below > > a certain > > > > > > > error norm > > > > > > > > whereas the inhouse solver would get almost exactly the > > > > > > > same answer as > > > > > > > > the direct solver solution. I'm thinking the petsc > > > > > solver should > > > > > > > > be able to solve this problem just as easily. It would be > > > > > > > great if anyone > > > > > > > > could help me figure out the appropriate settings I should > > > > > > > use in the petsc solver. > > > > > > > > > > > > > > > > Thanks, > > > > > > > > Julian. > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > > >
