Re: [petsc-dev] Inquiry Regarding Arc-Length Method Implementation in SNESNEWTONAL Solver

[Jed Brown]

Hi, David. I'm copying Zach, who developed this implementation and will be able 
to more directly answer your questions.

"David Jiawei LUO LIANG"        <12431...@mail.sustech.edu.cn> writes:

> Dear PETSc Developers,
>
>
> I am currently using the SNESNEWTONAL nonlinear solver and its arc-length 
> method module. I noticed from the tutorial that, besides the standard 
> residual and Jacobian definitions: &nbsp;
>
>
> SNESSetFunction(snes, r, FormFunction2nd, &amp;user); &nbsp;
> SNESSetJacobian(snes, J, J, FormJacobian2nd, &amp;user); &nbsp;
>
>
> it also requires: &nbsp;
>
>
> SNESNewtonALSetFunction(snes, ComputeTangentLoad, &amp;user); &nbsp;
>
>
> which is used to define the arc-length constraint function.
>
>
> I have several questions regarding the arc-length implementation, and I would 
> greatly appreciate your clarification:
>
>
> 1. According to the tutorial, the arc-length constraint function 
> `ComputeTangentLoad` is used to compute ( Q = dF^{ext}/d\lambda ), the 
> derivative of the external load with respect to lambda. In the case of a 
> linear load, i.e., ( F^{ext} = \lambda F^{ext} ), ( Q = F^{ext} ). In this 
> scenario, is it correct to simply return ( F^{ext} ) inside 
> `ComputeTangentLoad`?
>
>
> 2. From my experiments, I observed that the `ComputeTangentLoad` function is 
> called once before assembling the residual at every iteration. Interestingly, 
> it seems that during the first iteration, ( Q ) contributes to the residual 
> calculation, but from the second iteration onward, it no longer does. Could 
> you help me understand why this happens?
>
>
> 3. I have attached my own buckling example of a Lee Frame structure where I 
> implemented the arc-length method as the nonlinear solver. In my 
> implementation, I divided the arc-length method into two steps: &nbsp;
> &nbsp; &nbsp;(1) A predictor step, which predicts the intersection of the 
> constraint surface and the tangent; &nbsp;
> &nbsp; &nbsp;(2) A corrector step, which uses Newton-Raphson iterations to 
> refine the intersection point between the constraint surface and the actual 
> equilibrium path. &nbsp;
> &nbsp; &nbsp;This is the standard arc-length algorithm principle. I would 
> like to know if `SNESNEWTONAL` follows a similar predictor-corrector 
> procedure internally, or if I need to implement an outer load step loop and 
> predictor step myself?
>
>
> 4. I remember that the tutorial describes solving two linear systems during 
> each corrector step to obtain ( \Delta x^F ) and ( \Delta x^Q ) separately. 
> Where in the SNESNEWTONAL workflow are these two solves performed? Do they 
> happen right after assembling the residual and Jacobian?
>
>
> 5. Additionally, is the Jacobian matrix used in SNESNEWTONAL the original 
> system Jacobian, or is it augmented to include the constraint surface 
> function? Similarly, is the input vector expanded to ( n+1 ) dimensions to 
> account for ( \lambda ), or does it remain the original size ( n ) as in the 
> standard Newton-Raphson implementation?
>
>
> 6. In general, I would like to understand more about the internal iteration 
> details of `SNESSetType(snes, SNESNEWTONAL);` compared to `SNESSetType(snes, 
> SNESNEWTONLS);` so that I can more accurately replicate and improve upon the 
> arc-length algorithm I have implemented in the attached code using SNES.
>
>
> Thank you very much for your time and assistance!
>
>
> Best regards, &nbsp;
>
>
>
> David
>
>
>
>
>
> David&nbsp;Jiawei&nbsp;LUO&nbsp;LIANG
>
>
>
> 南方科技大学/学生/研究生/2024
>
>
>
> 广东省深圳市南山区学苑大道1088号
>
>
>
>
> &nbsp;