Skip to content
Snippets Groups Projects
Commit a09e8d5f authored by François Henrotte's avatar François Henrotte
Browse files

mma and scipy as options on the same footing

parent f4f4c7cd
No related branches found
No related tags found
No related merge requests found
...@@ -387,7 +387,7 @@ PostOperation { ...@@ -387,7 +387,7 @@ PostOperation {
CreateDir["res"]; CreateDir["res"];
//Print[js, OnElementsOf Vol_S_Mag, File "res/js.pos"]; //Print[js, OnElementsOf Vol_S_Mag, File "res/js.pos"];
Print[b, OnElementsOf Vol_Mag, File "res/b.pos"]; //Print[b, OnElementsOf Vol_Mag, File "res/b.pos"];
Print[az, OnElementsOf Vol_Mag, File "res/az.pos"]; Print[az, OnElementsOf Vol_Mag, File "res/az.pos"];
Echo[ Str["k= PostProcessing.NbViews-1;", Echo[ Str["k= PostProcessing.NbViews-1;",
...@@ -501,10 +501,10 @@ PostOperation{ ...@@ -501,10 +501,10 @@ PostOperation{
File "tmp.geo", LastTimeStepOnly] ; File "tmp.geo", LastTimeStepOnly] ;
//Print[az, OnElementsOf Vol_Mag, File "res/azstar.pos"]; //Print[az, OnElementsOf Vol_Mag, File "res/azstar.pos"];
Echo[ Str["k= PostProcessing.NbViews-1;", // Echo[ Str["k= PostProcessing.NbViews-1;",
"View[k].IntervalsType = 3;", // "View[k].IntervalsType = 3;",
"View[0].NbIso = 30;"], // "View[0].NbIso = 30;"],
File "tmp.geo", LastTimeStepOnly] ; // File "tmp.geo", LastTimeStepOnly] ;
Print[ Lie_az, OnElementsOf Vol_Mag, File "res/Lie_az.pos" ]; Print[ Lie_az, OnElementsOf Vol_Mag, File "res/Lie_az.pos" ];
Echo[ Str["k= PostProcessing.NbViews-1;", Echo[ Str["k= PostProcessing.NbViews-1;",
......
...@@ -6,6 +6,10 @@ ParamIndex = ...@@ -6,6 +6,10 @@ ParamIndex =
Choices {0="Optimize on D and E", Choices {0="Optimize on D and E",
1="D: air gap length", 1="D: air gap length",
2="E: C-core thickness"} ]; 2="E: C-core thickness"} ];
Optimizer =
DefineNumber[0, Name "Optimization/Optimizer",
Choices {0="scipy",
1="MMA"}];
StepsInRange = StepsInRange =
DefineNumber[20, Name "Optimization/Steps in range", Visible ParamIndex]; DefineNumber[20, Name "Optimization/Steps in range", Visible ParamIndex];
......
...@@ -8,6 +8,7 @@ import conveks ...@@ -8,6 +8,7 @@ import conveks
import onelab import onelab
from scipy.optimize import minimize from scipy.optimize import minimize
modelName = 'ccore'
c = onelab.client(__file__) c = onelab.client(__file__)
...@@ -27,8 +28,6 @@ if(not len(mygetdp)): ...@@ -27,8 +28,6 @@ if(not len(mygetdp)):
c.sendInfo('Will use gmsh={0} and getdp={1}'.format(mygmsh, mygetdp)) c.sendInfo('Will use gmsh={0} and getdp={1}'.format(mygmsh, mygetdp))
modelName = 'ccore'
# append correct path to model file names # append correct path to model file names
file_geo = c.getPath(modelName+'.geo') file_geo = c.getPath(modelName+'.geo')
file_msh = c.getPath(modelName+'.msh') file_msh = c.getPath(modelName+'.msh')
...@@ -46,10 +45,6 @@ c.sendCommand('Solver.AutoMesh = -1;') ...@@ -46,10 +45,6 @@ c.sendCommand('Solver.AutoMesh = -1;')
# dry getdp run (without -solve or -pos option) to get model parameters in the GUI # dry getdp run (without -solve or -pos option) to get model parameters in the GUI
c.runSubClient('myGetDP', mygetdp) c.runSubClient('myGetDP', mygetdp)
# define optimization parameters as Onelab parameter (editable in the GUI)
maxIter = c.defineNumber('Optimization/00Max iterations', value=100)
maxChange = c.defineNumber('Optimization/01Max change', value=1e-0)
# end of check phase. We're done if we do not start the optimization # end of check phase. We're done if we do not start the optimization
if c.action == 'check' : if c.action == 'check' :
exit(0) exit(0)
...@@ -145,13 +140,22 @@ def test(n): ...@@ -145,13 +140,22 @@ def test(n):
# List of design parameters with format # List of design parameters with format
# [ "Onelab name", value, Min, Max, "variable name" ] # [ "Onelab name", value, Min, Max, "variable name" ]
designParameters={ designParameters={
0 : ['Optimization/Parameters/D', 0.01, 0.001, 0.050, 'D'], 0 : ['Optimization/Parameters/D', 0.02, 0.001, 0.050, 'D'],
1 : ['Optimization/Parameters/E', 0.01 , 0.01, 0.15 , 'E'] 1 : ['Optimization/Parameters/E', 0.01 , 0.01, 0.15 , 'E']
} }
Optimizers={
0 : 'scipy',
1 : 'mma'
}
x = {} x = {}
index = int(c.getNumber('Optimization/Sweep on parameter'))-1 index = int(c.getNumber('Optimization/Sweep on parameter'))-1
Optimizer = Optimizers[int(c.getNumber('Optimization/Optimizer'))]
# define optimization parameters as Onelab parameter (editable in the GUI)
maxIter = c.defineNumber('Optimization/00Max iterations', value=100, visible=(Optimizer=='mma'))
maxChange = c.defineNumber('Optimization/01Max change', value=1e-1, visible=(Optimizer=='mma'))
if index >= 0: # parameter sweep and optimisation over one selected design parameter if index >= 0: # parameter sweep and optimisation over one selected design parameter
x[0] = designParameters[index]; x[0] = designParameters[index];
...@@ -162,40 +166,28 @@ if index >= 0: # parameter sweep and optimisation over one selected design param ...@@ -162,40 +166,28 @@ if index >= 0: # parameter sweep and optimisation over one selected design param
print res print res
else: # optimisation over full parameter space else: # optimisation over full parameter space
if Optimizer=='scipy':
for label, var in designParameters.iteritems(): for label, var in designParameters.iteritems():
x[label] = var x[label] = var
values = [designParameters[i][1] for i in range(len(designParameters))] values = [designParameters[i][1] for i in range(len(designParameters))]
bounds = [(designParameters[i][2], designParameters[i][3]) for i in range(len(designParameters))] bounds = [(designParameters[i][2], designParameters[i][3]) for i in range(len(designParameters))]
## print x
## print values, bounds
## print fobj(values)
## print grad(values)
## exit(1)
res = minimize(fobj, values, jac=grad, bounds=bounds, callback=callbackF,\ res = minimize(fobj, values, jac=grad, bounds=bounds, callback=callbackF,\
method = 'L-BFGS-B', tol=None, options={'ftol': 1e-5, 'gtol': 1e-3, 'disp':True} ) method = 'L-BFGS-B', tol=None, options={'ftol': 1e-5, 'gtol': 1e-3, 'disp':True} )
print res print res
elif Optimizer=='mma':
numVariables = len(designParameters.keys())
exit(1)
# number of design variables and number of constraints
numVariables = len(x.keys())
# initial value of the variables in the design space,
# the lower bound for design variables,
# as well as the upper bound for design variables.
initialPoint = np.zeros(numVariables) initialPoint = np.zeros(numVariables)
lowerBound = np.zeros(numVariables) lowerBound = np.zeros(numVariables)
upperBound = np.zeros(numVariables) upperBound = np.zeros(numVariables)
for label, var in x.iteritems(): for label, var in designParameters.iteritems():
x[label] = var
initialPoint[label] = var[1] initialPoint[label] = var[1]
lowerBound[label] = var[2] lowerBound[label] = var[2]
upperBound[label] = var[3] upperBound[label] = var[3]
# Initialize the MMA optimizer # Initialize the MMA optimizer
conveks.mma.initialize(initialPoint, lowerBound, upperBound) conveks.mma.initialize(initialPoint, lowerBound, upperBound)
...@@ -212,17 +204,11 @@ conveks.mma.option.setNumber('SubProblem.asymptotesRmin', 0.001) # def 0.001 ...@@ -212,17 +204,11 @@ conveks.mma.option.setNumber('SubProblem.asymptotesRmin', 0.001) # def 0.001
# Get iteration count (here it will be 1 - could be different in case of restart) # Get iteration count (here it will be 1 - could be different in case of restart)
it = conveks.mma.getOuterIteration() it = conveks.mma.getOuterIteration()
change = 1.0 change = 1.0
while it <= maxIter and c.getString(c.name+'/Action') != 'stop' and change > 1e-4 :
while it <= maxIter and c.getString(c.name+'/Action') != 'stop': #and change > 1e-6 :
c.setNumber('Optimization/01Current iteration', value=it, readOnly=1, c.setNumber('Optimization/01Current iteration', value=it, readOnly=1,
attributes={'Highlight':'LightYellow'}) attributes={'Highlight':'LightYellow'})
# get (copy of) current point
xFromMMA = np.array(conveks.mma.getCurrentPoint()) xFromMMA = np.array(conveks.mma.getCurrentPoint())
print 'xFromMMA = ', xFromMMA
# get the value of the objective function and of the constraints
# as well as their sensitivity with respect to design variables at `x'
objective = fobj(xFromMMA) objective = fobj(xFromMMA)
grad_objective = grad(xFromMMA) grad_objective = grad(xFromMMA)
constraints = np.array([np.sum(xFromMMA)/100.0-1.0]) constraints = np.array([np.sum(xFromMMA)/100.0-1.0])
...@@ -232,7 +218,7 @@ while it <= maxIter and c.getString(c.name+'/Action') != 'stop': #and change > 1 ...@@ -232,7 +218,7 @@ while it <= maxIter and c.getString(c.name+'/Action') != 'stop': #and change > 1
objective *= fscale objective *= fscale
grad_objective *= fscale grad_objective *= fscale
print it, xFromMMA, objective, grad_objective print it, change, xFromMMA, objective, grad_objective
c.sendInfo('Optimization: it. {}, obj. {}, constr. {}'.format(it,objective,constraints[0])) c.sendInfo('Optimization: it. {}, obj. {}, constr. {}'.format(it,objective,constraints[0]))
# call MMA and update the current point # call MMA and update the current point
...@@ -244,3 +230,4 @@ while it <= maxIter and c.getString(c.name+'/Action') != 'stop': #and change > 1 ...@@ -244,3 +230,4 @@ while it <= maxIter and c.getString(c.name+'/Action') != 'stop': #and change > 1
# This should be called at the end # This should be called at the end
conveks.mma.finalize() conveks.mma.finalize()
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment