Reimplementing constraint solver

[DeepSOIC]

What i guess i am asking is, would it be possible to demonstrate in code, one assembly relation (pre-alpha version).

The funniest thing about it, the demonstration itself should happen in py console... because where else?!

[triplus]

The funniest thing about it, the demonstration itself should happen in py console... because where else?!

For the very first (pre-alpha) version, i could easily settle for a button, lets say Assembly workbench.

[DeepSOIC]

I have just opened gitter chat: https://gitter.im/FreeCAD-ellipse/community

[DeepSOIC]

A small report on progress.
Base classes for geometry and constraints are up, with just one concrete geometric thing and one constraint implemented (point and distance). This should be enough for initial experiments on solver.
Now I'm busy with the solver itself. So far, I've only reimplemented jacobi matrix and gradient calculations, and lineSearch (seems to be the first part of solve). Untested. There are a lot of things left, such as diagnose (which is a big monster do-just-about-everything procedure), and 4 solve algorithms.
Once these are done (and I for sure will find something else that absolutely needs to be implemented to test), a long-term expansion work is to follow, adding all geometries and constraints.

[ickby]

nice progress! I really looking forward for this and hope you succeed in the improvements!

[triplus]

A small report on progress.
Base classes for geometry and constraints are up, with just one concrete geometric thing and one constraint implemented (point and distance). This should be enough for initial experiments on solver.

OK, taking a quick look at the code, i guess already with Python API included. Therefore i guess all, @kbwbe, @realthunder or @Zolko could use it, as a back-end, and add a front-end, a button, in their assembly module, or lets say to Assembly workbench, just thinking out loud. To determine things like tree structure, datum features as an assembly interfaces ... After i guess @Zolko could start asking for another experimental relation (to superimpose placement in between LCS). Too bad, the solver doesn't know anything about placements.

P.S. @ickby likely has some notes left, from times when he made that video. There had to be some math involved, and things like that. Third relation? No, that would be too much, for an experimental solver. And that would i guess already start the FreeCAD 1.0 development cycle, we are still months away from releasing FreeCAD 0.19.

[DeepSOIC]

i guess already with Python API included. Therefore i guess all, @kbwbe, @realthunder or @Zolko could use it,

yes, Py API is being made almost in sync with c++, but it's nowhere near to being useful for anything, even planning, as I might want to do big changes as I dive more and more into the sketcher solver.

(to superimpose placement in between LCS).

not sure what you're talking about... but, with some pointers from wmayer, I've recently implemented placement interpolation (method sclerp of App.Placement).


One more report of progress, I've implemented dogleg solver, and now I'm busy with helper stuff to give it a try from py console.

[DeepSOIC]

Yay! DogLeg solver works

Code: Select all

import ConstraintSolver as CS
ps = CS.ParameterStore()
p1 = CS.G2D.ParaPoint(ps)
p2 = CS.G2D.ParaPoint(ps)
p2.x.Value = 3
p2.y.Value = 4

c = CS.G2D.ConstraintDistance(
    p1 = p1,
    p2 = p2,
    store = ps,
    Label = "Constraint1"
)
c.dist.Value = 3

c.update()
c.NetError

sys = CS.SubSystem()
sys.addUnknown(p1.Parameters)
sys.addUnknown(p2.Parameters)
sys.addConstraint(c)

vs = CS.ValueSet(CS.ParameterSubset(ps.allFree()))
slv = CS.SketchSolver()
slv.solveDogLeg(sys, vs)

this is the iteration log:

Code: Select all

Begin Dogleg solving
  new iteration
    err = 2.000000
    gauss-newton step and gradient step exceed trust region. 
    -> using gradient step constrained into trust region.
    trust region = 0.100000
    GN step norm = 2.500000
    SD step norm = 1.414214
    dogleg step norm = 0.100000
    err_new = 1.727157
    linearity factor = 1.000000
    expanding trust region
  new iteration
    err = 1.727157
    gauss-newton step and gradient step exceed trust region. 
    -> using gradient step constrained into trust region.
    trust region = 0.300000
    GN step norm = 2.323223
    SD step norm = 1.314214
    dogleg step norm = 0.300000
    err_new = 1.028629
    linearity factor = 1.000000
    expanding trust region
  new iteration
    err = 1.028629
    gauss-newton step and gradient step exceed trust region. 
    -> using gradient step constrained into trust region.
    trust region = 0.900000
    GN step norm = 1.792893
    SD step norm = 1.014214
    dogleg step norm = 0.900000
    err_new = 0.013045
    linearity factor = 1.000000
    expanding trust region
  new iteration
    err = 0.013045
    gauss-newton step is within trust region -> use it as the step
    trust region = 2.700000
    GN step norm = 0.201903
    SD step norm = 0.114214
    dogleg step norm = 0.201903
    err_new = 0.000003
    linearity factor = 1.000229
    (trust region unchanged)
  new iteration
    err = 0.000003
    gauss-newton step is within trust region -> use it as the step
    trust region = 2.700000
    GN step norm = 0.003154
    SD step norm = 0.001729
    dogleg step norm = 0.003154
    err_new = 0.000000
    linearity factor = 1.000000
    (trust region unchanged)
  new iteration
    err = 0.000000
    gauss-newton step is within trust region -> use it as the step
    trust region = 2.700000
    GN step norm = 0.000001
    SD step norm = 0.000000
    dogleg step norm = 0.000001
    err_new = 0.000000
    linearity factor = 1.000000
    (trust region unchanged)
  new iteration
    err = 0.000000
solved

[DeepSOIC]

Analyzing the iteration log, there are some things I don't quite understand.

It seems there is an error in the implementation. On https://optimization.mccormick.northwes ... on_methods , in the formula for rho, the actual error change is in nominator. However, in the implementation, it is upside-down:

Code: Select all

        double dL = err - 0.5*(fx + Jx*h_dl).squaredNorm();
        double dF = err - err_new;
        double rho = dL/dF;

The rest of the math seems similar. So I think I will turn it upside down in my implementation and see how it goes.. The page I linked, unfortunately, is quite poor, as it doesn't explain most of the symbols. It is however remarkably similar in symbols used to the dogleg implementation we have in sketcher, so it was useful for me to understand it.

[DeepSOIC]

Abdullah, the master of sketcher. You are the most knowledgeable person about the solver. Any input from you is highly welcome!