October 2010 - petsc-users - mailman3.cels.anl.gov

CHKERRQ in Fortran
by Gaetan Kenway 06 Oct '10

06 Oct '10

Hello I use PETSc with fortran. I was wondering if the CHKERRQ(ierr) command is supposed to work in Fortran? My compiler (mpif90 with ifort). If I do something like this: call VecCreate(WARP_COMM_WORLD,globalSurfForce,ierr) CHKERRQ(ierr) ifort complains there is a syntax error. I also tried: call VecCreate(WARP_COMM_WORLD,globalSurfForce,ierr) call CHKERRQ(ierr) But then it complains that it can't find the chkerrq function while linking. I'm using PETSc-3.1-p3 which was compiled with the following options: --with-shared --download-superlu_dist=yes --download-spooles=yes --download-parmetis=yes --with-fortran-interfaces=1 The subroutine includes: #include "include/finclude/petsc.h" Am I missing something obvious? Gaetan

3 2

petsc and meshless type method: Forming parallel coefficient matrix
by M. Scot Breitenfeld 05 Oct '10

05 Oct '10

Hi, I'm working on implementing petsc in a meshfree type method. I have the serial version working with Petsc, but I have some questions about the parallel implementation. I guess the easiest way to explain the problem is with an example, lets take the 1D problem, 1 2 3 4 5 6 7 8 o o o o | o o o o Nodal partitioning: proc 0: 1-4 proc 1: 5-8 Now, Proc 0 contributes values to 5 and 6, so they are included in my "local numbering" and the total number of local nodes is 6 on proc 0 Proc 1 contributes values to 3 and 4, so they are included in my "local numbering" and the total number of local nodes is 6 on proc 1 1 2 3 4 5 6 o o o o o o Proc 0 o o o o o o Proc 1 3 4 5 6 7 8 Each node has 3 dof, so the global coefficient matrix, A, would be 24x24: Processor 0 has rows (1-12) and Processor 1 has rows (13-24) When each processor loops over its nodes (Proc 0: 1-6, Proc 1: 3-4) it adds its contribution into A: CALL MatSetValues(A, 1, "global row id", "number of column entries", "column entries", "values", ADD_VALUES,ierr) I'm unclear what the best way to create the A matrix is. For each processor I have its number of global nodes (proc 0: 1-6, proc 1: 3-8) so I can get the corresponding rows in the global A that it contributes (note that some of the rows are not local to the processor, those values need to be sent). Also, for each processor I have a list of the actual global nodes the processors owns (proc 0: 1-4, proc 1:5-8). I can set the mapping with: CALL ISLocalToGlobalMappingCreate(PETSC_COMM_WORLD, "number of local elements", "global index of elements", mapping, ierr) where I'm assuming the "number of local elements" includes the ghost nodes. This is essentially what global row this processor will contribute. How do I tell Petsc what nodes the processor actually owns so it knows what nodes need to be sent? Is there another mapping command to use for mapping the owned global nodes to the processor and then is there a command that takes both maps and figures out what nodes need to be sent. Or is the correct procedure to use VecScatterCreate and create these arrays myself. Can I still use MatSetValues or do I have to use MatSetValuesLoc (together with VecSetLocalToGlobalMapping) to assemble the matrix, or can I still do this all using global numbering. If there is a similar example to what I'm doing that you could point me to it would be helpful. Thanks for your time, Scot

4 5

preallocation for FE matrix
by Tabrez Ali 04 Oct '10

04 Oct '10

I am trying to assemble a FE matrix (adding one element at a time from the local stiffness matrix) and cant seem to get the preallocation right. I am correctly calculating the number of non zeros per row and storing the value in the array [nnzpr]. Here is part of the relevant code ... call MatCreateSeqBAIJ(Petsc_Comm_Self, 1, m, n, petsc_null_integer, nnzpr, Mat_A, ierr) ... ! Assume a bilinear quad (2 dof per node) do j1=1,8 do j2=1,8 call MatSetValues(Mat_A, 1, indx(j1)-1, 1, indx(j2)-1, k(j1,j2), Add_Values, ierr) end do end do ... call MatAssemblyBegin(Mat_A,Mat_Final_Assembly,ierr) call MatAssemblyEnd(Mat_A,Mat_Final_Assembly,ierr) ... On running it I get -bash-3.00$ ./a.out <inp -info -mat_view_info ... Total non-zero elements estimated i.e., sum(nnzpr) = 112 [0] PetscCommDuplicate(): Duplicating a communicator 1140850689 -2080374784 max tags = 2147483647 [0] PetscCommDuplicate(): returning tag 2147483647 [0] MatAssemblyEnd_SeqBAIJ(): Matrix size: 12 X 12, block size 1; storage space: 42 unneeded, 112 used [0] MatAssemblyEnd_SeqBAIJ(): Number of mallocs during MatSetValues is 9 [0] MatAssemblyEnd_SeqBAIJ(): Most nonzeros blocks in any row is 12 ... Matrix Object: type=seqbaij, rows=12, cols=12 total: nonzeros=112, allocated nonzeros=154 block size is 1 I dont understand why -mat_view_info shows allocated nonzeros to be 154 when sum(nnzpr) is 112. Thanks in advance.

2 4

Project: Matlab binding for PETSc
by Barry Smith 03 Oct '10

03 Oct '10

PETSc-users, I've started a little project to develop Matlab bindings for PETSc. It is not difficult but requires manually translating the API so is a bit tedious. If anyone is interesting in helping out the code is bin/matlab/classes; I have only the most rudimentary parts done and it is not ready for users but with help it could be. It you are interested then get petsc-dev http://www.mcs.anl.gov/petsc/petsc-as/developers/index.html and join petsc-dev(a)mcs.anl.gov http://www.mcs.anl.gov/petsc/petsc-as/miscellaneous/mailing-lists.html and take a look at the comments in bin/matlab/classes/PetscInitialize.m We need to add the rest of the bindings to matlabheader.h and all the classes. Plus, of course, simple tests. Happy computing, Barry

1 0

petsc4py Vec: how to copy from one index to multiple indices of the same vector?
by Amal Alghamdi 03 Oct '10

03 Oct '10

Dear all, I would like to ask about the proper way to copy a value in vector x, let us say x[0] into places in the same vector, let us say x[1], x[100]. I want to do this using petsc4py. Actually I have tried the method getValue in Vec class. But I received an error telling me that this method "Can only get local values", so I tried to follow example 10 in the manual, which is posted below, but when trying to create sequential vector I get the error " [1] Invalid argument [1] Cannot create VECSEQ on more than one process " Vec p, x; /* initial vector, destination vector */ VecScatter scatter; /* scatter context */ IS from, to; /* index sets that define the scatter */ PetscScalar *values; int idx_from[] = {100,200}, idx_to[] = {0,1}; VecCreateSeq(PETSC COMM SELF,2,&x); ISCreateGeneral(PETSC COMM SELF,2,idx from,&from); ISCreateGeneral(PETSC COMM SELF,2,idx to,&to); VecScatterCreate(p,from,x,to,&scatter); VecScatterBegin(scatter,p,x,INSERT VALUES,SCATTER FORWARD); VecScatterEnd(scatter,p,x,INSERT VALUES,SCATTER FORWARD); VecGetArray(x,&values); ISDestroy(from); ISDestroy(to); VecScatterDestroy(scatter); Thank you very much Amal

2 2

does petsc support schur complement
by Gong Ding 02 Oct '10

02 Oct '10

Dear all, I had take a look with the manual of petsc.It has refered schur complement. However, I can not find more useful information. Anyway, I am new to this method and I would like to get any information about it. Yours Gong Ding

3 2

[Fortran] subroutines inside modules?
by Leo van Kampenhout 01 Oct '10

01 Oct '10

Hi all, since it is mandatory to declare all subroutines as "external" in Fortran, is it possible for Modules to have subroutines? I'm unable to declare the subroutine external inside the module itself, nor in the program which is using it. Not declaring it external at all results in the following compilation error: /net/users/csg/csg4035/master/workdir/src/main.F:97: undefined reference to `__grid_MOD_readgrid' (the module is here is named "grid", the subroutine "readgrid" ) Thanks, Leo

3 5