[petsc-users] Tao iterations
Justin Chang
jychang48 at gmail.com
Mon Jun 8 05:16:53 CDT 2015
Hi Jason,
One more question about BLMVM... if it only uses gradient information and
does not require the definitition of a Hessian Matrix, can this method be
applied to solve problems that are nonsymmetric by nature? (e.g.,
advection-diffusion equations). If I had wanted to solve the same equations
using TRON I would need to rewrite the problem as "normal equations" (i.e.,
min 1/2*||K*u - f||^2 or min 1/2*u^T*K^T*K*u - u^T*K^T*f) so that I get a
symmetric Hessian (but this method results in extremely high condition
numbers and renders the numerical solution less reliable and accurate).
Thanks,
Justin
On Wed, Apr 29, 2015 at 5:53 PM, Justin Chang <jchang27 at uh.edu> wrote:
> Okay that's what I figured, thanks you very much
>
>
>
> On Wed, Apr 29, 2015 at 10:39 AM, Jason Sarich <jason.sarich at gmail.com>
> wrote:
>
>> Hi Justin,
>>
>> This expected behavior due to the accumulation of numerical round-offs.
>> If this is a problem or if you just want to confirm that this is the cause,
>> you can try configuring PETSc for quad precision
>> (--with-precision=__float128, works with GNU compilers) and the results
>> should match better.
>>
>> Jason
>>
>>
>> On Tue, Apr 28, 2015 at 10:19 PM, Justin Chang <jchang27 at uh.edu> wrote:
>>
>>> Jason (or anyone),
>>>
>>> I am noticing that the iteration numbers reported by
>>> TaoGetSolutionStatus() for blmvm differ whenever I change the number of
>>> processes. The solution seems to remain the same though. Is there a reason
>>> why this could be happening?
>>>
>>> Thanks,
>>>
>>> On Tue, Apr 21, 2015 at 10:40 AM, Jason Sarich <jason.sarich at gmail.com>
>>> wrote:
>>>
>>>> Justin,
>>>>
>>>> 1) The big difference between TRON and BLMVM is that TRON requires
>>>> hessian information, BLMVM only uses gradient information. Thus TRON will
>>>> usually converge faster, but requires more information, memory, and a KSP
>>>> solver. GPCG (gradient projected conjugate gradient) is another
>>>> gradient-only option, but usually performs worse than BLMVM.
>>>>
>>>> 2) TaoGetLinearSolveIterations() will get the total number of KSP
>>>> iterations per solve
>>>>
>>>> Jason
>>>>
>>>>
>>>> On Tue, Apr 21, 2015 at 10:33 AM, Justin Chang <jychang48 at gmail.com>
>>>> wrote:
>>>>
>>>>> Jason,
>>>>>
>>>>> Tightening the tolerances did the trick. Thanks. Though I do have a
>>>>> couple more related questions:
>>>>>
>>>>> 1) Is there a general guideline for choosing tron over blmvm or vice
>>>>> versa? Also is there another tao type that is also suitable given only
>>>>> bounded constraints?
>>>>>
>>>>> 2) Is it possible to obtain the total number of KSP and/or PG
>>>>> iterations from tron?
>>>>>
>>>>> Thanks,
>>>>> Justin
>>>>>
>>>>> On Tue, Apr 21, 2015 at 9:52 AM, Jason Sarich <jason.sarich at gmail.com>
>>>>> wrote:
>>>>>
>>>>>> Hi Justin,
>>>>>>
>>>>>> blmvm believes that it is already sufficiently close to a minimum,
>>>>>> so it doesn't do anything. You may need to tighten some of the tolerance to
>>>>>> force an iteration.
>>>>>>
>>>>>> Jason
>>>>>>
>>>>>>
>>>>>> On Tue, Apr 21, 2015 at 9:48 AM, Justin Chang <jychang48 at gmail.com>
>>>>>> wrote:
>>>>>>
>>>>>>> Time step 1:
>>>>>>>
>>>>>>> Tao Object: 1 MPI processes
>>>>>>> type: blmvm
>>>>>>> Gradient steps: 0
>>>>>>> TaoLineSearch Object: 1 MPI processes
>>>>>>> type: more-thuente
>>>>>>> Active Set subset type: subvec
>>>>>>> convergence tolerances: fatol=0.0001, frtol=0.0001
>>>>>>> convergence tolerances: gatol=0, steptol=0, gttol=0
>>>>>>> Residual in Function/Gradient:=0.0663148
>>>>>>> Objective value=-55.5945
>>>>>>> total number of iterations=35, (max: 2000)
>>>>>>> total number of function/gradient evaluations=37, (max: 4000)
>>>>>>> Solution converged: estimated |f(x)-f(X*)|/|f(X*)| <= frtol
>>>>>>>
>>>>>>> Time step 2:
>>>>>>>
>>>>>>> Tao Object: 1 MPI processes
>>>>>>> type: blmvm
>>>>>>> Gradient steps: 0
>>>>>>> TaoLineSearch Object: 1 MPI processes
>>>>>>> type: more-thuente
>>>>>>> Active Set subset type: subvec
>>>>>>> convergence tolerances: fatol=0.0001, frtol=0.0001
>>>>>>> convergence tolerances: gatol=0, steptol=0, gttol=0
>>>>>>> Residual in Function/Gradient:=0.0682307
>>>>>>> Objective value=-66.9675
>>>>>>> total number of iterations=23, (max: 2000)
>>>>>>> total number of function/gradient evaluations=25, (max: 4000)
>>>>>>> Solution converged: estimated |f(x)-f(X*)|/|f(X*)| <= frtol
>>>>>>>
>>>>>>> Time step 3:
>>>>>>>
>>>>>>> Tao Object: 1 MPI processes
>>>>>>> type: blmvm
>>>>>>> Gradient steps: 0
>>>>>>> TaoLineSearch Object: 1 MPI processes
>>>>>>> type: more-thuente
>>>>>>> Active Set subset type: subvec
>>>>>>> convergence tolerances: fatol=0.0001, frtol=0.0001
>>>>>>> convergence tolerances: gatol=0, steptol=0, gttol=0
>>>>>>> Residual in Function/Gradient:=0.0680522
>>>>>>> Objective value=-71.8211
>>>>>>> total number of iterations=19, (max: 2000)
>>>>>>> total number of function/gradient evaluations=22, (max: 4000)
>>>>>>> Solution converged: estimated |f(x)-f(X*)|/|f(X*)| <= frtol
>>>>>>>
>>>>>>> Time step 4:
>>>>>>>
>>>>>>> Tao Object: 1 MPI processes
>>>>>>> type: blmvm
>>>>>>> Gradient steps: 0
>>>>>>> TaoLineSearch Object: 1 MPI processes
>>>>>>> type: more-thuente
>>>>>>> Active Set subset type: subvec
>>>>>>> convergence tolerances: fatol=0.0001, frtol=0.0001
>>>>>>> convergence tolerances: gatol=0, steptol=0, gttol=0
>>>>>>> Residual in Function/Gradient:=0.0551556
>>>>>>> Objective value=-75.1252
>>>>>>> total number of iterations=18, (max: 2000)
>>>>>>> total number of function/gradient evaluations=20, (max: 4000)
>>>>>>> Solution converged: estimated |f(x)-f(X*)|/|f(X*)| <= frtol
>>>>>>>
>>>>>>> Time step 5:
>>>>>>>
>>>>>>> Tao Object: 1 MPI processes
>>>>>>> type: blmvm
>>>>>>> Gradient steps: 0
>>>>>>> TaoLineSearch Object: 1 MPI processes
>>>>>>> type: more-thuente
>>>>>>> Active Set subset type: subvec
>>>>>>> convergence tolerances: fatol=0.0001, frtol=0.0001
>>>>>>> convergence tolerances: gatol=0, steptol=0, gttol=0
>>>>>>> Residual in Function/Gradient:=0.0675667
>>>>>>> Objective value=-77.4414
>>>>>>> total number of iterations=6, (max: 2000)
>>>>>>> total number of function/gradient evaluations=8, (max: 4000)
>>>>>>> Solution converged: estimated |f(x)-f(X*)|/|f(X*)| <= frtol
>>>>>>>
>>>>>>> Time step 6:
>>>>>>>
>>>>>>> Tao Object: 1 MPI processes
>>>>>>> type: blmvm
>>>>>>> Gradient steps: 0
>>>>>>> TaoLineSearch Object: 1 MPI processes
>>>>>>> type: more-thuente
>>>>>>> Active Set subset type: subvec
>>>>>>> convergence tolerances: fatol=0.0001, frtol=0.0001
>>>>>>> convergence tolerances: gatol=0, steptol=0, gttol=0
>>>>>>> Residual in Function/Gradient:=0.059143
>>>>>>> Objective value=-79.5007
>>>>>>> total number of iterations=3, (max: 2000)
>>>>>>> total number of function/gradient evaluations=5, (max: 4000)
>>>>>>> Solution converged: estimated |f(x)-f(X*)|/|f(X*)| <= frtol
>>>>>>>
>>>>>>> Time step 7:
>>>>>>>
>>>>>>> Tao Object: 1 MPI processes
>>>>>>> type: blmvm
>>>>>>> Gradient steps: 0
>>>>>>> TaoLineSearch Object: 1 MPI processes
>>>>>>> type: more-thuente
>>>>>>> Active Set subset type: subvec
>>>>>>> convergence tolerances: fatol=0.0001, frtol=0.0001
>>>>>>> convergence tolerances: gatol=0, steptol=0, gttol=0
>>>>>>> Residual in Function/Gradient:=0.0433683
>>>>>>> Objective value=-81.3546
>>>>>>> total number of iterations=5, (max: 2000)
>>>>>>> total number of function/gradient evaluations=8, (max: 4000)
>>>>>>> Solution converged: estimated |f(x)-f(X*)|/|f(X*)| <= frtol
>>>>>>>
>>>>>>> Time step 8:
>>>>>>>
>>>>>>> Tao Object: 1 MPI processes
>>>>>>> type: blmvm
>>>>>>> Gradient steps: 0
>>>>>>> TaoLineSearch Object: 1 MPI processes
>>>>>>> type: more-thuente
>>>>>>> Active Set subset type: subvec
>>>>>>> convergence tolerances: fatol=0.0001, frtol=0.0001
>>>>>>> convergence tolerances: gatol=0, steptol=0, gttol=0
>>>>>>> Residual in Function/Gradient:=0.0840676
>>>>>>> Objective value=-82.9382
>>>>>>> total number of iterations=0, (max: 2000)
>>>>>>> total number of function/gradient evaluations=1, (max: 4000)
>>>>>>> Solution converged: estimated |f(x)-f(X*)|/|f(X*)| <= frtol
>>>>>>>
>>>>>>> Time step 9:
>>>>>>>
>>>>>>> Tao Object: 1 MPI processes
>>>>>>> type: blmvm
>>>>>>> Gradient steps: 0
>>>>>>> TaoLineSearch Object: 1 MPI processes
>>>>>>> type: more-thuente
>>>>>>> Active Set subset type: subvec
>>>>>>> convergence tolerances: fatol=0.0001, frtol=0.0001
>>>>>>> convergence tolerances: gatol=0, steptol=0, gttol=0
>>>>>>> Residual in Function/Gradient:=0.0840676
>>>>>>> Objective value=-82.9382
>>>>>>> total number of iterations=0, (max: 2000)
>>>>>>> total number of function/gradient evaluations=1, (max: 4000)
>>>>>>> Solution converged: estimated |f(x)-f(X*)|/|f(X*)| <= frtol
>>>>>>>
>>>>>>> Time step 10:
>>>>>>>
>>>>>>> Tao Object: 1 MPI processes
>>>>>>> type: blmvm
>>>>>>> Gradient steps: 0
>>>>>>> TaoLineSearch Object: 1 MPI processes
>>>>>>> type: more-thuente
>>>>>>> Active Set subset type: subvec
>>>>>>> convergence tolerances: fatol=0.0001, frtol=0.0001
>>>>>>> convergence tolerances: gatol=0, steptol=0, gttol=0
>>>>>>> Residual in Function/Gradient:=0.0840676
>>>>>>> Objective value=-82.9382
>>>>>>> total number of iterations=0, (max: 2000)
>>>>>>> total number of function/gradient evaluations=1, (max: 4000)
>>>>>>> Solution converged: estimated |f(x)-f(X*)|/|f(X*)| <= frtol
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> On Tue, Apr 21, 2015 at 9:28 AM, Jason Sarich <
>>>>>>> jason.sarich at gmail.com> wrote:
>>>>>>>
>>>>>>>> Hi Justin,
>>>>>>>>
>>>>>>>> what reason is blmvm giving for stopping the solve? (you can use
>>>>>>>> -tao_view or -tao_converged_reason to get this)
>>>>>>>>
>>>>>>>> Jason
>>>>>>>>
>>>>>>>> On Mon, Apr 20, 2015 at 6:32 PM, Justin Chang <jychang48 at gmail.com>
>>>>>>>> wrote:
>>>>>>>>
>>>>>>>>> Jason,
>>>>>>>>>
>>>>>>>>> I am using TaoGetSolutionStatus(tao,&its, ...) and it gives me
>>>>>>>>> exactly what I want. However, I seem to be having an issue with blmvm
>>>>>>>>>
>>>>>>>>> I wrote my own backward euler code for a transient linear
>>>>>>>>> diffusion problem with lower bounds >= 0 and upper bounds <= 1. For the
>>>>>>>>> first several time steps I am getting its > 0, and it decreases over time
>>>>>>>>> due to the nature of the discrete maximum principles. However, at some
>>>>>>>>> point my its become 0 and the solution does not "update", which seems to me
>>>>>>>>> that TaoSolve is not doing anything after that. This doesn't happen if I
>>>>>>>>> were to use tron (my KSP and PC are cg and jacobi respectively).
>>>>>>>>>
>>>>>>>>> Do you know why this behavior may occur?
>>>>>>>>>
>>>>>>>>> Thanks,
>>>>>>>>>
>>>>>>>>> On Tue, Apr 14, 2015 at 9:35 AM, Jason Sarich <
>>>>>>>>> jason.sarich at gmail.com> wrote:
>>>>>>>>>
>>>>>>>>>> Hi Justin,
>>>>>>>>>>
>>>>>>>>>> I have pushed these changes to the "next" branch, your code
>>>>>>>>>> snippet should work fine there.
>>>>>>>>>>
>>>>>>>>>> Note that there is also available (since version 3.5.0) the
>>>>>>>>>> routine TaoGetSolutionStatus(tao,&its,NULL,NULL,NULL,NULL,NULL) which will
>>>>>>>>>> provide the
>>>>>>>>>> same information
>>>>>>>>>>
>>>>>>>>>> Jason
>>>>>>>>>>
>>>>>>>>>> On Fri, Apr 10, 2015 at 6:28 PM, Justin Chang <
>>>>>>>>>> jychang48 at gmail.com> wrote:
>>>>>>>>>>
>>>>>>>>>>> Whatever is convenient and/or follow the "PETSc" standards.
>>>>>>>>>>> Something similar to SNESGetIterationNumber() or KSPGetIterationNumber()
>>>>>>>>>>> would be nice. Ideally I want my code to look like this:
>>>>>>>>>>>
>>>>>>>>>>> ierr = TaoGetIterationNumber(tao,&its);CHKERRQ(ierr);
>>>>>>>>>>> ierr = PetscPrintf(PETSC_COMM_WORLD, "Number of Tao iterations
>>>>>>>>>>> = %D\n", its);
>>>>>>>>>>>
>>>>>>>>>>> Thanks :)
>>>>>>>>>>>
>>>>>>>>>>> On Fri, Apr 10, 2015 at 5:53 PM, Jason Sarich <
>>>>>>>>>>> jason.sarich at gmail.com> wrote:
>>>>>>>>>>>
>>>>>>>>>>>> Hi Justin, I'll get this in. I assume that displaying the
>>>>>>>>>>>> number of iterations with tao_converged_reason is what you are asking for
>>>>>>>>>>>> in particular? Or did you have something else in mind?
>>>>>>>>>>>>
>>>>>>>>>>>> Jason
>>>>>>>>>>>> On Apr 10, 2015 16:42, "Smith, Barry F." <bsmith at mcs.anl.gov>
>>>>>>>>>>>> wrote:
>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> Justin,
>>>>>>>>>>>>>
>>>>>>>>>>>>> Sorry TAO simply doesn't even collect this information
>>>>>>>>>>>>> currently. But yes we should definitely make it available!
>>>>>>>>>>>>>
>>>>>>>>>>>>> Jason,
>>>>>>>>>>>>>
>>>>>>>>>>>>> Could you please add this; almost all the TaoSolve_xxx()
>>>>>>>>>>>>> have a local variable iter; change that to tao->niter (I'm guess this is
>>>>>>>>>>>>> suppose to capture this information) and add a TaoGetIterationNumber() and
>>>>>>>>>>>>> the uses can access this. Also modify at the end of TaoSolve()
>>>>>>>>>>>>> -tao_converged_reason to also print the iteration count. At the same time
>>>>>>>>>>>>> since you add this you can add a tao->totalits which would accumulate all
>>>>>>>>>>>>> iterations over all the solves for that Tao object and the routine
>>>>>>>>>>>>> TaoGetTotalIterations() to access this. Note that TaoSolve() would
>>>>>>>>>>>>> initialize tao->niter = 0 at the top.
>>>>>>>>>>>>>
>>>>>>>>>>>>> Thanks
>>>>>>>>>>>>>
>>>>>>>>>>>>> Barry
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> > On Apr 10, 2015, at 4:16 PM, Justin Chang <jchang27 at uh.edu>
>>>>>>>>>>>>> wrote:
>>>>>>>>>>>>> >
>>>>>>>>>>>>> > Hi all,
>>>>>>>>>>>>> >
>>>>>>>>>>>>> > Is there a way to generically obtain the number of Tao
>>>>>>>>>>>>> iterations? I am looking through the -help options for Tao and I don't see
>>>>>>>>>>>>> any metric where you can output this quantity in the manner that you could
>>>>>>>>>>>>> for SNES or KSP solves. I am currently using blmvm and tron, and the only
>>>>>>>>>>>>> way I can see getting this metric is by outputting -tao_view and/or
>>>>>>>>>>>>> -tao_monitor and manually finding this number. I find this cumbersome
>>>>>>>>>>>>> especially for transient problems where I would like to simply have this
>>>>>>>>>>>>> number printed for each step instead of ending up with unnecessary info.
>>>>>>>>>>>>> >
>>>>>>>>>>>>> > Thanks,
>>>>>>>>>>>>> >
>>>>>>>>>>>>> >
>>>>>>>>>>>>> > --
>>>>>>>>>>>>> > Justin Chang
>>>>>>>>>>>>> > PhD Candidate, Civil Engineering - Computational Sciences
>>>>>>>>>>>>> > University of Houston, Department of Civil and Environmental
>>>>>>>>>>>>> Engineering
>>>>>>>>>>>>> > Houston, TX 77004
>>>>>>>>>>>>> > (512) 963-3262
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>
>>>>>>>>
>>>>>>>
>>>>>>
>>>>>
>>>>
>>>
>>>
>>> --
>>> Justin Chang
>>> PhD Candidate, Civil Engineering - Computational Sciences
>>> University of Houston, Department of Civil and Environmental Engineering
>>> Houston, TX 77004
>>> (512) 963-3262
>>>
>>
>>
>
>
> --
> Justin Chang
> PhD Candidate, Civil Engineering - Computational Sciences
> University of Houston, Department of Civil and Environmental Engineering
> Houston, TX 77004
> (512) 963-3262
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20150608/452df7ef/attachment-0001.html>
More information about the petsc-users
mailing list