[petsc-users] Why am I getting bad scalability for least squares problems?
Matthew Knepley
knepley at gmail.com
Tue Jan 29 19:57:30 CST 2013
On Tue, Jan 29, 2013 at 8:11 PM, Gaurish Telang <gaurish108 at gmail.com>wrote:
> Hello,
>
> I am trying to solve some over-determined / under-determined least-squares
> systems in parallel using the LSQR routine of PETSc.
> In spite of being able to solve the least squares problem || Ax-b ||
> _correctly_ with this Krylov method on multiple processors, I don't seem
> to be acheiving scalability with respect to wall-clock time
> for the set of matrices that I am interested in.
>
> Below, I have listed the matrix sizes, the number of non-zeros in each
> matrix and the wall-clock times to required for 10 iterations of the LSQR
> krylov method.
> The timings have been amortized over 50 iterations of the solver i.e. I
> solve the least-square problem 50 times ( where each time 10 iterations of
> LSQR are carried out)
> and obtain the arithmetic mean of wall-clock times so recorded (i.e
> (sum-of-clock-times)/50).
>
Never ask performance questions without sending the output of -log_summary.
Also:
http://www.mcs.anl.gov/petsc/documentation/faq.html#computers
http://www.mcs.anl.gov/petsc/documentation/faq.html#slowerparallel
Matt
> Matrix size: 28407 x 19899
> Non-zeros: 725363
> Wall clock times
> 1 proc: 0.094 s
> 2 proc: 0.087 s
> 3 proc: 0.089 s
> 4 proc: 0.093 s
>
> Matrix size: 95194 x 155402
> Non-zeros: 3877164
> Wall clock times :
> 1 proc: 0.23 s
> 2 proc: 0.21 s
> 3 proc: 0.22 s
> 4 proc: 0.23 s
>
> Matrix size: 125878 x 207045
> Non-zeros: 3717995
> Wall clock times
> 1 proc: 0.24 s
> 2 proc: 0.22 s
> 3 proc: 0.23 s
> 4 proc: 0.24 s
>
> I have other matrices which show similar bad scaling as I increase the
> processor count.
>
> Please let me know what I can do to improve the performance after I
> increase the processor count.
> I feel it is a data-layout problem i.e. I need to chose some other
> data-structure with PETSc to represent the matrix and vector of my
> least-squares problem.
>
> Currently my Matrix type is "mpiaij" and Vec type is "mpi" which I set at
> the command-line while running the executable.
>
> I use the terminal command :
> mpirun -np <num-procs> ./test_parallel
> -fmat <binary-mat-file>
> -frhs <binary-vec-file>
> -vec_type mpi
> -mat_type mpiaij
> -ksp_type lsqr
> -ksp_max_it 10
>
>
> I have pasted the code I am using below for your reference.
>
> Thank you,
> Gaurish.
>
>
>
>
>
>
>
> ==============================================================================================================================
>
> #include <petscksp.h>
> #include <petscsys.h>
> #include <petsctime.h>
> #include <stdio.h>
> #undef __FUNCT__
> #define __FUNCT__ "main"
> int main(int argc,char **args)
> {
> Vec x, b, residue; /* approx solution, RHS, residual */
> Mat A; /* linear system matrix */
> KSP ksp; /* linear solver context */
> PC pc; /* preconditioner context */
> PetscErrorCode ierr;
> PetscInt m,n ; /* # number of rows and columns of the
> matrix read in*/
> PetscViewer fd ;
> PetscInt size;
> //tscInt its;
> //PetscScalar norm, tol=1.e-5;
> PetscBool flg;
> char fmat[PETSC_MAX_PATH_LEN]; /* input file names */
> char frhs[PETSC_MAX_PATH_LEN]; /* input file names */
>
> PetscInitialize(&argc,&args,(char *)0,help);
> ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr);
>
> ierr =
> PetscOptionsGetString(PETSC_NULL,"-fmat",fmat,PETSC_MAX_PATH_LEN,&flg);
> ierr =
> PetscOptionsGetString(PETSC_NULL,"-frhs",frhs,PETSC_MAX_PATH_LEN,&flg);
>
>
> /* Read in the matrix */
> ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD, fmat, FILE_MODE_READ, &fd
> ); CHKERRQ(ierr);
> ierr = MatCreate(PETSC_COMM_WORLD,&A);
> // ierr = MatSetType(A,MATSEQAIJ);
> ierr = MatSetFromOptions(A);CHKERRQ(ierr);
> ierr = MatLoad(A,fd);
> ierr = PetscViewerDestroy(&fd);
>
>
> /*Read in the right hand side*/
> ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD, frhs, FILE_MODE_READ, &fd
> ); CHKERRQ(ierr);
> ierr = VecCreate(PETSC_COMM_WORLD,&b);
> ierr = VecSetFromOptions(b);CHKERRQ(ierr);
> ierr = VecLoad(b,fd);
> ierr = PetscViewerDestroy(&fd);
>
>
> /* Get the matrix size. Used for setting the vector sizes*/
> ierr = MatGetSize(A , &m, &n);
> //printf("The size of the matrix read in is %d x %d\n", m , n);
>
>
> VecCreate(PETSC_COMM_WORLD, &x);
> VecSetSizes(x, PETSC_DECIDE, n);
> VecSetFromOptions(x);
>
> VecCreate(PETSC_COMM_WORLD, &residue);
> VecSetSizes(residue, PETSC_DECIDE, m);
> VecSetFromOptions(residue);
>
>
> /* Set the solver type at the command-line */
> ierr = KSPCreate(PETSC_COMM_WORLD,&ksp);CHKERRQ(ierr);
> ierr = KSPSetOperators(ksp,A,A,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr);
> ierr = KSPGetPC(ksp,&pc);CHKERRQ(ierr);
> ierr = PCSetType(pc,PCNONE);CHKERRQ(ierr);
> ierr =
> KSPSetTolerances(ksp,1.e-5,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);CHKERRQ(ierr);
> /* Set runtime options, e.g.,
> -ksp_type <type> -pc_type <type> -ksp_monitor -ksp_rtol <rtol>
> These options will override those specified above as long as
> KSPSetFromOptions() is called _after_ any other customization
> routines */
> ierr = KSPSetFromOptions(ksp);CHKERRQ(ierr);
>
> /* Initial guess for the krylov method set to zero*/
> PetscScalar p = 0;
> ierr = VecSet(x,p);CHKERRQ(ierr);
> ierr = KSPSetInitialGuessNonzero(ksp,PETSC_TRUE);CHKERRQ(ierr);
> /* Solve linear system */
>
> PetscLogDouble v1,v2,elapsed_time[50];
> int i;
> for (i = 0; i < 50; ++i)
> {
> ierr = PetscGetTime(&v1);CHKERRQ(ierr);
> ierr = KSPSolve(ksp,b,x);CHKERRQ(ierr);
> ierr = PetscGetTime(&v2);CHKERRQ(ierr);
> elapsed_time[i] = v2 - v1;
> PetscPrintf(PETSC_COMM_WORLD,"[%d] Time for the solve: %g s\n", i
> , elapsed_time[i]);
> }
>
> PetscLogDouble sum=0,amortized_time ;
> for ( i = 0; i < 50; ++i)
> {
> sum += elapsed_time[i];
> }
> amortized_time = sum/50;
> PetscPrintf(PETSC_COMM_WORLD,"\n\n***********************\nAmortized
> Time for the solve: %g s\n***********************\n\n", amortized_time);
>
>
> /* View solver info; we could instead use the option -ksp_view to print
> this info to the screen at the conclusion of KSPSolve().*/
> ierr = KSPView(ksp,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
>
> /* Clean up */
> ierr = VecDestroy(&x);CHKERRQ(ierr);
> ierr = VecDestroy(&residue);CHKERRQ(ierr);
> ierr = VecDestroy(&b);CHKERRQ(ierr);
> ierr = MatDestroy(&A);CHKERRQ(ierr);
> ierr = KSPDestroy(&ksp);CHKERRQ(ierr);
>
> ierr = PetscFinalize();
> return 0;
> }
>
>
>
>
--
What most experimenters take for granted before they begin their
experiments is infinitely more interesting than any results to which their
experiments lead.
-- Norbert Wiener
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20130129/3306108d/attachment-0001.html>
More information about the petsc-users
mailing list