> On Aug 9, 2020, at 1:14 AM, Nidish <[email protected]> wrote:
> 
> Hello,
> 
> My question is related to solving a generic nonlinear algebraic system in 
> parallel. I'm using the word generic to point to the fact that the system may 
> or may not be from a finite element/mesh type object - it's just a residual 
> function and a Jacobian function that I can provide.
> 
> I followed the serial example in $PETSC_DIR/src/snes/tutorials/ex1.c pretty 
> much verbatim (albeit with a slightly different objective function whose size 
> can be controlled), and have been trying to extend it to a parallel 
> application (the given example is strictly serial).
> 
> When I run my code in serial everything works well and I get the solution I 
> expect. However, even if I use just two cores, the code starts returning 
> garbage. I had a feeling it has something to do with index ordering, but I 
> want to get this code working without using a DM object. I'm attaching my 
> code with this email.
> 
> TL;DR Version: What modifications must be made in 
> "$PETSC_DIR/src/snes/tutorials/ex1.c" to make it work in Parallel? Is it 
> possible to make minimal modifications to make it run in parallel too?

  This is wrong:

  VecGetSize(x, &N);
  VecGetOwnershipRange(x, &is, &ie);
  /* VecGetArrayRead(x, &xx); */
  VecGetArray(x, &xx);  
  for (int i=is; i<ie; ++i)

  You appear to be trying to mimic a DM example but without a DM.

   VecGetArray() always returns an array with indices from 0 to 
VecGetLocalSize()

   Thus you need to make your loops over local indices, not global indices.

   In addition you are likely to need "ghost" values of the input vectors. If 
you are not using DM then you need to have code that determines what ghost 
values are needed and create a VecScatter to manage communicating those ghost 
values. Look for examples that use VecScatterBegin in the formfunction routine 
for how this may be done.

  Barry



> 
> Thank You,
> Nidish
> <sness.c>

Reply via email to