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

[David Jiawei LUO LIANG]

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:  


SNESSetFunction(snes, r, FormFunction2nd, &user);  
SNESSetJacobian(snes, J, J, FormJacobian2nd, &user);  


it also requires:  


SNESNewtonALSetFunction(snes, ComputeTangentLoad, &user);  


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:  
   (1) A predictor step, which predicts the intersection of the constraint surface and the tangent;  
   (2) A corrector step, which uses Newton-Raphson iterations to refine the intersection point between the constraint surface and the actual equilibrium path.  
   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,  



David





David Jiawei LUO LIANG



南方科技大学/学生/研究生/2024



广东省深圳市南山区学苑大道1088号




 
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20250317/17befbe1/attachment-0001.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: LeeFrame.zip
Type: application/octet-stream
Size: 20644 bytes
Desc: not available
URL: <http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20250317/17befbe1/attachment-0001.obj>