[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Re: SNES Problem

To: petsc-users@mcs.anl.gov
Subject: Re: SNES Problem
From: Barry Smith <bsmith@mcs.anl.gov>
Date: Mon, 27 Feb 2006 09:07:18 -0600 (CST)
In-reply-to: <82f936c20602270607t775d7e03o6a37afa4c0b3d2a8@mail.gmail.com>
References: <82f936c20602270607t775d7e03o6a37afa4c0b3d2a8@mail.gmail.com>
Reply-to: petsc-users@mcs.anl.gov
Sender: owner-petsc-users@mcs.anl.gov

  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@hotmail.com
Skype: schweingaard

References:
- SNES Problem
  - From: "Nils Erik Svangård" <nilserik@gmail.com>

Prev by Date: SNES Problem
Next by Date: Re: VecGetValues?
Previous by thread: SNES Problem
Next by thread: attaching pathdb
Index(es):
- Date
- Thread