[petsc-users] unexpected ordering when VecSetValues set multiples values

Matteo Parsani parsani.matteo at gmail.com
Fri Jun 7 14:03:37 CDT 2013


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 <bsmith at mcs.anl.gov> 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 <parsani.matteo at gmail.com>
> 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 <bsmith at mcs.anl.gov> wrote:
> >
> > On Apr 30, 2013, at 7:21 AM, Matteo Parsani <parsani.matteo at gmail.com>
> 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 <bsmith at mcs.anl.gov>
> 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 <parsani.matteo at gmail.com>
> 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 <bsmith at mcs.anl.gov>
> wrote:
> > > >
> > > > On Apr 26, 2013, at 10:20 AM, Matteo Parsani <
> parsani.matteo at gmail.com> 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
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20130607/515e6964/attachment.html>


More information about the petsc-users mailing list