Candidate

resources/project/EEtUlUb-dLAdf0KpMVivaUlztwA/oXIaZ30PV4AuIRRG1GZvJ_SIyRMd.xml

@@ -0,0 +1,2 @@
+<?xml version='1.0' encoding='UTF-8'?>
+<Info Ref="src/core/TreatmentOptimization/costTerms" Type="Relative"/>

resources/project/EEtUlUb-dLAdf0KpMVivaUlztwA/oXIaZ30PV4AuIRRG1GZvJ_SIyRMp.xml

@@ -0,0 +1,2 @@
+<?xml version='1.0' encoding='UTF-8'?>
+<Info location="cc1969ad-586b-4b5c-8496-892291913346" type="Reference"/>

resources/project/NDPaR6DsgCrsBoJWF55J1vGaXlA/8fabh-WBmBGbDz7_h2NqYd4ovLcd.xml

@@ -0,0 +1,2 @@
+<?xml version='1.0' encoding='UTF-8'?>
+<Info/>

resources/project/NDPaR6DsgCrsBoJWF55J1vGaXlA/8fabh-WBmBGbDz7_h2NqYd4ovLcp.xml

@@ -0,0 +1,2 @@
+<?xml version='1.0' encoding='UTF-8'?>
+<Info location="1" type="DIR_SIGNIFIER"/>

resources/project/D7EGYpvIiO_nGXQuTOYjsSP_JQw/70qGohwUKdJepml9XrUSQBg8iXAd.xml → resources/project/NDPaR6DsgCrsBoJWF55J1vGaXlA/liHb9kbr54OOKjgYDq6G1l2UsfId.xml

@@ -1,6 +1,6 @@
 <?xml version='1.0' encoding='UTF-8'?>
 <Info>
     <Category UUID="FileClassCategory">
-        <Label UUID="derived"/>
+        <Label UUID="design"/>
     </Category>
 </Info>

resources/project/NDPaR6DsgCrsBoJWF55J1vGaXlA/liHb9kbr54OOKjgYDq6G1l2UsfIp.xml

@@ -0,0 +1,2 @@
+<?xml version='1.0' encoding='UTF-8'?>
+<Info location="isTrackingCostTerm.m" type="File"/>

resources/project/_2w7KV6EHc1SMHPJ5sv5BJ_TE_4/NDPaR6DsgCrsBoJWF55J1vGaXlAd.xml

@@ -0,0 +1,2 @@
+<?xml version='1.0' encoding='UTF-8'?>
+<Info/>

resources/project/_2w7KV6EHc1SMHPJ5sv5BJ_TE_4/NDPaR6DsgCrsBoJWF55J1vGaXlAp.xml

@@ -0,0 +1,2 @@
+<?xml version='1.0' encoding='UTF-8'?>
+<Info location="costTerms" type="File"/>

resources/project/D7EGYpvIiO_nGXQuTOYjsSP_JQw/jwQJUlUQLZW-X5QWOUQ_Ed6LSTUd.xml → resources/project/jwn0mMf3vUOgeAWeaJY8jqgb0nM/g0w2BX3hbHutyPMvNe1ZC1Hym4sd.xml

@@ -1,6 +1,6 @@
 <?xml version='1.0' encoding='UTF-8'?>
 <Info>
     <Category UUID="FileClassCategory">
-        <Label UUID="derived"/>
+        <Label UUID="design"/>
     </Category>
 </Info>

resources/project/jwn0mMf3vUOgeAWeaJY8jqgb0nM/g0w2BX3hbHutyPMvNe1ZC1Hym4sp.xml

@@ -0,0 +1,2 @@
+<?xml version='1.0' encoding='UTF-8'?>
+<Info location="writeNeuralControlPersonalizationOsimxFile.m" type="File"/>

src/DesignOptimization/DesignOptimization.m

@@ -1,7 +1,7 @@
 % This function is part of the NMSM Pipeline, see file for full license.
 %
 % () -> ()
-% 
+%
 
 % ----------------------------------------------------------------------- %
 % The NMSM Pipeline is a toolkit for model personalization and treatment  %
@@ -11,7 +11,7 @@
 % National Institutes of Health (R01 EB030520).                           %
 %                                                                         %
 % Copyright (c) 2021 Rice University and the Authors                      %
-% Author(s): Marleny Vega                                                 %
+% Author(s): Marleny Vega, Claire V. Hammond                              %
 %                                                                         %
 % Licensed under the Apache License, Version 2.0 (the "License");         %
 % you may not use this file except in compliance with the License.        %
@@ -26,19 +26,9 @@
 % ----------------------------------------------------------------------- %
 
 function [output, inputs] = DesignOptimization(inputs, params)
-pointKinematics(inputs.mexModel);
-inverseDynamics(inputs.mexModel);
-inputs = getStateDerivatives(inputs);
-inputs = setupGroundContact(inputs);
-inputs = getSplines(inputs);
-inputs = checkStateGuess(inputs);
-inputs = checkControlGuess(inputs);
-inputs = checkParameterGuess(inputs);
-inputs = getIntegralBounds(inputs);
-inputs = getPathConstraintBounds(inputs);
-inputs = getTerminalConstraintBounds(inputs); 
-inputs = getDesignVariableInputBounds(inputs);
-if strcmp(inputs.controllerType, 'synergy_driven') 
+inputs = makeTreatmentOptimizationInputs(inputs, params);
+initializeMexOrMatlabParallelFunctions(inputs.mexModel);
+if strcmp(inputs.controllerType, 'synergy_driven')
     inputs = setupMuscleSynergies(inputs);
 end
 output = computeDesignOptimizationMainFunction(inputs, params);
@@ -47,4 +37,4 @@
 inputs.splineSynergyActivations = spaps(inputs.initialGuess.time, ...
     inputs.initialGuess.control(:, inputs.numCoordinates + 1:end)', 0.0000001);
 inputs.synergyLabels = inputs.initialGuess.controlLabels(:, inputs.numCoordinates + 1:end);
-end
+end

src/DesignOptimization/DesignOptimizationTool.m

@@ -33,6 +33,6 @@ function DesignOptimizationTool(settingsFileName)
 settingsTree = xml2struct(settingsFileName);
 [inputs, params] = parseDesignOptimizationSettingsTree(settingsTree);
 [outputs, inputs] = DesignOptimization(inputs, params);
-reportDesignOptimizationResults(outputs, inputs);
+reportTreatmentOptimizationResults(outputs, inputs);
 saveDesignOptimizationResults(outputs, inputs);
-end
+end

src/DesignOptimization/calcDesignOptimizationDiscreteObjective.m

@@ -25,34 +25,11 @@
 % permissions and limitations under the License.                          %
 % ----------------------------------------------------------------------- %
 
-function discrete = calcDesignOptimizationDiscreteObjective(values, params)
+function discrete = calcDesignOptimizationDiscreteObjective(values, ...
+    modeledValues, auxdata)
 
-discrete = [];
-if isfield(params, 'discrete')
-    for i = 1:length(params.discrete.tracking)
-        costTerm = params.discrete.tracking{i};
-        if costTerm.isEnabled
-            switch costTerm.type
-                case "synergy_vectors"
-                    discrete = cat(2, discrete, ...
-                        calcTrackingCoordinateIntegrand(params, ...
-                        values.time, values.statePositions, ...
-                        costTerm.coordinate));   
-                otherwise
-                    throw(MException('', ['Cost term type ' costTerm.type ...
-                        ' does not exist for this tool.']))   
-            end
-        end
-    end
-    for i = 1:length(params.discrete.minimizing)
-        costTerm = params.discrete.minimizing{i};
-        if costTerm.isEnabled
-            switch costTerm.type
-                otherwise
-                    throw(MException('', ['Cost term type ' costTerm.type ...
-                        ' does not exist for this tool.']))   
-            end
-        end
-    end
-end
+[costTermCalculations, allowedTypes] = ...
+    generateCostTermStruct("discrete", "DesignOptimization");
+discrete = calcTreatmentOptimizationCost( ...
+    costTermCalculations, allowedTypes, values, modeledValues, auxdata);
 end

src/DesignOptimization/calcDesignOptimizationIntegrand.m

@@ -1,7 +1,7 @@
 % This function is part of the NMSM Pipeline, see file for full license.
 %
 % () -> ()
-% 
+%
 
 % ----------------------------------------------------------------------- %
 % The NMSM Pipeline is a toolkit for model personalization and treatment  %
@@ -11,7 +11,7 @@
 % National Institutes of Health (R01 EB030520).                           %
 %                                                                         %
 % Copyright (c) 2021 Rice University and the Authors                      %
-% Author(s): Marleny Vega                                                 %
+% Author(s): Marleny Vega, Claire V. Hammond                              %
 %                                                                         %
 % Licensed under the Apache License, Version 2.0 (the "License");         %
 % you may not use this file except in compliance with the License.        %
@@ -25,42 +25,12 @@
 % permissions and limitations under the License.                          %
 % ----------------------------------------------------------------------- %
 
-function integrand = calcDesignOptimizationIntegrand(values, params, ...
-    phaseout)
-integrand = [];
-for i = 1:length(params.integral.tracking)
-    costTerm = params.integral.tracking{i};
-    if costTerm.isEnabled
-        switch costTerm.type
-            case "coordinate"
-                integrand = cat(2, integrand, ...
-                    calcTrackingCoordinateIntegrand(params, ...
-                    values.time, values.statePositions, ...
-                    costTerm.coordinate));
-            case "controller"
-                integrand = cat(2, integrand, ...
-                    calcTrackingControllerIntegrand(params, values, ...
-                    costTerm.controller));      
-            otherwise
-                throw(MException('', ['Cost term type ' costTerm.type ...
-                    ' does not exist for this tool.']))   
-        end
-    end
-end
-for i = 1:length(params.integral.minimizing)
-    costTerm = params.integral.minimizing{i};
-    if costTerm.isEnabled
-        switch costTerm.type
-            case "joint_jerk"
-                integrand = cat(2, integrand, ...
-                    calcMinimizingJointJerkIntegrand(values.controlJerks, ...
-                    params, costTerm.coordinate));
-            otherwise
-                throw(MException('', ['Cost term type ' costTerm.type ...
-                    ' does not exist for this tool.']))   
-        end
-    end
-end
-integrand = scaleToBounds(integrand, params.maxIntegral, params.minIntegral);
+function integrand = calcDesignOptimizationIntegrand(values, ...
+    modeledValues, auxdata)
+[costTermCalculations, allowedTypes] = ...
+    generateCostTermStruct("continuous", "DesignOptimization");
+integrand = calcTreatmentOptimizationCost( ...
+    costTermCalculations, allowedTypes, values, modeledValues, auxdata);
+integrand = integrand ./ (auxdata.maxIntegral - auxdata.minIntegral);
 integrand = integrand .^ 2;
 end

src/DesignOptimization/calcDesignOptimizationObjective.m

@@ -25,8 +25,12 @@
 % permissions and limitations under the License.                          %
 % ----------------------------------------------------------------------- %
 
-function objective = calcDesignOptimizationObjective(discrete, continuous)
+function objective = calcDesignOptimizationObjective(discrete, ...
+    continuous, finalTime, inputs)
 continuousObjective = sum(continuous) / length(continuous);
+if isfield(inputs, "finalTimeRange")
+    continuousObjective = continuousObjective / finalTime;
+end
 discreteObjective = sum(discrete) / length(discrete);
 if isnan(discreteObjective); discreteObjective = 0; end
 objective = continuousObjective + discreteObjective;

src/DesignOptimization/computeDesignOptimizationContinuousFunction.m

@@ -1,7 +1,7 @@
 % This function is part of the NMSM Pipeline, see file for full license.
 %
 % () -> ()
-% 
+%
 
 % ----------------------------------------------------------------------- %
 % The NMSM Pipeline is a toolkit for model personalization and treatment  %
@@ -25,13 +25,20 @@
 % permissions and limitations under the License.                          %
 % ----------------------------------------------------------------------- %
 
-function phaseout = computeDesignOptimizationContinuousFunction(inputs)
+function modeledValues = computeDesignOptimizationContinuousFunction(inputs)
 
 values = getDesignOptimizationValueStruct(inputs.phase, inputs.auxdata);
-phaseout = calcTorqueBasedModeledValues(values, inputs.auxdata);
-phaseout = calcSynergyBasedModeledValues(values, inputs.auxdata, phaseout);
-phaseout.dynamics = calcDesignOptimizationDynamicsConstraint(values, inputs.auxdata);
-phaseout.path = calcDesignOptimizationPathConstraint(values, phaseout, inputs.auxdata);
-phaseout.integrand = calcDesignOptimizationIntegrand(values, inputs.auxdata, ...
-    phaseout);
-end
+inputs = updateSystemFromUserDefinedFunctions(inputs, values);
+modeledValues = calcTorqueBasedModeledValues(values, inputs.auxdata);
+modeledValues = calcSynergyBasedModeledValues(values, inputs.auxdata, ...
+    modeledValues);
+modeledValues.dynamics = calcDesignOptimizationDynamicsConstraint(values, ...
+    inputs.auxdata);
+if ~isempty(inputs.auxdata.path)
+    modeledValues.path = calcDesignOptimizationPathConstraint(values, ...
+        modeledValues, inputs.auxdata);
+end
+modeledValues.integrand = calcDesignOptimizationIntegrand(values, ...
+    modeledValues, inputs.auxdata);
+end
+

src/DesignOptimization/computeDesignOptimizationEndpointFunction.m

@@ -1,7 +1,7 @@
 % This function is part of the NMSM Pipeline, see file for full license.
 %
 % () -> ()
-% 
+%
 
 % ----------------------------------------------------------------------- %
 % The NMSM Pipeline is a toolkit for model personalization and treatment  %
@@ -11,7 +11,7 @@
 % National Institutes of Health (R01 EB030520).                           %
 %                                                                         %
 % Copyright (c) 2021 Rice University and the Authors                      %
-% Author(s): Marleny Vega                                                 %
+% Author(s): Marleny Vega, Claire V. Hammond                              %
 %                                                                         %
 % Licensed under the Apache License, Version 2.0 (the "License");         %
 % you may not use this file except in compliance with the License.        %
@@ -31,14 +31,35 @@
 inputs.phase.time = [inputs.phase.initialtime; inputs.phase.finaltime];
 inputs.phase.control = ones(size(inputs.phase.time,1), ...
     length(inputs.auxdata.minControl));
-values = getDesignOptimizationValueStruct(inputs.phase, inputs.auxdata);
+phase = inputs.phase;
+if isfield(inputs, "parameter")
+    phase.parameter = inputs.parameter;
+end
+values = getDesignOptimizationValueStruct(phase, inputs.auxdata);
+inputs = updateSystemFromUserDefinedFunctions(inputs, values);
 modeledValues = calcTorqueBasedModeledValues(values, inputs.auxdata);
 
 if ~isempty(inputs.auxdata.terminal)
-output.eventgroup.event = calcDesignOptimizationTerminalConstraint( ...
-   values, modeledValues, inputs.auxdata);
+    output.eventgroup.event = calcDesignOptimizationTerminalConstraint( ...
+        values, modeledValues, inputs.auxdata);
 end
-discrete = calcDesignOptimizationDiscreteObjective(values, inputs.auxdata);
+discrete = calcDesignOptimizationDiscreteObjective(values, ...
+    modeledValues, inputs.auxdata);
+% discrete = computeStaticParameterCost(inputs);
 output.objective = calcDesignOptimizationObjective(discrete, ...
-    inputs.phase.integral);
+    inputs.phase.integral, values.time(end), inputs.auxdata);
+end
+
+function cost = computeStaticParameterCost(inputs)
+costTerms = inputs.auxdata.costTerms;
+cost = 0;
+for i = 1:length(costTerms)
+    costTerm = costTerms{i};
+    if strcmp(costTerm.type, "user_defined") && ...
+            strcmp(costTerm.cost_term_type, "discrete")
+        func = str2func(costTerm.function_name);
+        cost = cost + func(inputs);
+    end
+end
 end
+