SNES Problem

Nils Erik Svangård nilserik at gmail.com
Thu Mar 2 07:29:35 CST 2006


Hi again,
I realised that my implementation was a bit flawed so I have rewritten
it again, the push and pop() functions dont seem to work, and they are
keeping me from compiling the program, this is the main timestep loop:

      call VecCreateSeq(PETSC_COMM_SELF,N,x,ierr)
      call VecDuplicate(x,r,ierr)
      call VecDuplicate(x,xb,ierr)
      call SNESCreate(PETSC_COMM_SELF,snes,ierr)
      call SNESSetFunction(snes,r,FormFunction,PETSC_NULL_OBJECT,ierr)
      call SNESSetFromOptions(snes,ierr)
      call FormInitialGuess(x,ierr)
      call VecCopy(x,xb,ierr)
      write(6,*)zero,dt,tmax,i1000
      step=0
      timestep=dt
C While loop
   10 if (step .le. 100) then
        call PetscExceptionPush(87,ierr)
        call SNESSolve(snes,PETSC_NULL_OBJECT,x,ierr)
        call PetscExceptionPop(87,ierr2)

        call SNESGetConvergedReason(snes,reason,ierr)
        IF ( reason .LT. 0 ) THEN
           timestep=timestep/2
           VecCopy(xb,x,ierr)
           write(6,*)"Solution diverged, trying with",
     .               " smaller timestep: ",timestep
        ELSE
         IF ( ierr .EQ. 87) THEN
           timestep=timestep/2
           VecCopy(xb,x,ierr)
           write(6,*) "No physical solution, trying",
     .                " with smaller timestep: ",timestep
         ELSE
           timestep=dt
           step=step+1
           VecCopy(x,xb,ierr)
           write(6,*) "Solve success, going for step: ",step
         ENDIF
        ENDIF
       goto 10
      endif
C End While loop

And in my Formfunction I have added the following check:
if (ET(L) .LT. 0) then
                ierr = 87
                return
              endif

PETSC_ERR_ARG_DOMAIN didnt work:  error on line 6500 of SNESG3D2.F:
data type is undefined for variable PETSC_ERR_ARG_DOMAIN
Even thought I had petsc.h included.

Can I do this check when I use TS also? I have the sameproblem with
the values being negative when using TS.


/nisse


On 2/27/06, Barry Smith <bsmith at mcs.anl.gov> wrote:
>
>    SNES works by computing p = -approxinv(J)*F(uold) and
> then does a line search on unew = uold + lambda*p to get the
> new u. First it uses a test value of 1 for lambda so it
> tries to compute F(uold + p). It is possible that uold + p
> has some "non-physical" values in it.  There are two ways
> you can try handling it:
>
> 1) Before doing the linesearch SNES calls a "precheck" function,
> that can change the step if it decides there is a problem
> with the step (like it is too long). You can provide your
> own precheck function with SNESLineSearchSetPreCheck() it could,
> for example, shrink the Newton direction to make it remain physical.
>
> or 2) when your form function detects an illegal value it calls
> SETERRQ(PETSC_ERR_ARG_DOMAIN,"Nonphysical function input"); Note:
> this is before it takes the square root of the negative number,
> so check the number and call the error before calling SETERRQ().
> (If using Fortran then simply set ierr = PETSC_ERR_ARG_DOMAIN and
> then immediately return). Then replace your call to SNESSolve with
>    ierr = PetscExceptionTry1(SNESSolve(snes,b,x),PETSC_ERR_ARG_DOMAIN);
>    if (PetscExceptionCaught(ierr,PETSC_ERR_ARG_DOMAIN)) {
>      /* this means your function found a non-physical value so
>         cut your time-step and continue through the loop again. */
>         Put code to do this here.
>    }
> If using fortran then do
> call PetscExceptionPush(PETSC_ERR_ARG_DOMAIN,ierr)
> call SNESSolve(snes,b,x,ierr)
> call PetscExceptionPop(PETSC_ERR_ARG_DOMAIN,anotherierr)
> if (ierr == PETSC_ERR_ARG_DOMAIN) then
>    non-physical value so cut the time-step and try again
> else ! everything is normal so take the next time-step
>
>    Barry
>
> The Fortran interface may be missing PetscExceptionPush() and Pop()
> if so let us know and we'll provide the patch.
>
> On Mon, 27 Feb 2006, Nils Erik Svangård wrote:
>
> > Hi all,
> >
> > I have problems solving some CFD problems using SNES and my custom
> > back-euler. I have 7 equations that I want to solve.
> >
> > First my FormFunction copies the values from the Vec that SNES uses to
> > the variables that my code use:
> >
> >      RO(L)=xx(1,L)
> >      RU(L)=xx(2,L)
> >      RV(L)=xx(3,L)
> >      RW(L)=xx(4,L)
> >      ET(L)=xx(5,L)
> >      RQ(L)= xx(6,L)
> >      REPS(L)=xx(7,L)
> >
> > Here I also print the values of ET(2) for debugging:
> >
> >      write(6,*)"ET(2) = xx(5,2) : ",ET(2)," = ",xx(5,2)
> >
> > Then I get the fluxes by calling custom functions
> >
> >      call AUXVR
> >      call VGRAD
> >      call FLUX
> >      call KESRC
> >
> > Then I perform back-euler  save the new value to so that I can use it
> > in the next iterationi (TSF(L) is Time step function, which is set by
> > hand and is the same for all L):
> >
> >      ff(1,L) = RO(L)-OLD(1,L)-TSF(L)*DRO(L)
> >      ff(2,L) = RU(L)-OLD(2,L)-TSF(L)*DRU(L)
> >      ff(3,L) = RV(L)-OLD(3,L)-TSF(L)*DRV(L)
> >      ff(4,L) = RW(L)-OLD(4,L)-TSF(L)*DRW(L)
> >      ff(5,L) = ET(L)-OLD(5,L)-TSF(L)*DET(L)
> >      ff(6,L) = RQ(L)-OLD(6,L)-TSF(L)*DRQ(L)
> >      ff(7,L) = REPS(L)-OLD(7,L)-TSF(L)*DREPS(L)
> >
> > And save the new value of RO-REPS to use in the next iteration of back-euler:
> >
> >      OLD(1,L)=RO(L)
> >      OLD(2,L)=RU(L)
> >      OLD(3,L)=RV(L)
> >      OLD(4,L)=RW(L)
> >      OLD(5,L)=ET(L)
> >      OLD(6,L)=RQ(L)
> >      OLD(7,L)=REPS(L)
> >
> > Here I print the values of L=2 as above for debugging:
> >
> >      write(6,*)ff(5,2)," = ",ET(2),"-",OLD(5,2),"-",TSF(2),"*",DET(2)
> >
> >
> >
> > The program exit abnormaly after 33 runs of the FormFunction. The
> > cause of this is that AUXVR tries to perform sqrt(ET(2)) when ET(2) is
> > negative, ET is the total energy and should never be negative. It
> > seems that all of a sudden the PETSc SNES solver decides to supply the
> > FormFunction with a negative ET value. I need help understanding why
> > and how to fix it (if it is fixable).
> >
> > Here is the output of:
> >      write(6,*)"ET(2) = xx(5,2) : ",ET(2)," = ",xx(5,2)
> > Which is in the beginning of my FormFunction.
> > And:
> >      write(6,*)ff(5,2)," = ",ET(2),"-",OLD(5,2),"-",TSF(2),"*",DET(2)
> > Which is in the end of my Formfunction.
> >
> > 1. ET(2) = xx(5,2) :   253250.  =   253250.000000000
> >   ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4126563437192  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.941266E+08
> > 2. ET(2) = xx(5,2) :   253250.  =   253250.000010259
> >   ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4126563437192  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.941266E+08
> > 3. ET(2) = xx(5,2) :   253250.  =   253249.999999988
> >   ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4132419437260  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.941324E+08
> > 4. ET(2) = xx(5,2) :   253250.  =   253250.000048432
> >   ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4126563437192  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.941266E+08
> > 5. ET(2) = xx(5,2) :   253250.  =   253249.999999872
> >   ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4177283437784  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.941773E+08
> > 6. ET(2) = xx(5,2) :   253250.  =   253250.000071113
> >   ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4126563437192  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.941266E+08
> > 7. ET(2) = xx(5,2) :   253250.  =   253249.999980115
> >   ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4453347441010  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.944533E+08
> > 8. ET(2) = xx(5,2) :   253250.  =   253250.000037298
> >   ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4082723436679  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.940827E+08
> > 9. ET(2) = xx(5,2) :   253250.  =   253249.999942706
> >   ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4082723436679  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.940827E+08
> > 10. ET(2) = xx(5,2) :   253250.  =   253249.999948703
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4126563437192  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.941266E+08
> > 11. ET(2) = xx(5,2) :   253250.  =   253250.000007808
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.7283811474088  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.972838E+08
> > 12. ET(2) = xx(5,2) :   253250.  =   253250.000018146
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4126563437192  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.941266E+08
> > 13. ET(2) = xx(5,2) :   253250.  =   253250.000001865
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   28.3742819315846  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.837428E+08
> > 14. ET(2) = xx(5,2) :   253250.  =   253249.999997296
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4126563437192  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.941266E+08
> > 15. ET(2) = xx(5,2) :   253250.  =   253250.000000049
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   30.6548963582361  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -3.065490E+08
> > 16. ET(2) = xx(5,2) :   253250.  =   253250.000022870
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4126563437192  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.941266E+08
> > 17. ET(2) = xx(5,2) :   253250.  =   253249.999999716
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.5326051451209  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.953260E+08
> > 18. ET(2) = xx(5,2) :   253250.  =   253250.000074328
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4165187437643  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.941652E+08
> > 19. ET(2) = xx(5,2) :   253250.  =   253250.000049809
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4085059436707  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.940851E+08
> > 20. ET(2) = xx(5,2) :   253250.  =   253250.000019253
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4126563437192  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.941266E+08
> > 21. ET(2) = xx(5,2) :   253250.  =   253250.000000189
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.8851075492403  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.988511E+08
> > 22. ET(2) = xx(5,2) :   253250.  =   253249.999997897
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4085059436707  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.940851E+08
> > 23. ET(2) = xx(5,2) :   253250.  =   253249.999975142
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.3961539435263  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.939615E+08
> > 24. ET(2) = xx(5,2) :   253250.  =   253250.000018324
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4085059436707  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.940851E+08
> > 25. ET(2) = xx(5,2) :   253250.  =   253249.999987073
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4165187437643  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.941652E+08
> > 26. ET(2) = xx(5,2) :   253250.  =   253249.999980982
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4165187437643  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.941652E+08
> > 27. ET(2) = xx(5,2) :   253250.  =   253249.999994645
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4126563437192  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.941266E+08
> > 28. ET(2) = xx(5,2) :   253250.  =   253249.999962364
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4165187437643  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.941652E+08
> > 29. ET(2) = xx(5,2) :   253250.  =   253249.999999921
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4126563437192  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.941266E+08
> > 30. ET(2) = xx(5,2) :   253250.  =   253250.000000058
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   28.4019971319085  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.840200E+08
> > 31. ET(2) = xx(5,2) :   253250.  =   253250.000040348
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4085059436707  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.940851E+08
> > 32. ET(2) = xx(5,2) :   253250.  =   253250.000005245
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4126563437192  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.941266E+08
> > 33. ET(2) = xx(5,2) :   253250.  =   253250.000005245
> >    ff(5,2) = ET(2)-OLD(5,2)-TSF(2)*DET(2) :   29.4126563437192  =
> > 253250. -  253250.000000000 -  1.000000E-07 * -2.941266E+08
> > 34. ET(2) = xx(5,2) :  -848141.  =  -848141.388090847
> > --------------------------------------------------------------------------
> > Petsc Release Version 2.3.0, Patch 32, April, 26, 2005
> > See docs/changes/index.html for recent updates.
> > See docs/faq.html for hints about trouble shooting.
> > See docs/index.html for manual pages.
> > -----------------------------------------------------------------------
> > ./SNESG3D2 on a linux-gnu named sethnx004.vac.com by yy26539 Mon Feb
> > 27 14:21:13 2006
> > Libraries linked from /home/yy26539/work/NISSE/petsc-2.3.0/lib/linux-gnu
> > Configure run at Thu Oct 13 08:23:56 2005
> > Configure options --with-cc=gcc --with-fc="f77 -N109"
> > --download-mpich=1 --download-mpich-pm=gforker
> > --download-f-blas-lapack=1 --with-shared=0
> > -----------------------------------------------------------------------
> > [0]PETSC ERROR: Caught signal number 8 FPE: Floating Point
> > Exception,probably divide by zero
> > [0]PETSC ERROR: Try option -start_in_debugger or -on_error_attach_debugger
> > [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] SNES user function line 0 unknownunknown
> > [0]PETSC ERROR: [0] SNESComputeFunction line 788 src/snes/interface/snes.c
> > [0]PETSC ERROR: [0] SNESLineSearchCubic line 514 src/snes/impls/ls/ls.c
> > [0]PETSC ERROR: --------------------------------------------
> > [0]PETSC ERROR: User provided function() line 0 in unknown directory
> > unknown file
> > [0]PETSC ERROR: Signal received!
> > [0]PETSC ERROR:  !
> >
> >
> > Suggestions and tips are very welcome!
> > /nisse
> >
> >
> > --
> > Nils-Erik Svang�rd
> > MSN: schweingaard at hotmail.com
> > Skype: schweingaard
> >
> >
>


--
Nils-Erik Svangård
MSN: schweingaard at hotmail.com
Skype: schweingaard




More information about the petsc-users mailing list