Get Started

Solving optimal control problems with mesh refinement

This guide will help you setup the solve of optimal control problems with mesh refinement.


Overall structure

The overall structure for the main control file (main_MeshRefinement.m) shares many similarities with that of a single solve (see HERE). First we will fetch the formulated optimal control problem and the corrsponding options and solver settings.

clear all;close all;format compact;

global sol;  
sol=[];                             % Initialize solution structure

[problem,guess]=myProblem;  		% Fetch the problem definition

%%%% Choice among the ones belwo %%%%
% options= settings_auto(40);                  
% options= settings_h(40);                
% options= settings_hp(5,4);                 
options= settings_hp([4 5 3],[-1 0.3 0.4 1]);    % Get options and solver settings 

Then we can declare the variables for storage

% Declare variables for storage
errorHistory=zeros(2,length(problem.states.x0));
npsegmentHistory=zeros(2,1);
ConstraintErrorHistory=zeros(2,length(problem.constraintErrorTol));
timeHistory=zeros(1,2);
solutionHistory=cell(1,2);

and assign inital values to iteration control variables

maxAbsError=1e9; % the initial (before start) maximum absolute local eeror
i=1; % starting iteration
imax=50; % maximum number of mesh refinement iterations

The next step is to iteratively solve the OCP, analyze the errors and refine the mesh accordingly, until the absolute local errors and local constraint violation errors are all within the user defined tolerances.

while (any(maxAbsError > problem.states.xErrorTol) || any(maxAbsConstraintError > problem.constraintErrorTol)) && i<=imax    
    [infoNLP,data,options]=transcribeOCP(problem,guess,options); % Format for NLP solver
    [solution,infoNLP1,data] = solveNLP(infoNLP,data);      % Solve the NLP
    [solution]=output(problem,solution,options,data,4);         % Output solutions
    
    maxAbsError=max(abs(solution.Error));
    maxAbsConstraintError=max(solution.ConstraintError);
    errorHistory(i,:)=maxAbsError;
    ConstraintErrorHistory(i,:)=maxAbsConstraintError;
    timeHistory(i)=solution.computation_time;
    solutionHistory{i}=solution;
    
    if (any(maxAbsError > problem.states.xErrorTol) || any(maxAbsConstraintError > problem.constraintErrorTol)) && i <=imax
        [ options, guess ] = doMeshRefinement( options, problem, guess, data, solution, i );
    end
    i=i+1;
end

If wantted, the mesh refinement history can be organized and grouped with

MeshRefinementHistory.errorHistory=errorHistory;
MeshRefinementHistory.timeHistory=timeHistory;
MeshRefinementHistory.ConstraintErrorHistory=ConstraintErrorHistory;