Page 2 of 3
Re: [Help] Goal seek macro
Posted: Sun Jan 24, 2021 8:29 am
by Chris_G
davidosterberg wrote: ↑Sat Jan 23, 2021 3:42 pm
Scipy is awesome, definitely the right way. You can even tell it to respect bounds on the parameters (i.e max and min value)
by changing method to e.g "L-BFGS-B" and setting the
bounds parameter to a list of (min, max) tuples.
Thanks.
I was interested with parameters bounds, but get lost by all available algorithms.
edwilliams16 wrote: ↑Sat Jan 23, 2021 10:49 pm
The minimization may work better using
Code: Select all
return (o2.ddGoal.Value - o2.ddResult.Value)**2
because the latter function is smooth at the minimum, and gradient descent methods will automatically take smaller steps as the minimum is approached.
Thanks, that's interesting. I'll give a try at that.
Re: [Help] Goal seek macro
Posted: Sun Jan 24, 2021 9:44 am
by kisolre
Chris_G wrote: ↑Sun Jan 24, 2021 8:29 am.
Sorry for the bump but can you suggest a way to exit the macro on failed Recompute? I tried generating an Exception but that just prints the message in report view and minimize happily continues until 100 iterations pass. Which takes time without doing any meaningful job.
This version (when the minimize starts from 1) causes the sketch in GoalSeekTest003 to fail and so does the Recompute:
Code: Select all
import FreeCAD
from scipy.optimize import minimize
sel = Gui.Selection.getSelectionEx()
if len(sel)!=2:
raise Exception("Select 2 DD objects - Parameter and Goal")
o1=sel[0].Object
o2=sel[1].Object
def compute(value):
o1.ddParameter = value[0]
print(o1.ddParameter)
if not(App.ActiveDocument.recompute()):
raise Exception("Recompuute failed, Please check model for errors (broken sketckes, toponaming, ... ")
return (o2.ddGoal.Value - o2.ddResult.Value)**2
#res = minimize(compute, o1.ddParameter, method='Nelder-Mead', options={"maxiter": 100, "disp": True})
res = minimize(compute, 1, method='Nelder-Mead', options={"maxiter": 100, "disp": True})
print("For Parameter = {}, Result = {}".format(o1.ddParameter, o2.ddResult))
Re: [Help] Goal seek macro
Posted: Sun Jan 24, 2021 10:00 am
by davidosterberg
kisolre wrote: ↑Sun Jan 24, 2021 9:44 am
Sorry for the bump but can you suggest a way to exit the macro on failed Recompute?
There is no proper way (cf
https://github.com/scipy/scipy/issues/7306).
I suggest
It should make the optimization stop, but I havent tried it.
Re: [Help] Goal seek macro
Posted: Sun Jan 24, 2021 10:37 am
by keithsloan52
Chris_G wrote: ↑Sat Jan 23, 2021 9:30 am
If you can install python.scipy, this little script seems to work :
Code: Select all
import FreeCAD
from scipy.optimize import minimize
doc1 = FreeCAD.getDocument('GoalSeekTest001')
o1 = doc1.getObject('dd')
o2 = doc1.getObject('dd001')
def compute(value):
o1.ddParameter = value[0]
doc1.recompute()
return abs(o2.ddGoal.Value - o2.ddResult.Value)
res = minimize(compute, 1, method='Nelder-Mead', options={"maxiter": 2000, "disp": True})
print("For Height = {}, Volume = {}".format(o1.ddParameter, o2.ddResult))
Minor query : Would it not br better to use o1.recompute() rather than a recompute of the whole document, which may contain a large number of objects?
Re: [Help] Goal seek macro
Posted: Sun Jan 24, 2021 11:42 am
by openBrain
keithsloan52 wrote: ↑Sun Jan 24, 2021 10:37 am
Minor query : Would it not br better to use o1.recompute() rather than a recompute of the whole document, which may contain a large number of objects?
Nope. 'o1' is just the parameter. You need to at least recompute the result in 'o2'.
Maybe doing 'o2.recompute(True)' -- 'True' is for recursive because you also need intermediate objects to be recomputed -- would be enough, but I'm not sure if recursive recomputing took into account the expression bindings.
Re: [Help] Goal seek macro
Posted: Sun Jan 24, 2021 1:00 pm
by kisolre
Code: Select all
if not(App.ActiveDocument.recompute()):
This actually does not work because Recompute returns number of objects recomputed. What I need to check is if the document is still touched after the recompute, which should be the state if recompute does not complete (fail). But I cant find how to do that.
The API docs say
Code: Select all
App::Document::isTouched()
◆ isTouched()
bool Document::isTouched ( void ) const
check if there is any touched object in this document
References App::DocumentP::objectArray.
but trying that in the console gives
Code: Select all
>>> App.ActiveDocument.isTouched()
Traceback (most recent call last):
File "<input>", line 1, in <module>
AttributeError: 'App.Document' object has no attribute 'isTouched'
>>>
I surely dont understand much Python but is this wrong?
Re: [Help] Goal seek macro
Posted: Sun Jan 24, 2021 1:19 pm
by openBrain
kisolre wrote: ↑Sun Jan 24, 2021 1:00 pm
I surely dont understand much Python but is this wrong?
I didn't find any trace of the function you pointed out.
There is 'DocumentObject::isTouched()' but as it says, it applies to objects, not directly to document.
You can use an iterator to test if there is still touched objects in the document (snippets returns 'True' if at least one object is touched, 'False' otherwise):
Code: Select all
next((True for o in App.ActiveDocument.Objects if o.State==['Touched']),False)
Another solution would be to recompute a second time the document. If output isn't 0, there is an error in the document.
HTH
Re: [Help] Goal seek macro
Posted: Sun Jan 24, 2021 2:08 pm
by kisolre
openBrain wrote: ↑Sun Jan 24, 2021 1:19 pm
I didn't find any trace of the function you pointed out.
Here is a link
https://freecad.github.io/SourceDoc/d8/ ... 6bdfca2af0
openBrain wrote: ↑Sun Jan 24, 2021 1:19 pm
Another solution would be to recompute a second time the document. If output isn't 0, there is an error in the document.
This actually works. and then the exception properly kills the minimize and the macro.
So here is the latest "complete" version:
Code: Select all
# Goal Seeek FreeCAD macro
# Author: Chris_G ,Kisolre
# License: LGPL-3.0-or-later
#
# Usage: This macro uses scipy and DynamicData WB to adjust a Parameter so that Result reaches Goal.
#
# An example is adjusting the height of water glass model so the internal volume reaches desired value (260ml for example)
#
# For entering desired values two DD objects are used.
# First one should have a property named "ddParameter" under it which should be used by expression in the model (the glass height in the above example).
# Second one should have properties named "ddResult" and "ddGoal" under it.
# ddResult should reference through expressions some value from the model dependent on ddParameter (volume of the glass (Body.Shape.Volume)).
# ddGoal is the desired value that ddResult should reach (260ml)
#
# Ctrl-select the dd objects in proper order, execute macro, wait for result. Each iteration recomputes the documend so depending on the model this can be slow.
import FreeCAD
from scipy.optimize import minimize
sel = Gui.Selection.getSelectionEx()
if len(sel)!=2:
raise Exception("Select 2 DD objects - Parameter and Goal")
o1=sel[0].Object
o2=sel[1].Object
def compute(value):
o1.ddParameter = value[0]
print(o1.ddParameter)
App.ActiveDocument.recompute()
# This will check if there are still objects to recompute - previous recompute failed.
if App.ActiveDocument.recompute()!=0:
raise Exception ("Recompuute failed, Please check model for errors (broken sketckes, toponaming, ... ")
return (o2.ddGoal.Value - o2.ddResult.Value)**2
res = minimize(compute, o1.ddParameter, method='Nelder-Mead', options={"maxiter": 100, "disp": True})
print("For Parameter = {}, Result = {}".format(o1.ddParameter, o2.ddResult))
Re: [Help] Goal seek macro
Posted: Sun Jan 24, 2021 2:20 pm
by openBrain
OK. I missed the point. This hasn't been exposed to Python API ATM, so available only in C++.
Re: [Help] Goal seek macro
Posted: Sun Jan 24, 2021 8:49 pm
by openBrain
As a hobby, below a quick and dirty seeking function using basic dichotomy
Code: Select all
SEEKMAX = 100 # max number of iterations to found a range where a solution exists
DICHOMAX = 100 # max number of iterations to found a solution inside the found range
SEEKFACTOR = 0.5 # range size around nominal value, shall be in ]0, 1[ ; 0.5 typical ; reducing may improve finding a solution on periodic functions
DICHOFACTOR = 0.5 # dichotomy factor, no valuable reason to set something else than 0.5
THRESHOLD = 1e-7 # deviation with goal at which a solution is considered found
from PySide2 import QtWidgets
import FreeCAD as App
import FreeCADGui as Gui
doc = App.ActiveDocument
ddp = doc.dd
ddg = doc.dd001
goal = ddg.ddGoal
sb = Gui.getMainWindow().findChild(QtWidgets.QProgressBar)
def recompute(val = ddp.ddParameter):
if not ddg.recompute(True):
raise Exception("Recompute error")
def goalseek(baseValue, currValue = None, itseek = 1, itdicho = 1, up = False, down = False):
global SEEKFACTOR
sb.setValue((itseek+itdicho)/(SEEKMAX+DICHOMAX)*100)
Gui.updateGui()
if itseek > SEEKMAX+1:
App.Console.PrintWarning("Maximum number of seek iterations reached without solution range\n")
return None
if itdicho > DICHOMAX+1:
App.Console.PrintWarning("Maximum number of dichotomy iterations reached without solution\n")
return None
if currValue:
ddp.ddParameter = currValue
else:
ddp.ddParameter = baseValue
recompute()
currResult = ddg.ddResult
## SEEK phase
if up == down:
upValue = baseValue + baseValue*SEEKFACTOR
ddp.ddParameter = upValue
recompute()
upResult = ddg.ddResult
downValue = baseValue - baseValue*SEEKFACTOR
ddp.ddParameter = downValue
recompute()
downResult = ddg.ddResult
if (currResult < goal) == (upResult > goal): # found a range where a solution exists upside
return goalseek(baseValue, baseValue, itseek, itdicho, True, False)
if (currResult > goal) == (downResult < goal): # found a range where a solution exists downside
return goalseek(baseValue, baseValue, itseek, itdicho, False, True)
if (upResult > currResult) == (currResult > downResult) == (downResult > goal): # found a linear range, solution shoud be downside
return goalseek(downValue, currValue, itseek+1, itdicho, up, down)
if (upResult > currResult) == (currResult > downResult) == (upResult < goal): # found a linear range, solution shoud be upside
return goalseek(upValue, currValue, itseek+1, itdicho, up, down)
## found a concave/convexe range
upDiff = upResult - currResult
downDiff = currResult - downResult
SEEKFACTOR = SEEKFACTOR / 2
if ((currResult < goal) == (downDiff < 0)) != (abs(downDiff) > abs(upDiff)):
return goalseek(baseValue, currValue, itseek+1, itdicho, up, down)
else:
return goalseek(baseValue, currValue, itseek+1, itdicho, up, down)
## DICHO phase
if up:
currValue = currValue + baseValue*SEEKFACTOR*DICHOFACTOR**itdicho
if down:
currValue = currValue - baseValue*SEEKFACTOR*DICHOFACTOR**itdicho
ddp.ddParameter = currValue
recompute()
if abs(ddg.ddResult - goal) <= THRESHOLD:
print(f"""Found a solution with value : {currValue}
Deviation : {abs(ddg.ddResult - goal)}
Seek iterations : {itseek}
Dicho iterations : {itdicho}""")
return currValue
if (currResult < goal) == (ddg.ddResult < goal):
return goalseek(baseValue, currValue, itseek, itdicho+1, up, down)
else:
return goalseek(baseValue, currValue, itseek, itdicho+1, not up, not down)
def runSeek():
App.Console.PrintMessage(f"Seeking for goal : {goal}\n")
sb.setValue(0)
sb.show()
startVal = ddp.ddParameter
res = None
try:
res = goalseek(startVal)
except Exception:
App.Console.PrintError(f"Error : can't recompute with value : {ddp.ddParameter}\n")
if not res:
ddp.ddParameter = startVal # Restore start value if no solution found
sb.hide()
runSeek()
I didn't test it extensively but it should deal quite well with any polynomial (and maybe even periodic) transfer function.