forum.alglib.net http://forum.alglib.net/ 

achieving ~e10 KKT feasibility for wellscaled problem http://forum.alglib.net/viewtopic.php?f=2&t=4606 
Page 1 of 1 
Author:  8699oriel [ Mon May 27, 2024 2:45 pm ] 
Post subject:  achieving ~e10 KKT feasibility for wellscaled problem 
Hi, Apologies for the lengthy context before asking my question. I have a simple/small subproblem that I will solve thousands of times within a framework to solve a larger problem; the subproblem maximizes the distance (for this purpose, Euclidean only) between two convex neighborhoods. For simplicity's sake, let them be two ellipsoids with explicitly fairly different radii with no additional linear constraints imposed. I have multiple means of solving this problem to global optimality including, now (thank you!!!!), alglib's options (filter SQP, AUL, etc) that I just implemented/tested; I also have tested SNOPT via dll, SNOPT, CONOPT, MINOS and COUENNE via GAMS, and a geometric parameter search method I wrote. I say 'globally optimal' on purpose because I have a geometric partitioning scheme that mimics what a branch/cut implementation would do here when used in conjunction with a local optimizer. I provide more context in this recent unanswered post: https://math.stackexchange.com/questions/4917464/analyticalsolutiontoatoysizedbutinterestingoptimizationproblem. Note that I've correctly modified the object function and gradient to solve as a minimization within ALGLIB and match the expected solution. I solve the problem in what I call unitcanonical space (the ellipsoids are transformed into axisaligned unit spheres centered at the Cartesian origin) because this provides numerical stability over many orders of magnitude difference in radii sizes (so my ALGLIB scaling is simply a unitvector). The KKT conditions are quite 'simple' in that I have only 6 variables (the xyzposition for each neighborhood) and 2 quadratic constraints, so I have 6 equations to evaluate: df/dx lambda1*dg1/dx lambda2*dg2/dx. A further simplification is that since the dg1 vector is nonzero for only ellipsoid 1 and the dg2 vector is nonzero for only ellipsoid 2, I end up with 3 equations for dflambda1*dg1 and 3 for dflambda2*dg2 (because the other lambda is irrelevant in each case). Here are two ALGLIB questions that are probably simple to you but not to me: 1. I have the analytical gradient (it's linear) and the analytical Hessian (it's a constant matrix) for f and g but there seem to be no algorithms making use of the exact Hessian  e.g. filter SQP algorithms like SNOPT talk about an approximate Hessian and a modified Lagrangian to deal with large and/or highly nonlinear problems; I see nothing in minnlcsetalgosqp or the other algorithms where I can use that. Is there a means of using this within ALGLIB's optimization schemes? Are there algorithms for simple problems like mine that utilize the analytical Hessian? 2. ALGLIB and SNOPT's optimal (termination status 1) objective function value matches mine but the KKT conditions that I calculate are infeasible at the 1e3 to 1e5 level whereas in my (slow) geometric algorithm, I can sometimes get to 1e11 to 1e13 even for highly ellipsoidal cases. Even when I hack the ALGLIB sqp code to force things like egap and edual to be < 1e14 instead of relative, the solution still is KKt feasible only at ~1e4. Where can I extract the options for KKT multipliers that ALGLIB computes so I can see what ALGLIB says df/dxlambda1*dg1/dxlambda2*dg2/dx are? All I can find are the box constraint and linear constraint KKT multipliers which are not part of my current problem? What can I do in ALGLIB to find a solution feasible to the KKT conditions at the ~e10 level? Thank you, William 
Author:  8699oriel [ Mon May 27, 2024 10:49 pm ] 
Post subject:  Re: achieving ~e10 KKT feasibility for wellscaled problem 
If it helps, I'm willing to post the solution found by alglib fsqp, the KKT conditions I derived, the multipliers and evaluations produced from the KKT point as well as my own solution. I'm convinced there is a better way to solve this very small problem and I'd like to get the KKT feasibility down to ~e12 in general. 
Author:  Sergey.Bochkanov [ Tue May 28, 2024 7:07 pm ] 
Post subject:  Re: achieving ~e10 KKT feasibility for wellscaled problem 
Hi! Would you like to run ALGLIB with logging enabled, i.e. call alglib::trace_file("sqp","trace.log") prior to running the optimizer? I would like to take a closer look at the optimization process. The fact that the objective value shows good agreement, but KKT error is not so good, suggests that there is some kind of degeneracy here. Sergey. 
Author:  8699oriel [ Wed May 29, 2024 11:09 am ]  
Post subject:  Re: achieving ~e10 KKT feasibility for wellscaled problem  
Thank you for reading my post and or the suggestion. I've attached the trace file as well as a more detailed one. FWIW, to show how similar the solution values are, here is the alglib solution. Please excuse the many decimal places. I generally ignore Visual Studio C++ numerical places beyond 1e13 but I value extreme precision in my problem. Alglib filter sqp: UnitCanonical Ellipsoid 1: (0.329176838420615, 0.94401746027854, 0.0217633576606146) UnitCanonical Ellipsoid 2: (0.930895021061298, 0.329232641233283, 0.158241359036905) In general space, this corresponds to General Ellipsoid 1: (10.5967582185939, 206.347913283321, 6.77869831250489) General Ellipsoid 2: (18.8893523034911, 49.9697297257856, 570.805604711329) The Euclidean norm is: 620.237467405261 In code, of the 6 KKT condition equations, the worst feasibility is 9.7502106655156240e06 and in my Excel "implementation" of my geometric convergence algorithm, the KKT feasibility evaluates as 0.0001294260728173 or roughly 1e4. In contrast, the point I find in my Excel "algorithm" differs from the ALGLIB solution by at most e9 in unitcanonical space and e8 in general space! My solution: UnitCanonical Ellipsoid 1: (0.32917684136241, 0.944017459191415, 0.0217633603207913) UnitCanonical Ellipsoid 2: (0.930895019622935, 0.329232643072075, 0.158241363672689) In general space, this corresponds to General Ellipsoid 1: (10.5967582026675, 206.347913213863, 6.77869828010731) General Ellipsoid 2: (18.8893522318179, 49.9697297427802, 570.805604707352) The Euclidean norm is: 620.237467405261 In code, of the 6 KKT condition equations, the worst feasibility is 9.0949470177292824e12 and in my Excel "implementation" of my geometric convergence algorithm, the KKT feasibility evaluates as 0.00000000003547029 or roughly 3.5e11. Does this give you any hints? I realize it may seem like polishing an already very nice apple, but I'd love to squeeze out the extra accuracy from alglib if at all possible...

Author:  8699oriel [ Wed May 29, 2024 11:13 am ] 
Post subject:  Re: achieving ~e10 KKT feasibility for wellscaled problem 
I forgot to give you an idea of the original problem. Ellipsoid 1: radii: (16.765, 55.8273, 1.3), origin: (20.0, 155.9, 19.99), rotations: (1.4, 236.83383, 12.999) degrees Ellipsoid 2: radii: (8.7, 10.1, 15.1617), origin: (20.09, 45.6, 562.94), rotations: (1.23456, 98.764523, 186.34) degrees So a 'large' difference in radii which I like because it shows robustness of the approach. More spherical ellipsoids solve to e13 KKT feasibility by any method I mentioned above. 
Author:  8699oriel [ Sat Jun 01, 2024 12:05 pm ] 
Post subject:  Re: achieving ~e10 KKT feasibility for wellscaled problem 
I just tested in GAMS (CONOPT and SNOPT) a "postprocessing" NLP that makes the KKT multiplier calculations actual variables and constraints and minimizes the infeasibility as the objective function while requiring my previous objective function value be maintained...A poor person's dual formulation in a way. I have to check for divide by zero bu It seems to only take a couple internal iterations (I guess that's the QP solver) to achieve what I want. However, it seems I should be able to get this from a primal solution. Recognizing ALGLIB advertises the "no tuning needed" for its excellent implementations, any thoughts on internal settings that might help here? Or parts of the code I could look into different from what I reference above? Thank you!! 
Author:  Sergey.Bochkanov [ Sun Jun 02, 2024 12:08 pm ] 
Post subject:  Re: achieving ~e10 KKT feasibility for wellscaled problem 
Hi! Actually, I thought for several days about using arbitrary precision to further polish solution obtained in a double precision. ALGLIB has legacy version (2.x branch) supporting MPFR, probably if you use run a 128bit Newton iteration straight from the solution returned by ALGLIB it will converge in several iterations even without merit functions, line search and other kinds of safeguards that make optimizer so difficult to implement. That legacy version does not support advanced optimization algos, but you can use its linear algebra to implement Newton method yourself. Your approach is also a good option, though. 
Author:  8699oriel [ Wed Jun 05, 2024 11:43 am ] 
Post subject:  Re: achieving ~e10 KKT feasibility for wellscaled problem 
Thank you! Yours is also an excellent idea. I will try to implement that if it doesn't take me too long to get it running. In the guts of the current alglib release, when I comment out the termination conditions that specify if the change in objective function is too small and if the trust region is smaller than the sqrt of machine epsilon (optimization.cpp ~line 63417 and line 63438), the solution will be perturbed until the trust region is smaller than the actual epsilon I input (rather than the sqrt of epsilon). Of course, since it's not optimizing anything really at that point, it's "random"....but perhaps you could add a "solution polishing" step under these conditions that used extra numerical precision, per your idea, or explicitly minimizes the infeasibility per an extended formulation. I'll let you know how your idea goes once I get to it  but may be a few weeks. Thank you for thinking about it! William 
Author:  8699oriel [ Sat Jun 08, 2024 9:38 pm ] 
Post subject:  Re: achieving ~e10 KKT feasibility for wellscaled problem 
Just an update  I think your NewtonRaphson suggestion is superior. Although my idea is certainly robust, it's far more complicated whereas this is, so far, 12 iterations of NewtonRaphson to get to 1e13 roughly...and that's simply prototyping it in Excel. I was unable to get the ALGLIB 2.6.0 mpr fully imported (I create dlls out of all 3rd party software I use and export methods I want but I can't figure out some of the unresolved externals I'm getting)...but given the success using regular 64bit Excel accuracy, I expect to be just fine without it. Hopefully, if I ever have to use it, I'll be able to just pull in some of the specific code that had no issues linking. So for now, I think this should be good enough. Now to figure out how to calculate KKT multipliers for additional linear constraints. Many thanks!! 
Page 1 of 1  All times are UTC 
Powered by phpBB © 2000, 2002, 2005, 2007 phpBB Group http://www.phpbb.com/ 