SNES Problem
Nils Erik Svangård
nilserik at gmail.com
Tue Mar 14 23:40:43 CST 2006
>
> Just get a copy of rk.c and link it in with your code; you don't have
> to muck with the actual file that goes into the PETSc library.
>
Ok.
> > In the back-euler implementation, I think I just have to add one line
> > to get the stability checking of the timestep working.
>
> If it is one line for backward Euler, why isn't one line for RK?
> You can simply call the code in fortran that computes the time-step,
> can't you then use MPI_Allreduce to get the min over all processors.
It is one line to add in the SNES implementation because there I set
timestep for each iteration. The RK is used via TS which manage the
timestep, and I am not sure how to manage the timestep myself in when
using TS,
>
> Here is the issue: the ONLY reason to use an implicit method is to
> increase the timestep which you cannot do with the explicit code due
> to stability constraints. But if you use the same time-step with the
> implicit code as with the explicit code you gain nothing over using
> the explicit code (in fact you do worse because the implicit code
> is more expensive per timestep).
So using the old timestep control function might make the code to work
but probably will reach the solution much slower. So if this work, I
probably have to rewrite the stability check to fit the implicit
timestepping?
/nisse
>
>
> Barry
>
> >
> > /nisse
> >> Barry
> >>
> >>
> >>
> >> On Tue, 14 Mar 2006, Nils Erik Svangård wrote:
> >>
> >>>>> 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.
> >>>>
> >>>
> >>> Sorry, I meant TS.
> >>>
> >>>>> 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?
> >>>
> >>> Sorry again, using TS I use the same size of the timestep on every
> >>> cell (both when trying use back-euler and RK).
> >>> In my old 3-stage RK I can use both variable timestep and
> >>> time-accurate timestep.
> >>> In my seperate back-euler implementation using SNES I use the same
> >>> timestep on every cell.
> >>>
> >>>>
> >>>>> 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".
> >>>>
> >>> Hmm.. Ok so I could use my back-euler SNES implementation with my old
> >>> timestep calculation function (using time-accurate timestep)? And it
> >>> might work better.
> >>>
> >>>
> >>> Thanks
> >>> /nisse
> >>>
> >>>>>
> >>>>> 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
> >>>>>
> >>>>>
> >>>>
> >>>
> >>>
> >>> --
> >>> 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