Hello Barry, I have tried to pass to VectSetValues two values *call VecSetValues(x_vec_in,2,index_list(1:2),x_vec_local(1:2),INSERT_VALUES,ierr_local) *
By looking at the values of index_list and x_vec_local with gdb I see the correct numbers: (gdb) print x_vec_local(1) $1 = 0.9394877605834675 (gdb) print x_vec_local(2) $2 = 0.55294328167143081 (gdb) print index_list(1) $3 = 0 (gdb) print index_list(2) $4 = 1 Now, when I ask gdb to print the values of ni, ix and y in *PetscErrorCode VecSetValues_Seq(Vec xin,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode m)* I get (gdb) print ni $5 = 2 (gdb) print ix[0] $6 = 0 (gdb) print ix[1] $7 = 0 print y[0] $13 = 0.9394877605834675 (gdb) print y[1] $14 = 0.55294328167143081 The values of y match with the values I am passing in, i.e. x_vec_local. However the index ix is wrong because the first two components of it are zero. Thus, when i = 1 in the VecSetValues_Seq for loop, the first element that was placed correctly is replaced by the second element. This is also waht I see if I use call VecView(x_vec_in,PETSC_VIEWER_STDOUT_WORLD,ierr_local) Now if I keep printing ix in gdb I see the following pattern: (gdb) print ix[0] $18 = 0 (gdb) print ix[1] $19 = 0 (gdb) print ix[2] $20 = 1 (gdb) print ix[3] $21 = 0 (gdb) print ix[4] $22 = 2 (gdb) print ix[5] $23 = 0 (gdb) print ix[6] $24 = 3 which I think is strange since I ix should be a vector of length 2 because I am calling VecSetValues as *call VecSetValues(x_vec_in,2,index_list(1:2),x_vec_local(1:2),INSERT_VALUES,ierr_local) * Am I doing something wrong? Thank you in advance. On Fri, Jun 7, 2013 at 3:03 PM, Matteo Parsani <[email protected]>wrote: > Yes, they are correct. When asked, the debugger prints the right values > for bot the vector of the index and the vector of the values. > > Sorry, how can I set a break point in VecSetValues_Seq() and > VecSetValues_MPI()? > If I type > > break VecSetValues_Seq() > > or > > break VecSetValues_MPI() > > or > > /home0/pmatteo/research/lib_src/petsc/src/vec/vec/impls/seq/bvec2.c:1050 > > I get > > Function "vecsetvalues_seq" not defined. > > or > > Function "vecsetvalues_mpi" not defined. > > > Thanks > > > > On Fri, Jun 7, 2013 at 2:32 PM, Barry Smith <[email protected]> wrote: > >> >> Are, seemingly the correct values getting into the vecsetvalues_ call? >> >> Put a break point also in VecSetValues_Seq() and VecSetValues_MPI() >> when it gets into either of those routines you can do step in the debugger >> and watch it take the values passed in and attempt to put them in the >> appropriate places in the vector's array. >> >> Barry >> >> On Jun 7, 2013, at 1:20 PM, Matteo Parsani <[email protected]> >> wrote: >> >> > Hello Barry, >> > I have tried to follow the debugger but I got lost. >> > >> > I have run my code in this way: >> > >> > ./NSE -start_in_debugger noxterm >> > >> > then after getting some output I get >> > >> > gdb() >> > >> > Then I type break vecsetvalues_ >> > >> > >> > and a break point is generated >> > >> > Breakpoint 1 at 0x2ae97eda5827: file >> /scratch/home0/pmatteo/research/lib_src/petsc/src/vec/vec/interface/ftn-auto/rvectorf.c, >> line 239. >> > >> > >> > I keep typing n and it goes on. It shows me every time vectsetvalues is >> called but it never goes inside vectsetvalues. >> > >> > I have tested the VecSetValues with just two elements and it does not >> work. It sets correctly the values both in serial and in parallel if I set >> one value at the time. >> > >> > Thanks >> > >> > >> > >> > On Tue, Apr 30, 2013 at 8:35 AM, Barry Smith <[email protected]> >> wrote: >> > >> > On Apr 30, 2013, at 7:21 AM, Matteo Parsani <[email protected]> >> wrote: >> > >> > > Hello Barry, >> > > it does not work also if m is not so large. >> > >> > Then you need to run in the debugger. You can use the command line >> options >> > >> > -start_in_debugger noxterm >> > >> > in the debugger put a break point in vecsetvalues_ (yes lower >> case followed by an underscore) then when it >> > gets to that point print the size of the array and the values in the >> integer and double array to see if they match, then step into the >> VecSetValues() function that is called by vecsetvalues_() and look at the >> array values ago, continue to step along and it will go into >> VecSetValues_Seq() and start actually putting the values into the Vec array. >> > >> > Barry >> > >> > > >> > > Thanks again. >> > > >> > > >> > > >> > > >> > > >> > > >> > > On Fri, Apr 26, 2013 at 3:23 PM, Barry Smith <[email protected]> >> wrote: >> > > >> > > Shouldn't matter that it is called from Fortran we do it all the >> time. >> > > >> > > Does it work if the final m is not very large? >> > > >> > > You may need to run in the debugger and follow the values through >> the code to see why they don't get to where they belong. >> > > >> > > >> > > Barry >> > > >> > > Could also be a buggy fortran compiler. >> > > >> > > >> > > On Apr 26, 2013, at 2:06 PM, Matteo Parsani <[email protected]> >> wrote: >> > > >> > > > Hello Barry, >> > > > sorry I modified few things just before to write the mail. >> > > > The correct loop with the correct name of the variables is the >> following >> > > > >> > > > ! Number of DOFs owned by this process >> > > > len_local = (elem_high-elem_low+1)*nodesperelem*nequations >> > > > >> > > > ! Allocate memory for x_vec_local and index_list >> > > > allocate(x_vec_local(len_local)) >> > > > allocate(index_list(len_local)) >> > > > >> > > > >> > > > m = 1 >> > > > ! Loop over all elements >> > > > do ielem = elem_low, elem_high >> > > > ! Loop over all nodes in the element >> > > > do inode = 1, nodesperelem >> > > > !Loop over all equations >> > > > do ieq = 1, nequations >> > > > ! Add element to x_vec_local >> > > > x_vec_local(m) = ug(ieq,inode,ielem) >> > > > ! Add element to index list >> > > > index_list(m) = >> (elem_low-1)*nodesperelem*nequations+m-1 >> > > > ! Update m index >> > > > m = m+1 >> > > > end do >> > > > end do >> > > > end do >> > > > >> > > > ! HERE I HAVE PRINTED x_vec_local, ug and index_list >> > > > >> > > > ! Set values in the portion of the vector owned by the >> process >> > > > call >> VecSetValues(x_vec_in,m-1,index_list,x_vec_local,INSERT_VALUES,& >> > > > & ierr_local) >> > > > >> > > > ! Assemble initial guess >> > > > call VecAssemblyBegin(x_vec_in,ierr_local) >> > > > call VecAssemblyEnd(x_vec_in,ierr_local) >> > > > >> > > > >> > > > >> > > > I have printed the values and the indices I want to pass to >> VecSetValues() and they are correct. >> > > > >> > > > I also printed the values after VecSetValues() has been called and >> they are wrong. >> > > > >> > > > The attachment shows that. >> > > > >> > > > >> > > > Could it be a problem of VecSetValues() + F90 when more than 1 >> elements is set? >> > > > >> > > > I order to debug my the code I am running just with 1 processor. >> Thus the process owns all the DOFs. >> > > > >> > > > Thank you. >> > > > >> > > > >> > > > >> > > > >> > > > >> > > > >> > > > >> > > > >> > > > >> > > > >> > > > >> > > > On Fri, Apr 26, 2013 at 2:39 PM, Barry Smith <[email protected]> >> wrote: >> > > > >> > > > On Apr 26, 2013, at 10:20 AM, Matteo Parsani < >> [email protected]> wrote: >> > > > >> > > > > Hello, >> > > > > I have some problem when I try to set multiple values to a PETSc >> vector that I will use later on with SNES. I am using Fortran 90. >> > > > > Here the problem and two fixes that however are not so good for >> performances reasons. The code is very simple. >> > > > > >> > > > > Standard approach that does not work correctly: (I am probably >> doing something wrong) >> > > > > >> > > > > m = 1 >> > > > > ! Loop over all elements >> > > > > do ielem = elem_low, elem_high >> > > > > ! Loop over all nodes in the element >> > > > > do inode = 1, nodesperelem >> > > > > !Loop over all equations >> > > > > do ieq = 1, nequations >> > > > > ! Add element to x_vec_local >> > > > > x_vec_local(m) = ug(ieq,inode,ielem) >> > > > > ! Add element to index list >> > > > > ind(m) = >> (elem_low-1)*nodesperelem*nequations+m-1 >> > > > > ! Update m index >> > > > > m = m+1 >> > > > > end do >> > > > > end do >> > > > > end do >> > > > > >> > > > > ! Set values in the portion of the vector owned by the >> process >> > > > > call >> VecSetValues(x_vec_in,len_local,index_list,x_vec_local,INSERT_VALUES,& >> > > > >> > > > What is len_local and index_list? They do not appear in the >> loop above. Shouldn't you be passing m-1 for the length and ind for the >> indices? >> > > > >> > > > I would first print out the all the values in your input to >> VecSetValues() and make sure they are correct. >> > > > >> > > > Barry >> > > > >> > > > > & ierr_local) >> > > > > >> > > > > ! Assemble initial guess >> > > > > call VecAssemblyBegin(x_vec_in,ierr_local) >> > > > > call VecAssemblyEnd(x_vec_in,ierr_local) >> > > > > >> > > > > Then I print my expected values and the values contained in the >> PETSc vector to a file. See attachment. I am running in serial for the >> moment BUT strangely if you look at the file I have attached the first 79 >> DOFs values have a wrong ordering and the remaining 80 are zero. >> > > > > >> > > > > >> > > > > 1st approach: set just one value at the time inside the loop. >> > > > > m = 1 >> > > > > ! Loop over all elements >> > > > > do ielem = elem_low, elem_high >> > > > > ! Loop over all nodes in the element >> > > > > do inode = 1, nodesperelem >> > > > > !Loop over all equations >> > > > > do ieq = 1, nequations >> > > > > ! Add element to x_vec_local >> > > > > value = ug(ieq,inode,ielem) >> > > > > ! Add element to index list >> > > > > ind = (elem_low-1)*nodesperelem*nequations+m-1 >> > > > > call >> VecSetValues(x_vec_in,1,ind,value,INSERT_VALUES,& >> > > > > & ierr_local) >> > > > > ! Update m index >> > > > > m = m+1 >> > > > > end do >> > > > > end do >> > > > > end do >> > > > > >> > > > > >> > > > > This works fine. As you can see I am using the same expression >> used in the previous loop to compute the index of the element that I have >> to add in the x_vec_in, i.e. >> > > > > ind = (elem_low-1)*nodesperelem*nequations+m-1 >> > > > > >> > > > > Thus I cannot see which is the problem. >> > > > > >> > > > > 2nd approach: get the pointer to the local part of the global >> vector and use it to set the values in the global vector >> > > > > >> > > > > m = 1 >> > > > > ! Loop over all elements >> > > > > do ielem = elem_low, elem_high >> > > > > ! Loop over all nodes in the element >> > > > > do inode = 1, nodesperelem >> > > > > !Loop over all equations >> > > > > do ieq = 1, nequations >> > > > > ! Add element to x_vec_local >> > > > > tmp(m) = ug(ieq,inode,ielem) >> > > > > ! Update m index >> > > > > m = m+1 >> > > > > end do >> > > > > end do >> > > > > end do >> > > > > >> > > > > >> > > > > This works fine too. >> > > > > >> > > > > >> > > > > Jut to be complete. I use the following two approaches to view >> the vector: >> > > > > >> > > > > call VecView(x_vec_in,PETSC_VIEWER_STDOUT_WORLD,ierr_local) >> > > > > >> > > > > >> > > > > and >> > > > > >> > > > > call VecGetArrayF90(x_vec_in,tmp,ierr_local) >> > > > > >> > > > > >> > > > > m = 1 >> > > > > ! Loop over all elements >> > > > > do ielem = elem_low, elem_high >> > > > > ! Loop over all nodes in the element >> > > > > do inode = 1, nodesperelem >> > > > > !Loop over all equations >> > > > > do ieq = 1, nequations >> > > > > write(*,*) >> m,index_list(m),x_vec_local(m),tmp(m) >> > > > > ! Update m index >> > > > > m = m+1 >> > > > > end do >> > > > > end do >> > > > > end do >> > > > > >> > > > > >> > > > > Thank you. >> > > > > >> > > > > >> > > > > -- >> > > > > Matteo >> > > > > <values.txt> >> > > > >> > > > >> > > > >> > > > >> > > > -- >> > > > Matteo >> > > > <comparison.txt> >> > > >> > > >> > > >> > > >> > > -- >> > > Matteo >> > >> > >> > >> > >> > -- >> > Matteo >> >> > > > -- > Matteo > -- Matteo
