May 2009 - petsc-users - mailman3.cels.anl.gov

Some basic questions about PETSC
by Ben Tay 19 Jan '11

19 Jan '11

Hi, I'm still relatively new to PETSC but I've finally managed to solve my own poisson eqn using it in about a week's time.Thank you for writing such a good solver. However, I have some questions still about PETSC: 1. Are petscscaler the same as real nos.? I am writing in double precision or real(8). So are these 2 interchangeable? (same as petscint and integer). 2. while linking in compaq visual fortran, I get the warning message about conflict library with LIBCMT. Is it important or can I just ignore it? 3. I am writing in f90 free form format. When I modify from fixed f77 format to f90, it gives me a lot of error, which seems to be the include file error. So, is it possible to write in f90 format? 4. I am writing a Navier stokes solver and I'm using PETSC to solve the poisson eqn. If only the RHS changes with each time step, do I need to call VecCreateSeq <insert values in RHS b> VecAssemblyBegin VecAssemblyEnd VecDestroy at each time step? I don't think I can (or I should) change values of b after VecAssemblyBegin,VecAssemblyEnd, or can I? Thank you very much and have a nice day Regards, Ben

4 8

Petsc matrices in Matlab
by Manav Bhatia 18 Sep '09

18 Sep '09

Hi, Is there a way to directly read a Petsc sequential sparse matrix in Matlab using a Matlab script? I vaguely remember doing this a couple of years ago, and the script possibly was a part of the Petsc distribution. Is that still the case? I am, otherwise, aware of the methods outlined in chapter 9 of the Petsc users guide. Thanks, Manav

8 15

Mismatch in explicit fortran interface for MatGetInfo
by Stephan Kramer 18 Jun '09

18 Jun '09

Hi all, First of all thanks of a lot for providing explicit fortran interfaces for most functions in Petsc 3. This is of great help. I do however run into a problem using MatGetInfo. The calling sequence for fortran (according to the manual) is: double precision info(MAT_INFO_SIZE) Mat <cid:[email protected]> A integer ierr call MatGetInfo <cid:[email protected]>(A,MAT_LOCAL,info,ierr) The interface however seems to indicate the info argument has to be a single double precision (i.e. a scalar not an array). I guess with implicit interfaces this sort of thing would work, but with the provided explicit interface, at least gfortran won't let me have it. Cheers Stephan

3 9

VecView behaviour
by Andreas Grassl 10 Jun '09

10 Jun '09

Hello, I'm working with the PCNN preconditioner and hence with ISLocalToGlobalMapping. After solving I want to write the solution to an ASCII-file where only the values belonging to the "external" global numbering are given and not followed by the zeros. Currently I'm giving this commands: ierr = PetscViewerSetFormat(viewer,PETSC_VIEWER_ASCII_SYMMODU);CHKERRQ(ierr); ierr = VecView(X,viewer);CHKERRQ(ierr); Does anybody have an idea, which option or function could help me? cheers, ando -- /"\ Grassl Andreas \ / ASCII Ribbon Campaign Uni Innsbruck Institut f. Mathematik X against HTML email Technikerstr. 13 Zi 709 / \ +43 (0)512 507 6091

3 17

MatLUFactor for Complex Matrices - Returns inverse of the diagonal in LU??
by Paul Dostert 27 May '09

27 May '09

I'm a beginner with PETSc, so please forgive me if this is obvious, but I couldn't seem to find any help in the archives. I'm trying to just the hang of the software, so I've been messing around with routines. I'm going to need complex matrices (Maxwell's with PML) so everything is configured for this. I'm messing around with some very simple test cases, and have a symmetric (but not Hermitian) complex matrix with 2-2i on the diagonal and -1+i on the upper and lower diagonals. I am reading this in from a petsc binary file (again, for testing purposes, eventually I'm going to be just reading in my matrix and RHS). I view the matrix, and it has been read in correctly. I perform LU factorization by doing the following (where ISCreateStride(PETSC_COMM_WORLD,m,0,1,&perm); has been called earlier): MatConvert(A,MATSAME,MAT_INITIAL_MATRIX,&LU); MatFactorInfoInitialize(&luinfo); MatLUFactor(LU,perm,perm,&luinfo); I get that LU is: (1,1) 0.2500 + 0.2500i (2,1) -0.5000 (1,2) -1.0000 + 1.0000i (2,2) 0.3333 + 0.3333i (3,2) -0.6667 (2,3) -1.0000 + 1.0000i (3,3) 0.3750 + 0.3750i (4,3) -0.7500 (3,4) -1.0000 + 1.0000i (4,4) 0.4000 + 0.4000i (5,4) -0.8000 (4,5) -1.0000 + 1.0000i (5,5) 0.4167 + 0.4167i (6,5) -0.8333 (5,6) -1.0000 + 1.0000i (6,6) 0.4286 + 0.4286i Now, I am interpreting this as L being unit on the diagonal and the lower diagonal portion of this "LU" matrix, while U being the diagona + upper of this "LU" matrix. I can interpret this the other way around as well, and it doesn't matter. However, knowing the LU factorization, it is VERY clear the the proper LU decomposition would have the inverse of the diagonal elements presented here. So I believe I should have LU as: (1,1) 2.0000 - 2.0000i (2,1) -0.5000 (1,2) -1.0000 + 1.0000i (2,2) 1.5000 - 1.5000i (3,2) -0.6667 (2,3) -1.0000 + 1.0000i (3,3) 1.3333 - 1.3333i (4,3) -0.7500 (3,4) -1.0000 + 1.0000i (4,4) 1.2500 - 1.2500i (5,4) -0.8000 (4,5) -1.0000 + 1.0000i (5,5) 1.2000 - 1.2000i (6,5) -0.8333 (5,6) -1.0000 + 1.0000i (6,6) 1.1667 - 1.1667i Is there some reason this returns the inverse of the diagonal entries, or am I completely missing something? Is returning the inverse something standard?? Since I'm new, I'm not quite sure where to look for actual source code. Is there a location where the LU factorization code is written and well commented? Thank you very much!

2 3

looking for opportunity to take training in HPC
by Tomasz Jankowski 27 May '09

27 May '09

I'm looking for opportunity to take training in HPC in one of European Union country. I'm thinking about fifth scheme - "Support for training and career development of researchers" . I don't know, maybe it could be in framework of ITN - Initial Training Networks. But maybe someone know about another possibilities (for example IAPP - Industry-Academia Partnerships and Pathways). I decided to disturb members of this group because I didn't found yet anything interesting in this area. Any help greatly appreciated. thanks, tomasz jankowski skills/capabilities: - c/c++/perl (multithrading with pthread,boost::thread,MPI) - I have hands-on computing experience with building sparse systems of equations with C/C++ code from raw data and solving them with PETSC on PC's cluster set up by me. - & much,much more ;-) About me: After I finished my Ph.D. studies in quantitative genetics field in 2007 (defended in April 2008) I have switched from science to IT. (It was possible because I did a lot of programming to obtain the results in my Ph.D.). Actually I'm system programmer in supercomputing and networking center which serves for the needs (especially HPC) of the European scientific community and provides advanced services for the benefit of regional administration. Because of this two field career I'm fluent in quantitative genetics, molecular genetics, genomics, biotechnology, c/c++ programming, linux os administration and also practiced in statistics. ######################################################## # tomjan(a)jay.au.poznan.pl # # jay.au.poznan.pl/~tomjan/ # ########################################################

1 0

prometheus related question,
by Ryan Yan 21 May '09

21 May '09

Hi all, I did a numerical experiment of ksp solve on a large sparse matrix. For the first step, I load the matrix from petsc readable files. If I load the matrix as a MPIAIJ matrix and then apply the prometheus preconditioner with GMRES, it will not converge in 10000 iterations. However, if I load the matrix as a MPIBAIJ matrix and then apply the prometheus preconditioner with GEMRES, it will converge with in 800 iterations Can anyone please provide some hints on how to explain this numerical result. Thanks a lot, Yan

2 1

Re: About the -pc_type tfs
by Ryan Yan 21 May '09

21 May '09

Hi Matt, May I make an inquiry about the file that I send out last week. If there is any thing missing, please let me know. I have other choices on the -pc_type, but it seems like that the 'tfs' is very promising, given it's excellent convergence trajectory when it was applied onto the smaller matrix. Thank you very much, Yan On Fri, May 15, 2009 at 10:46 PM, Ryan Yan <vyan2000(a)gmail.com> wrote: > Hi Matt, > Since the bin file is quit big, I did not post this email to the group. You > can reply the next response to the group. > > The matrix bin file is in the attachment. > petsc_matrix_coef.bin is the bin file to generate the matrix, which have > the format aij. > > petsc_vec_knownsolu.bin is the bin file to generate the exact solution x. > > > petsc_vec_rhs.bin is the bin file to generate the right hand side b. > > A c script named "tfs_binarymatrix_verify.c" is also attached for your > convenience to check Ax-b =0. > > If you need any more information, please let me know. > > Thank you very much, > > Yan > > > > > > On Fri, May 15, 2009 at 1:34 PM, Matthew Knepley <knepley(a)gmail.com>wrote: > >> If you send the matrix in PETSc binary format we can check this. >> >> Matt >> >> >> On Fri, May 15, 2009 at 12:20 PM, Ryan Yan <vyan2000(a)gmail.com> wrote: >> >>> Hi all, >>> I am tring to use the tfs preconditioner to solve a large sparse mpiaij >>> matrix. >>> >>> 11111111111111111111111111111111111111111 >>> It works very well with a small matrix 45*45(Actually a 9*9 block matrix >>> with blocksize 5) on 2 processors; Out put is as follows: >>> >>> 0 KSP preconditioned resid norm 3.014544557924e+04 true resid norm >>> 2.219812091849e+04 ||Ae||/||Ax|| 1.000000000000e+00 >>> 1 KSP preconditioned resid norm 3.679021546908e-03 true resid norm >>> 1.502747104104e-03 ||Ae||/||Ax|| 6.769704109737e-08 >>> 2 KSP preconditioned resid norm 2.331909907779e-09 true resid norm >>> 8.737892755044e-10 ||Ae||/||Ax|| 3.936320910733e-14 >>> KSP Object: >>> type: gmres >>> GMRES: restart=30, using Classical (unmodified) Gram-Schmidt >>> Orthogonalization with no iterative refinement >>> GMRES: happy breakdown tolerance 1e-30 >>> maximum iterations=10000, initial guess is zero >>> tolerances: relative=1e-10, absolute=1e-50, divergence=10000 >>> left preconditioning >>> PC Object: >>> type: tfs >>> linear system matrix = precond matrix: >>> Matrix Object: >>> type=mpiaij, rows=45, cols=45 >>> total: nonzeros=825, allocated nonzeros=1350 >>> using I-node (on process 0) routines: found 5 nodes, limit used is >>> 5 >>> Norm of error 2.33234e-09, Iterations 2 >>> >>> 2222222222222222222222222222222222222222 >>> >>> However, when I use the same code for a larger sparse matrix, a 18656 * >>> 18656 block matrix with blocksize 5); it encounters the followins >>> error.(Same error message for using 1 and 2 processors, seperately) >>> >>> [0]PETSC ERROR: >>> ------------------------------------------------------------------------ >>> [0]PETSC ERROR: Caught signal number 11 SEGV: Segmentation Violation, >>> probably memory access out of range >>> [0]PETSC ERROR: Try option -start_in_debugger or >>> -on_error_attach_debugger >>> [0]PETSC ERROR: or see >>> http://www.mcs.anl.gov/petsc/petsc-as/documentation/troubleshooting.html#Si…<http://www.mcs.anl.gov/petsc/petsc-as/documentation/troubleshooting.html#Si…>ERROR: or try >>> http://valgrind.org on linux or man libgmalloc on Apple to find memory >>> corruption errors >>> [0]PETSC ERROR: likely location of problem given in stack below >>> [0]PETSC ERROR: --------------------- Stack Frames >>> ------------------------------------ >>> [0]PETSC ERROR: Note: The EXACT line numbers in the stack are not >>> available, >>> [0]PETSC ERROR: INSTEAD the line number of the start of the >>> function >>> [0]PETSC ERROR: is given. >>> [0]PETSC ERROR: [0] PCSetUp_TFS line 116 src/ksp/pc/impls/tfs/tfs.c >>> [0]PETSC ERROR: [0] PCSetUp line 764 src/ksp/pc/interface/precon.c >>> [0]PETSC ERROR: [0] KSPSetUp line 183 src/ksp/ksp/interface/itfunc.c >>> [0]PETSC ERROR: [0] KSPSolve line 305 src/ksp/ksp/interface/itfunc.c >>> [0]PETSC ERROR: --------------------- Error Message >>> ------------------------------------ >>> [0]PETSC ERROR: Signal received! >>> [0]PETSC ERROR: >>> ------------------------------------------------------------------------ >>> [0]PETSC ERROR: Petsc Release Version 2.3.3, Patch 15, Tue Sep 23 >>> 10:02:49 CDT 2008 HG revision: 31306062cd1a6f6a2496fccb4878f485c9b91760 >>> [0]PETSC ERROR: See docs/changes/index.html for recent updates. >>> [0]PETSC ERROR: See docs/faq.html for hints about trouble shooting. >>> [0]PETSC ERROR: See docs/index.html for manual pages. >>> [0]PETSC ERROR: >>> ------------------------------------------------------------------------ >>> [0]PETSC ERROR: ./kspex1reader_binmpiaij on a linux-gnu named >>> vyan2000-linux by vyan2000 Fri May 15 01:06:12 2009 >>> [0]PETSC ERROR: Libraries linked from >>> /home/vyan2000/local/PPETSc/petsc-2.3.3-p15//lib/linux-gnu-c-debug >>> [0]PETSC ERROR: Configure run at Mon May 4 00:59:41 2009 >>> [0]PETSC ERROR: Configure options >>> --with-mpi-dir=/home/vyan2000/local/mpich2-1.0.8p1/ --with-debugger=gdb >>> --with-shared=0 --download-hypre=1 --download-parmetis=1 >>> [0]PETSC ERROR: >>> ------------------------------------------------------------------------ >>> [0]PETSC ERROR: User provided function() line 0 in unknown directory >>> unknown file >>> application called MPI_Abort(MPI_COMM_WORLD, 59) - process 0[cli_0]: >>> aborting job: >>> application called MPI_Abort(MPI_COMM_WORLD, 59) - process 0 >>> >>> >>> 3333333333333333333333333333333333333333333333 >>> >>> I have the exact solution x in hands, so before I push the matrix into >>> the ksp solver, I did check the PETSC loaded matrix A and rhs vector b, by >>> verifying Ax-b=0, in both cases of 1 processor and 2 processors. >>> >>> Any sugeestions? >>> >>> Thank you very much, >>> >>> Yan >>> >> >> >> >> -- >> 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 >> > >

1 0

Parallelization of KSPSolve() in multidimensions
by Rolf Kuiper 20 May '09

20 May '09

Hi PETSc-Users, when solving an implicit equation with KSPSolve() in 3D (communication with 7-point-stencil) I experienced the following: Parallelization of the e.g. 64 x 64 x 64 domain on n cpus in the last direction (every cpu has a 64 x 64 x 64/n subdomain) leads to a parallel efficiency of approximately 90%, which is fine for us. Parallelization of the e.g. 64 x 64 x 64 domain on n cpus in more than one direction (every cpu has e.g. a 64 x 64/sqrt(n) x 64/sqrt(n) subdomain) leads to a parallel efficiency of approximately 10%, which is absolutely unusable. Is this behavior generally true for this kind of solver? If so, why? If not: What did I do wrong most presumably? Has anybody made the same experience and/or could help me? Thanks in advance, Rolf

2 1

Time for MatAssembly
by tribur＠vision.ee.ethz.ch 19 May '09

19 May '09

Distinguished PETSc experts, Assuming processor k has defined N entries of a parallel matrix using MatSetValues. The half of the entries are in matrix rows belonging to this processor, but the other half are situated in rows of other processors. My question: When does MatAssemblyBegin+MatAssemblyEnd take longer, if the rows where the second half of the entries are situated belong all to one single other processor, e.g. processor k+1, or if these rows are distributed across several, let's say 4, other processors? Is there a significant difference? Thanks in advance for your answer, Trini

3 3