SNES Problem
Barry Smith
bsmith at mcs.anl.gov
Tue Mar 14 13:20:17 CST 2006
On Tue, 14 Mar 2006, Nils Erik Svangård wrote:
> Hi,
> Yes, there is a difference between the old 3-stage RK-solver and the
> SNESsolver using RK. In the old a timestep was calculatade for each
^^^^^^^^^^^^^^^
Are you using PETSc's explicit RK? If so that does not use SNES.
> cell based on a stability criteria and if chosen it could run in
> time-accurate mode eg. the timestep is the same for each cell, in SNES
^^^^^^^^^
Do you mean for backward Euler for the explicit RK? Or both?
> I use a timestep that is the same for all cells and
> changes only if it returns non-physical answers or the solution diverges.
I think you need to put the "proper" computation of the time-step
based on the stability criteria that you have in the "old" code.
Unfortunately this is currently hardwired in src/ts/impls/explicit/rk/rk.c
(search for "Computing next stepsize"). I think you need to use the minimum
of the timestep computed here and what is needed for the "stability criteria".
>
> Another interesting observation is that when using back euler at the
> moment when it bugs out petsc has assigned neigboring cells pressure
> values that differs about 1 ATM, which is a very non physical value,
Likely some stability problem.
Barry
> but I have no clue as to why petsc want such a value in thesolution.
>
> /nisse
>
> On 3/3/06, Barry Smith <bsmith at mcs.anl.gov> wrote:
>>
>> There must be a difference between the two RKs. The initial
>> timestep it uses, or algorithm for changing time step or slight
>> difference in algorithms?
>>
>> Barry
>>
>> On Fri, 3 Mar 2006, Nils Erik Svangård wrote:
>>
>>> A interesting thing is when I change TS_BEULER to TS_RUNGE_KUTTA the
>>> solver exits fast due to the non physical values. My original solver
>>> is a 3-step Runge-Kutta which works flawlessly.
>>>
>>> The error when running RK:
>>>
>>> ./ImplicitG3D on a linux-gnu named sethnx004.vac.com by yy26539 Fri
>>> Mar 3 10:31:27 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] TS user right-hand-side function line 0 unknownunknown
>>> [0]PETSC ERROR: [0] TSComputeRHSFunction line 282 src/ts/interface/ts.c
>>> [0]PETSC ERROR: [0] TSRkqs line 256 src/ts/impls/explicit/rk/rk.c
>>> [0]PETSC ERROR: [0] TSStep_Rk line 322 src/ts/impls/explicit/rk/rk.c
>>> [0]PETSC ERROR: [0] TSStep line 1520 src/ts/interface/ts.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: !
>>> aborting job:
>>> application called MPI_Abort(MPI_COMM_WORLD, 59) - process 0
>>>
>>> On 3/3/06, Nils Erik Svangård <nilserik at gmail.com> wrote:
>>>> Hello again,
>>>>
>>>>> Is it just not compiling because PETSC_ERR_ARG_DOMAIN is not
>>>>> defined? Or does it still not work even when you have put the number in?
>>>>> What is the error? Or did you get it all working with SNES?
>>>>>
>>>>
>>>>
>>>> It is not compiling because PETSC_ERR_ARG_DOMAIN is not defined:
>>>>
>>>> error on line 585 of SNESG3D2.F: data type is undefined for variable
>>>> PETSC_ERR_ARG_DOMAIN
>>>>
>>>> I tried putting a number in, then it compiles but complains while
>>>> linking that the Push() and Pop() functions is not defined. I made a
>>>> workaround by adding a variable in a common block and setting it
>>>> before returning from a non physical value in the FormFunction. Now I
>>>> can compile and when there is a non-physical value it returns and
>>>> resumes with a smaller timestep.
>>>> SNES still give non physical values which are hard to add checks for
>>>> in the program (e.g sqrt(RO(L)-RO(L+1))).
>>>>
>>>> Did I understand correctly that SNES does a linesearch and probably
>>>> goes for a solution which isnt the one I looking for? Does SNES use
>>>> the linesearch algorithm when I use -snes_mf (which I use)?
>>>>
>>>>> It may work automatically with TS; you would just set TS with a
>>>>> smaller timestep and call TSSolve again?
>>>>>
>>>>
>>>> OK, I'll add a check and put TSSolve in a while loop.
>>>>
>>>> My SNESSolve loop looks like this right now:
>>>> 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
>>>> bang=0
>>>> C While loop
>>>> 10 if (step .le. 100) then
>>>> C call PetscExceptionPush(87,ierr)
>>>> C call PetscExceptionPush(PETSC_ERR_ARG_DOMAIN,ierr)
>>>> call SNESSolve(snes,PETSC_NULL_OBJECT,x,ierr)
>>>> ierr=bang
>>>> C call PetscExceptionPop(87,ierr2)
>>>> IF ( ierr .EQ. 87) THEN
>>>> timestep=timestep/2
>>>> call VecCopy(xb,x,ierr)
>>>> write(6,*) "No physical solution, trying",
>>>> . " with smaller timestep: ",timestep
>>>> bang=0
>>>> ELSE
>>>> call SNESGetConvergedReason(snes,reason,ierr)
>>>> IF ( reason .LT. 0 ) THEN
>>>> timestep=timestep/2
>>>> call VecCopy(xb,x,ierr)
>>>> write(6,*)"Solution diverged, trying with",
>>>> . " smaller timestep: ",timestep
>>>> ELSE
>>>> timestep=dt
>>>> step=step+1
>>>> call VecCopy(x,xb,ierr)
>>>> write(6,*) "Solve success, going for step: ",step
>>>> write(6,*) "Reason: ",reason," ierr: ",ierr
>>>> reason=0
>>>> ierr=0
>>>> ENDIF
>>>> ENDIF
>>>> goto 10
>>>> endif
>>>> C End While loop
>>>>
>>>> And the check in the FormFunction:
>>>>
>>>> if (ET(L) .LT. 0) then
>>>> ierr = 87
>>>> bang = 87
>>>> C SETERRQ(87,"Non physical value!",ierr)
>>>> return
>>>> elseif ( RO(L) .LT. 0 ) then
>>>> ierr = 87
>>>> bang = 87
>>>> return
>>>> endif
>>>>
>>>>
>>>> /nisse
>>>>
>>>>
>>>>
>>>>> Barry
>>>>>
>>>>> I'm adding a petscerror.h for Fortran that includes these values
>>>>> Sorry I forgot we hadn't given them for Fortran.
>>>>>
>>>>> On Thu, 2 Mar 2006, Nils Erik Svangård wrote:
>>>>>
>>>>>> 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
>>>>>>
>>>>>>
>>>>>
>>>>
>>>>
>>>> --
>>>> Nils-Erik Svangård
>>>> E-Mail: nilserik at gmail.com
>>>> MSN: schweingaard at hotmail.com
>>>> Skype: schweingaard
>>>> Mobil: +46-(0)70-3612178
>>>>
>>>
>>>
>>> --
>>> Nils-Erik Svangård
>>> E-Mail: nilserik at gmail.com
>>> MSN: schweingaard at hotmail.com
>>> Skype: schweingaard
>>> Mobil: +46-(0)70-3612178
>>>
>>>
>>
>
>
> --
> Nils-Erik Svangård
> E-Mail: nilserik at gmail.com
> MSN: schweingaard at hotmail.com
> Skype: schweingaard
> Mobil: +46-(0)70-3612178
>
>
More information about the petsc-users
mailing list