Discrete-Valued Variables in Response Optimization (Code)
This example shows how to use response optimization to tune discrete-valued variables. Discrete-valued variables represent model parameters that are restricted to a finite set of values, instead of continuously varying. To use discrete-valued variables for response optimization or parameter estimation, you must use the surrogateopt solver for mixed-integer optimization.
Open the Model
In the sdoMotorPosition model, a PI controller enables the angular position of a DC motor to match a desired reference value. The load on the motor is subject to disturbances, and the controller must reject these disturbances.
open_system('sdoMotorPosition')
Within the Controller subsystem, the PI gains are set by a 12 V source divided down by resistors R1, R2, R3, and R4, as shown in the following diagram.
The proportional gain, Kp, is given by 12*R2/(R1+R2) and the integral gain, Ki, is given by 12*R4/(R3+R4). The initial values of the resistances are R1 = 47 kΩ, R2 = 180 kΩ, R3 = 10 kΩ, and R4 = 10 kΩ.
Specify Discrete Design Variables
Those four resistor values are the model parameters to tune for the optimization. Because resistors are only available in discrete values, use discrete parameters to represent them. To do so, use third input argument to sdo.getParameterFromModel to specify discrete parameters.
DesignVars = sdo.getParameterFromModel('sdoMotorPosition',[],{'R1','R2','R3','R4'});
sdo.getParameterFromModel returns an array of param.Discrete parameter objects that correspond to the model parameters in the order you specified. Set the ValueSet property of each parameter to the set of discrete values allowed for the parameter. Set the Value property to the current value, which must be present in ValueSet.
% R1 values DesignVars(1).ValueSet = [39 43 47 51 56] * 1e3; DesignVars(1).Value = 47*1e3; % R2 values DesignVars(2).ValueSet = [150 160 180 200 220] * 1e3; DesignVars(2).Value = 180*1e3; % R3 values DesignVars(3).ValueSet = [8.2 9.1 10 11 12] * 1e3; DesignVars(3).Value = 10*1e3; % R4 values DesignVars(4).ValueSet = [8.2 9.1 10 11 12] * 1e3; DesignVars(4).Value = 10*1e3;
Specify Design Requirements and Signals
The model applies a step disturbance at 1 second. With the initial resistance values specified in the model, the disturbance causes the motor to deviate by about 20°. The response then settles back to within ±5° of the reference position by 4 seconds after the disturbance. For this example, find new resistor values to improve this specification by 10%, so that the motor deviates no more than 18°, and settles back to within ±4.5° degrees of the reference position by 4 seconds after the disturbance. First, specify these new design requirements using sdo.requirements.SignalBound.
Requirements = struct; Requirements.UpperBound = sdo.requirements.SignalBound(... 'BoundMagnitudes', [18 18 ; 4.5 4.5], ... 'BoundTimes', [1 5 ; 5 15]); Requirements.LowerBound = sdo.requirements.SignalBound(... 'BoundMagnitudes', [-4.5 -4.5], ... 'BoundTimes', [5 15], ... 'Type', '>=');
To visualize the initial performance against the requirement, first specify signals to log during simulation. In particular, log the output of the Integrator block, which is the angular position of the model, the signal you want to optimize. To do so, create a simulation scenario with sdo.SimulationTest and configure its LoggingInfo property. Also, to prepare the scenario for later use in the optimization, assign the design variables you created as its Parameters property.
Simulator = sdo.SimulationTest('sdoMotorPosition'); Sig_Info = Simulink.SimulationData.SignalLoggingInfo; Sig_Info.BlockPath = 'sdoMotorPosition/Integrator'; Sig_Info.LoggingInfo.LoggingName = 'Sig'; Sig_Info.LoggingInfo.NameMode = 1; Simulator.LoggingInfo.Signals = Sig_Info; Simulator.Parameters = DesignVars;
Now create a plot of the disturbance-rejection requirements. Then, simulate the model and add the logged signal to the plot.
hRequirement = plot( ... Requirements.LowerBound.BoundTimes', ... Requirements.LowerBound.BoundMagnitudes', ... 'color','black'); hold on plot( ... Requirements.UpperBound.BoundTimes', ... Requirements.UpperBound.BoundMagnitudes', ... 'color','black'); % Simulator.Parameters = DesignVars; Simulator = sim(Simulator); SimLog = find(Simulator.LoggedData, ... get_param('sdoMotorPosition','SignalLoggingName')); Sig_Log = find(SimLog,'Sig'); clrs = colororder; hOriginal = plot(Sig_Log.Values,'Color',clrs(2,:)); legend([hRequirement hOriginal],'Requirement','Original Variables');
The plot shows that the current resistor values do not quite satisfy the more stringent requirement.
Set Up Optimization
The optimization process runs the model many times with different values for the resistor variables. To expedite this, configure the simulator for fast restart.
Simulator = setup(Simulator,FastRestart='on');Optimization requires an objective or cost function that runs the model and determines whether the disturbance rejection requirements are satisfied. For this example, use the function sdoMotorPosition_optFcn, which is included at the end of the example. This function is called many times by the optimization solver. To pass the function to the optimizer, use an anonymous function with one argument that calls sdoMotorPosition_optFcn.
optimfcn = @(P) sdoMotorPosition_optFcn(P,Simulator,Requirements);
Specify optimization options. Set the optimization method to surrogateopt solver, which is the method that supports tuning of discrete parameters.
Options = sdo.OptimizeOptions;
Options.Method = 'surrogateopt';
Options.MethodOptions.MaxFunctionEvaluations = 100;
Options.MethodOptions.ObjectiveLimit = 0.001;
Options.OptimizedModel = Simulator;Optimize the Design
To find resistance values optimized for the requirements, call sdo.optimize with the cost function handle, parameters to optimize, and options.
rng('default'); % for reproducibility [Optimized_DesignVars,Info] = sdo.optimize(optimfcn,DesignVars,Options);
Optimization started 2025-Aug-12, 18:41:23
Current
F-count max constraint max constraint Trial type
1 0.0759805 0.0759805 initial
2 0.0751287 0.0751287 random
3 0.0751287 0.179179 random
4 0.0447169 0.0447169 random
5 0.0317818 0.0317818 random
6 0.0317818 0.261353 random
7 0.0257711 0.0257711 random
8 0.0257711 0.0930365 random
9 0.0257711 0.0971938 random
10 0.0257711 0.200359 random
11 0.0062997 0.0062997 random
12 0.0062997 0.206374 random
13 0.0062997 0.0818387 random
14 0.0062997 0.0585416 random
15 0.0062997 0.114002 random
16 0.0062997 0.112485 random
17 0.0062997 0.0759805 random
18 0.0062997 0.0709489 random
19 0.0062997 0.280136 random
20 0.0062997 0.099985 random
21 -0.00759693 -0.00759693 adaptive
22 -0.00759693 0.0691728 random
23 -0.00759693 0.0347164 random
24 -0.00759693 0.0977682 random
25 -0.00759693 0.0589627 random
26 -0.00759693 0.131778 random
27 -0.00759693 0.147347 random
28 -0.00759693 0.104488 random
29 -0.00759693 0.0534262 random
30 -0.00759693 0.1643 random
31 -0.00759693 0.0419731 random
32 -0.00759693 0.0396579 random
33 -0.00759693 0.077945 random
34 -0.00759693 0.0883493 random
35 -0.00762259 -0.00762259 random
36 -0.00762259 0.0455853 random
37 -0.00762259 0.13976 random
38 -0.00762259 0.140634 random
39 -0.00762259 0.0662368 random
40 -0.00762259 0.256738 random
41 -0.00762259 0.0348764 random
42 -0.00762259 0.18995 random
43 -0.00762259 0.0615952 random
44 -0.00762259 0.0597358 random
45 -0.00762259 0.015983 random
46 -0.00762259 0.0770592 random
47 -0.00762259 0.174666 random
48 -0.00762259 0.081758 random
49 -0.00762259 0.0240828 random
50 -0.00762259 0.173292 random
Current
F-count max constraint max constraint Trial type
51 -0.00762259 0.116519 random
52 -0.00762259 0.217099 random
53 -0.00762259 0.104371 random
54 -0.00762259 0.0937562 random
55 -0.00762259 0.104225 random
56 -0.00762259 0.0380129 random
57 -0.00762259 0.150819 random
58 -0.00762259 0.0636666 random
59 -0.00762259 0.0320499 random
60 -0.00762259 0.0273158 random
61 -0.00762259 0.0798937 random
62 -0.00762259 0.174146 random
63 -0.00762259 0.0166127 random
64 -0.00762259 0.00433873 adaptive
65 -0.00762259 0.294106 random
66 -0.00518218 -0.00518218 adaptive
67 -0.00518218 0.0319864 random
68 -0.00518218 0.266843 random
69 -0.00518218 0.0342618 random
70 -0.00518218 0.0982987 random
71 -0.00518218 0.0419731 random
72 -0.00518218 0.0396579 random
73 -0.00518218 0.299815 random
74 -0.00518218 0.0712971 random
75 -0.00518218 0.0709489 random
76 -0.00518218 0.0589627 random
77 -0.00518218 0.131778 random
78 -0.00518218 0.00653156 random
79 -0.00518218 0.211877 random
80 -0.00518218 0.0967155 random
81 -0.00518218 0.0689487 random
82 -0.00518218 0.0257711 random
83 -0.00518218 0.0930365 random
84 -0.00518218 0.114002 random
85 -0.00518218 0.112485 random
86 -0.00518218 0.06706 random
87 -0.00518218 0.0751287 random
88 -0.00737024 -0.00737024 adaptive
89 -0.00737024 0.0759805 random
90 -0.00737024 0.10447 random
91 -0.00737024 0.190644 random
92 -0.00737024 0.142755 random
93 -0.00737024 0.104193 random
94 -0.00737024 0.100237 random
95 -0.00737024 0.097787 random
96 -0.00737024 0.0735855 random
97 -0.00737024 0.137557 random
98 -0.00737024 0.0783709 random
99 -0.00737024 0.0875674 random
100 -0.00737024 0.0794652 random
101 -0.00737024 -0.00737024 best value
The current solution is feasible. surrogateopt stopped because it exceeded the function evaluation limit set by the 'MethodOptions.MaxFunctionEvaluations' property in the sdo.OptimizeOptions object.
If the solution needs to be improved, you could try increasing the function evaluation limit.
Evaluate Optimized Design
Plot the model response after optimization.
Simulator.Parameters = Optimized_DesignVars; Simulator = sim(Simulator); SimLog = find(Simulator.LoggedData, ... get_param('sdoMotorPosition','SignalLoggingName')); Sig_Log = find(SimLog,'Sig'); clrs = colororder; hOptimized = plot(Sig_Log.Values, 'Color', clrs(1,:)); legend([hRequirement hOriginal hOptimized], ... 'Requirement','Original Variables','Optimized Variables');
Now the motor position is within the spec of the disturbance rejection requirement.
Update Model with Optimized Parameter Values
sdo.Optimize returns tuned versions of the parameters in the array Optimized_DesignVars. Each entry in this array is a param.Discrete parameter object whose Value property is set to the optimized value. For instance, examine the new value of R2.
Optimized_DesignVars(2).Value
ans = 220000
Update the model with the optimized parameter values.
sdo.setValueInModel('sdoMotorPosition',Optimized_DesignVars);Restore the simulator fast restart settings.
Simulator = restore(Simulator);
Conclusion
In this example you tuned variables to improve the disturbance rejection characteristics of a motor controller. The variables being tuned were electrical component parts that could only take on discrete values, rather than continuous values. You used sdo.getParameterFromModel to specify the variables as discrete, and used the surrogateopt solver to tune the parameters.
Objective Function
The function sdoMotorPosition_optFcn is called at each iteration of the optimization problem with a set of parameter values, P. The function returns the objective value and constraint violations, Vals, to the optimization solver. See sdo.optimize for a more detailed description of the function signature.
function Vals = sdoMotorPosition_optFcn(P,Simulator,Requirements) %SDOMOTORPOSITION_OPTFCN % % Simulate the model. Simulator.Parameters = P; Simulator = sim(Simulator); % Retrieve logged signal data. SimLog = find(Simulator.LoggedData, ... get_param('sdoMotorPosition','SignalLoggingName')); Sig_Log = find(SimLog,'Sig'); % Evaluate the design requirements. Cleq_UpperBound = evalRequirement( ... Requirements.UpperBound,Sig_Log.Values); Cleq_LowerBound = evalRequirement( ... Requirements.LowerBound,Sig_Log.Values); % Collect the evaluated design requirement values in a structure to % return to the optimization solver. Vals.Cleq = [... Cleq_UpperBound(:); ... Cleq_LowerBound(:)]; end
See Also
sdo.optimize | param.Discrete | sdo.getParameterFromModel | sdo.OptimizeOptions
