[petsc-users] TSSetIJacobian is never called
Matthew Knepley
knepley at gmail.com
Tue Oct 28 09:25:18 CDT 2014
On Tue, Oct 28, 2014 at 9:22 AM, Sharp Stone <thronesf at gmail.com> wrote:
> Hi Matt and Jed,
>
> Thank you for your replies. Now I see Petsc will calculate functions'
> residue in SNES before form the Jacobian. In order to get a initial
> non-zero residues of the functions, we cannot initialize some of the
> discretized spatial nodes of the functions to zero, am I right? In my
> problem, initially only a small region will be initialized with some
> values, while other spatial nodes remain zero. I don't know how the SNES
> calculate the initial residue? If it is a sum residue on all discretized
> nodes, then I should not get a zero residue, or at least, I should get a
> small residue which is not zero. But I still get the exact zero residue. So
> I wonder how the SNES calculate the residue?
>
You are calculating the residual in your callback function. Print out the
vector F which you pass back.
Matt
> Thanks in advance!
>
> TS Object: 1 MPI processes
>>>>> type: theta
>>>>> maximum steps=10
>>>>> maximum time=2e-11
>>>>> total number of nonlinear solver iterations=0
>>>>> total number of nonlinear solve failures=0
>>>>> total number of linear solver iterations=0
>>>>> total number of rejected steps=5
>>>>> Theta=1
>>>>> Extrapolation=no
>>>>> SNES Object: 1 MPI processes
>>>>> type: newtonls
>>>>> maximum iterations=50, maximum function evaluations=10000
>>>>> tolerances: relative=1e-08, absolute=1e-50, solution=1e-08
>>>>> total number of linear solver iterations=0
>>>>> total number of function evaluations=1
>>>>> SNESLineSearch Object: 1 MPI processes
>>>>> type: bt
>>>>> interpolation: cubic
>>>>> alpha=1.000000e-04
>>>>> maxstep=1.000000e+08, minlambda=1.000000e-12
>>>>> tolerances: relative=1.000000e-08, absolute=1.000000e-15,
>>>>> lambda=1.000000e-08
>>>>> maximum iterations=40
>>>>> KSP Object: 1 MPI processes
>>>>> type: gmres
>>>>> GMRES: restart=30, using Classical (unmodified) Gram-Schmidt
>>>>> Orthogonalization with no iterative refinement
>>>>> GMRES: happy breakdown tolerance 1e-30
>>>>> maximum iterations=10000, initial guess is zero
>>>>> tolerances: relative=1e-05, absolute=1e-50, divergence=10000
>>>>> left preconditioning
>>>>> using DEFAULT norm type for convergence test
>>>>> PC Object: 1 MPI processes
>>>>> type: ilu
>>>>> PC has not been set up so information may be incomplete
>>>>> ILU: out-of-place factorization
>>>>> 0 levels of fill
>>>>> tolerance for zero pivot 2.22045e-14
>>>>> using diagonal shift on blocks to prevent zero pivot [INBLOCKS]
>>>>> matrix ordering: natural
>>>>> linear system matrix = precond matrix:
>>>>> Mat Object: 1 MPI processes
>>>>> type: seqaij
>>>>> rows=20904, cols=20904, bs=4
>>>>> total: nonzeros=410816, allocated nonzeros=410816
>>>>> total number of mallocs used during MatSetValues calls =0
>>>>> using I-node routines: found 5226 nodes, limit used is 5
>>>>>
>>>>> On Mon, Oct 27, 2014 at 3:09 PM, Matthew Knepley <knepley at gmail.com>
>>>>> wrote:
>>>>>
>>>>>> On Mon, Oct 27, 2014 at 2:05 PM, Sharp Stone <thronesf at gmail.com>
>>>>>> wrote:
>>>>>>
>>>>>>> Hi Jed and Matt,
>>>>>>>
>>>>>>> Thank you very much for your replies.
>>>>>>>
>>>>>>> This time I copied the whole output from my screen, from compile to
>>>>>>> execution. This time I put the TSView() routine just after the TSSolve()
>>>>>>> routine to make the tsview output. And I intentionally to do only 5 time
>>>>>>> steps to make the test. Thank you in advance!
>>>>>>>
>>>>>>
>>>>>> You have to make the command line arguments work. Quit calling
>>>>>> TSView() in the code. We need many more options
>>>>>> for debugging. First track down why you are not processing options.
>>>>>> Do the examples take options for you?
>>>>>>
>>>>>> Once that works, run with -ts_view -ts_monitor -snes_monitor
>>>>>> -ksp_monitor_true_residual
>>>>>>
>>>>>> Matt
>>>>>>
>>>>>>
>>>>>>> host206-47:streamer_Implicit_SG GPL$ make pStreamer
>>>>>>> /Users/GPL/local/PETSc/petsc-3.5.2/arch-darwin-cxx-debug/bin/mpicxx
>>>>>>> -o pStreamer.o -c -Wall -Wwrite-strings -Wno-strict-aliasing
>>>>>>> -Wno-unknown-pragmas -g -O0 -fPIC
>>>>>>> -I/Users/GPL/local/PETSc/petsc-3.5.2/include
>>>>>>> -I/Users/GPL/local/PETSc/petsc-3.5.2/arch-darwin-cxx-debug/include
>>>>>>> -I/opt/X11/include -Wno-unused `pwd`/pStreamer.cc
>>>>>>> /Users/GPL/local/PETSc/petsc-3.5.2/arch-darwin-cxx-debug/bin/mpicc
>>>>>>> -Wl,-multiply_defined,suppress -Wl,-multiply_defined -Wl,suppress
>>>>>>> -Wl,-commons,use_dylibs -Wl,-search_paths_first
>>>>>>> -Wl,-multiply_defined,suppress -Wl,-multiply_defined -Wl,suppress
>>>>>>> -Wl,-commons,use_dylibs -Wl,-search_paths_first -fPIC -Wall
>>>>>>> -Wwrite-strings -Wno-strict-aliasing -Wno-unknown-pragmas -g3 -O0
>>>>>>> -Wno-unused -o pStreamer flux.o field.o glob_streamer.o input.o
>>>>>>> lookuptable.o nrutil.o pStreamer.o someutil.o photoionization.o heating.o
>>>>>>> -L/Users/GPL/local/PETSc/petsc-3.5.2/arch-darwin-cxx-debug/lib
>>>>>>> -L/Users/GPL/local/PETSc/petsc-3.5.2/arch-darwin-cxx-debug/lib -lpetsc
>>>>>>> -lsuperlu_4.3 -llapack -lblas -lparmetis -lmetis -L/opt/X11/lib -lX11
>>>>>>> -lpthread -lssl -lcrypto
>>>>>>> -L/Applications/Xcode.app/Contents/Developer/Toolchains/XcodeDefault.xctoolchain/usr/lib/clang/5.1/lib/darwin
>>>>>>> -lmpichf90 -lgfortran
>>>>>>> -L/usr/local/gfortran/lib/gcc/x86_64-apple-darwin13/4.8.2
>>>>>>> -L/usr/local/gfortran/lib -lgfortran -lgcc_ext.10.5 -lquadmath -lm
>>>>>>> -lclang_rt.osx -lmpichcxx -lc++
>>>>>>> -L/Applications/Xcode.app/Contents/Developer/Toolchains/XcodeDefault.xctoolchain/usr/bin/../lib/clang/5.1/lib/darwin
>>>>>>> -lclang_rt.osx
>>>>>>> -L/Users/GPL/local/PETSc/petsc-3.5.2/arch-darwin-cxx-debug/lib -ldl
>>>>>>> -lpmpich -lmpich -lopa -lmpl -lpthread -lSystem
>>>>>>> -L/Applications/Xcode.app/Contents/Developer/Toolchains/XcodeDefault.xctoolchain/usr/bin/../lib/clang/5.1/lib/darwin
>>>>>>> -L/Applications/Xcode.app/Contents/Developer/Toolchains/XcodeDefault.xctoolchain/usr/bin/../lib/clang/5.1/lib/darwin
>>>>>>> -lclang_rt.osx -ldl
>>>>>>> host206-47:streamer_Implicit_SG GPL$ mpiexec -n 1 ./pStreamer
>>>>>>>
>>>>>>> time=0.000000 step=0
>>>>>>> dt = 0.001000
>>>>>>>
>>>>>>> time=0.002000 step=1
>>>>>>> dt = 0.002000
>>>>>>>
>>>>>>> time=0.005000 step=2
>>>>>>> dt = 0.003000
>>>>>>>
>>>>>>> time=0.009000 step=3
>>>>>>> dt = 0.004000
>>>>>>>
>>>>>>> time=0.014000 step=4
>>>>>>> dt = 0.005000
>>>>>>>
>>>>>>> time=0.020000 step=5
>>>>>>> dt = 0.006000
>>>>>>>
>>>>>>> Step time elapse = 0.854144 s Total time elapse = 5.276643 s
>>>>>>> TS Object: 1 MPI processes
>>>>>>> type: theta
>>>>>>> maximum steps=10
>>>>>>> maximum time=2e-11
>>>>>>> total number of nonlinear solver iterations=0
>>>>>>> total number of nonlinear solve failures=0
>>>>>>> total number of linear solver iterations=0
>>>>>>> total number of rejected steps=5
>>>>>>> Theta=1
>>>>>>> Extrapolation=no
>>>>>>> SNES Object: 1 MPI processes
>>>>>>> type: newtonls
>>>>>>> maximum iterations=50, maximum function evaluations=10000
>>>>>>> tolerances: relative=1e-08, absolute=1e-50, solution=1e-08
>>>>>>> total number of linear solver iterations=0
>>>>>>> total number of function evaluations=1
>>>>>>> SNESLineSearch Object: 1 MPI processes
>>>>>>> type: bt
>>>>>>> interpolation: cubic
>>>>>>> alpha=1.000000e-04
>>>>>>> maxstep=1.000000e+08, minlambda=1.000000e-12
>>>>>>> tolerances: relative=1.000000e-08, absolute=1.000000e-15,
>>>>>>> lambda=1.000000e-08
>>>>>>> maximum iterations=40
>>>>>>> KSP Object: 1 MPI processes
>>>>>>> type: gmres
>>>>>>> GMRES: restart=30, using Classical (unmodified) Gram-Schmidt
>>>>>>> Orthogonalization with no iterative refinement
>>>>>>> GMRES: happy breakdown tolerance 1e-30
>>>>>>> maximum iterations=10000, initial guess is zero
>>>>>>> tolerances: relative=1e-05, absolute=1e-50, divergence=10000
>>>>>>> left preconditioning
>>>>>>> using DEFAULT norm type for convergence test
>>>>>>> PC Object: 1 MPI processes
>>>>>>> type: ilu
>>>>>>> PC has not been set up so information may be incomplete
>>>>>>> ILU: out-of-place factorization
>>>>>>> 0 levels of fill
>>>>>>> tolerance for zero pivot 2.22045e-14
>>>>>>> using diagonal shift on blocks to prevent zero pivot
>>>>>>> [INBLOCKS]
>>>>>>> matrix ordering: natural
>>>>>>> linear system matrix = precond matrix:
>>>>>>> Mat Object: 1 MPI processes
>>>>>>> type: seqaij
>>>>>>> rows=20904, cols=20904, bs=4
>>>>>>> total: nonzeros=410816, allocated nonzeros=410816
>>>>>>> total number of mallocs used during MatSetValues calls =0
>>>>>>> using I-node routines: found 5226 nodes, limit used is 5
>>>>>>>
>>>>>>>
>>>>>>> On Mon, Oct 27, 2014 at 2:46 PM, Jed Brown <jed at jedbrown.org> wrote:
>>>>>>>
>>>>>>>> Sharp Stone <thronesf at gmail.com> writes:
>>>>>>>>
>>>>>>>> > Hi Matt,
>>>>>>>> >
>>>>>>>> > Thank you very much for your reply.
>>>>>>>> >
>>>>>>>> > The ts_view output is attached below. I found my code results say
>>>>>>>> "SNES has
>>>>>>>> > not been set up", and PC "has not been set up". Does this cause
>>>>>>>> the
>>>>>>>> > problem? If so I do not see example ex17 explicitly set up the
>>>>>>>> snes object?
>>>>>>>>
>>>>>>>> Those objects are set up on the first time step (inside TSSolve).
>>>>>>>>
>>>>>>>> > PS: I don't know why the -ts_view option does not work in my
>>>>>>>> command line,
>>>>>>>>
>>>>>>>> Are you sure it doesn't output after the TSSolve? The example does.
>>>>>>>>
>>>>>>>> > so I use TSView() routine just before the TSSolve() routine.
>>>>>>>>
>>>>>>>> You should call it *after* the solve so we can see what happened
>>>>>>>> during
>>>>>>>> the solve.
>>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> --
>>>>>>> Best regards,
>>>>>>>
>>>>>>> Feng
>>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>> --
>>>>>> What most experimenters take for granted before they begin their
>>>>>> experiments is infinitely more interesting than any results to which their
>>>>>> experiments lead.
>>>>>> -- Norbert Wiener
>>>>>>
>>>>>
>>>>>
>>>>>
>>>>> --
>>>>> Best regards,
>>>>>
>>>>> Feng
>>>>>
>>>>
>>>>
>>>>
>>>> --
>>>> What most experimenters take for granted before they begin their
>>>> experiments is infinitely more interesting than any results to which their
>>>> experiments lead.
>>>> -- Norbert Wiener
>>>>
>>>
>>>
>>>
>>> --
>>> Best regards,
>>>
>>> Feng
>>>
>>
>>
>>
>> --
>> What most experimenters take for granted before they begin their
>> experiments is infinitely more interesting than any results to which their
>> experiments lead.
>> -- Norbert Wiener
>>
>
>
>
> --
> Best regards,
>
> Feng
>
--
What most experimenters take for granted before they begin their
experiments is infinitely more interesting than any results to which their
experiments lead.
-- Norbert Wiener
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20141028/363534b2/attachment.html>
More information about the petsc-users
mailing list