arbitrary objective function

[rgcoe]

I'm working through this... what I've done so far works but it's a bit clunky

closes #70

[rgcoe]

test_results.pdf

[rgcoe]

Having to create an RM3Device object is pretty clunky...

[H0R5E]

Ryan, I'm really sorry I haven't done this quickly enough, I'm happy for you to take it over if you want. It would be really helpful, in the future, if you could give me some rough deadlines for when you would like things like this done and then I can plan my upcoming work to suit (if possible). I think this will save us from duplicating effort and ensure the most efficient use of resources.

Regarding this implementation, I think you're sort of on the right track, but I would have tried to avoid exposing the implementation of the device model. My vision for this is that the user would build an "evaluation" class, with a single method, where the class has access to the "data transfer objects" (DTOs), so something like:

classdef MyCostFunc < WecOptTool.AbsEvaluation
   methods
      function cost = main(obj)
         % Motions is one of the DTOs (Struct like)
         cost = obj.motions.power 
      end
   end
end

Then once this is defined the user would add it to the study:

study.addObjective(MyCostFunc)

A lot of this rolls into what I was talking about yesterday, which I will write up shortly.

[rgcoe]

No worries (hope we weren't duplicating too much). This didn't take me too long and I am happy to hand it off if you have a good vision of how to be implement things.

On your suggestion above, I'm reluctant to force the user to implement OOP to accomplish this (I believe that's what you're suggesting?)... is there something easier?

[H0R5E]

Ha, I've just been duplicating grey hairs so far. :)

Another option is for the user to provide some kind of templated function of the form:

function result = myCost(envDev, loadings, motions)
    result = motions.power
end

In particular, the arguments of this function would have to be just right, and there would have to be some magic underneath to apply it, but it's doable and probably fits the basic interface of WecOptTool. This would be easier for the user to write, but it would be hard to debug if it were written incorrectly.

[rgcoe]

What if all calls to the performance simulation functions (PS, D, CC) return a single structure/object with a fixed set of fields? The idea would be that the returned structure/object has all the data needed to compute objective function values. If we want to be fancy, it could be an object that has methods for finding averages, maximums, even fatigue equivalent loadings. For example:

>> designRes = -1 * RM3Device.getPower(...); % <- probably want to rename, as we're not just getting power
>> disp(simRes)
       pow: [257×1 double]
       f_pto: [257×1 double]
       pos: [257×2 double]
       t: [257×1 double]
       cntrl: 'PS'
       m: 12345
       sa: 567
       designVars: [4x1 double]
fval = mObjectiveFun(desginRes)

We'd inevitably wind up adding fields to simRes, but at least we'd have a fixed argument structure when passing between functions.

[H0R5E]

Yes, this is basically the last DTO in #126. It needs some formalization around the structures (they must be well documented and reliable for the user), but it's essentially the same concept. You have also described abstracting the controllers, so that they are act in a polymorphic manner. Does this shed some light on #126 for you?

[H0R5E]

SosimRes is motions in #126 and in my comment about a format for defining a cost function above:

function result = myCost(envDev, loadings, motions)
result = motions.power
end

The problem is we only have one consistent output between the 3 controllers currently, so the user couldn't reliably build a cost function for that would work with any controller.

[rgcoe]

So simRes is motions in #126 and in my comment about a format for defining a cost function above:

Yes, I think you're right (motions ~= simRes)!

The problem is we only have one consistent output between the 3 controllers currently, so the user couldn't reliably build a cost function for that would work with any controller.

This is something we can resolve and need to resolve. Way back in the day (prior to 1c955cc), this was the case.

[H0R5E]

Yeah, I think there was some decision made that we would concentrate on power (@ssolon can confirm this) for the 0.1 release, and I don't think I ever saw any code that didn't do that.

I'm sure this is easy to hack together a solution for any cost function (which is fine for our paper), but I don't think it's that easy to offer an arbitrary, extendable solution without something like #126, in order to standardise what the user could select.

So, yeah, I guess the question is what does "arbitrary" mean?

[rgcoe]

Here's some pseudo-code to illustrate what I'm envisioning for the structure/workflow. It would seem to me that we would do the bulk of our work on WecOptTool.getDevice. I believe this basically represents some (minor?) restructuring of the existing code and more or less removing WecOptTool.run. The structure would be something like:

% same as existing
mySpectra = WecOptLib.tests.data.example8Spectra();

% for each device, you will have to create options where
% you define design variables and return an object of a given
% class. The device should inherit method from a superclass for
% evaluating performance with the 'P', 'CC', and 'PS' controlers
myDevice = WecOptTool.getDevice('RM3','CC','parametric')

% define a function handle
myObjFun = @(x) objFunWithStandardFormat(x, myDevice, mySpectra)

%% Set up some optimisation options for fmincon
options = optimoptions('fmincon');
options.MaxFunctionEvaluations = 5; % set artificial low for fast running
options.UseParallel = true;

% set constraints
lowerBounds = [...]
upperBounds = [...]
A = [...]
b = [...]
Aeq = [...]
beq = [...]
x0 = [...]

% solve problem
[sol,fval,exitflag,output] = fmincon(myObjFun,...
                					 x0,         	...
                					 A,          	...
                					 b,          	...
                					 Aeq,        	...
                					 beq,        	...
                					 lowerBounds,	...
                					 upperBounds,	...
                					 [],         	...
                					 options);

I've never made one of these before, but here's an attempt at a UML diagram.

[H0R5E]

EDIT: I fixed using the term "controller class" twice!

OK, I don't think we are too far apart, but I see something more general. I wouldn't abstract at device level and I see two possible approaches, one more verbose than the other.

Approach 1 (without a controller class):

cost = verboseCostFunction(param_values) <- takes the geometry values from the optimiser
    
    % S is the environmental data half of the environment/device DTO in #126
    S = WecOptTool.environment.fromTimeSeries("data.dat") <- User comes with some time series data

    % D is the device data half of the environment/device DTO in #126 
    geometry = WecOptTool.geometry.Gmsh("gmshfile", list_of_param)
    D = geometry.run(param_values)
    
    % F is the loadings DTO in #126
    solver = WecOptTool.solvers.NEMOH()  <- This an implementation of the abstract behaviour class in #126
    F = solver.run(S, D)
    
    % M is the motions DTO in #126
    controller = WecOptTool.control.PseudoSpectral() <- This an implementation of the abstract motions/control class in #126
    controller.addConstraints(constraints)
    M = controller.run(S, D, F)
    
    % Define a cost function
    cost = M.power / D.surface_area <- I just made these up, but the members would be documented for all the DTOs
    
end

% Call the optimiser
sol = Optimiser(verboseCostFunction,
                param_ranges)

% We would extra methods to recover results beyond the cost function
M_hist = controller.recall(sol)

Approach 2 (with an overlord class):

study = Study(some_options) <- There is an overarching overlord of the process
study.addTimeSeries("data.dat") <- User adds their time series data
study.addGeometry("gmshfile", list_of_param) <- User adds their geometry description
study.addSolver("NEMOH") <- Maybe could collect the class automagically

% Or we might need to build the class ourselves
controller = WecOptTool.control.PseudoSpectral()
controller.addConstraints(constraints)

study.addController(controller)

% User defines a cost function with this format
cost = shortCostFunction(S, D, F, M) <- Same DTOs as the first example
    cost = M.power / D.surface_area
end

study.addCost(shortCostFunction)

% Call the optimiser
Optimiser(study.run, <- the run method of the study joins everything together and record the history of the DTOs.
          param_ranges)

So, in some ways approach one is nice as it shows what is going on very clearly, on the other hand some of the passing of the DTOs is a little verbose once you get down to the controller. We could use the overlord class to store and recover the history of the outputs without an extra function. Also, if there were global options we wanted to set, this is easier with the overlord class, although you could have finer grained control of options with the first approach. Swings and roundabouts.

The second approach is an extension of the first, so perhaps just set the first approach up and then we can see if the second approach is necessary. What do you think?

[H0R5E]

Right, I optimised the first approach!

%% Setup

% S is the environmental data half of the environment/device DTO in #126
S = WecOptTool.environment.fromTimeSeries("data.dat") <- User comes with some time series data

geometry = WecOptTool.geometry.fromGmsh("gmshfile", list_of_param) <- Parametrise the geometry
solver = WecOptTool.solvers.NEMOH()  <- This an implementation of the abstract behaviour class in #126

controller = WecOptTool.control.PseudoSpectral() <- This an implementation of the abstract motions/control class in #126
controller.addConstraints(constraints)

%% Cost Function

cost = verboseCostFunction(param_values) <- takes the geometry values from the optimiser
    
    % D is the device data half of the environment/device DTO in #126 
    D = geometry.run(param_values)
    
    % F is the loadings DTO in #126
    F = solver.run(S, D)
    
    % M is the motions DTO in #126
    M = controller.run(S, D, F)
    
    % Define a cost function
    cost = M.power / D.surface_area <- I just made these up, but the members would be documented for all the DTOs
    
end

%% Optimise

sol = Optimiser(verboseCostFunction,
                param_ranges)

% We would need extra methods to recover results beyond the cost function
M_hist = controller.recall(sol)

[rgcoe]

Yes, I do think we are starting to converge to the same concept! As you might guess, I prefer "Approach 1," because it give the users a clearer view of what's happening (and perhaps more straightforward access to edit).

I think one challenge I have with your DTO structure in #126 and the code above is that it I don't like to think of running the geometry, solver and controller separately. Maybe its just me, but this somehow implies that they can exist independently of each other. I think I do get the intent (geometry.run produces a mesh and maybe volume, surface area, etc.).

Here's an updated diagram. This mostly aligns with your pseudo-code, with maybe some differences. Again, I'm sure my UML diagram is laughable, but hopefully gets the point across.

At this point, I'd say we either (A) make an actually correct UML diagram (not sure how to mix functions and classes...?) or (B) go ahead and try to write it in MATLAB (and return to update the UML diagram once we're done). I suppose I lean towards B at this point...

EDIT: add zip of StarUML file in case you want to edit
WecOptTool_Approach1.mdj.zip

[H0R5E]

Generally, I've not had great success with using UML class diagrams at this stage of a design because it just doesn't last into the not too distant future. There is a large burden in maintaining them as the design changes and they are not that useful on their own, I think. Nonetheless, I can try and formalise my concept for you this way, if it helps.

The "run" methods in my classes are because you can't call a class in MATLAB. In your geometry class "getNemoh" is your "run", but "getNemoh" doesn't work if were calling WAMIT, for instance. I think this is a key difference - I would have another class for using WAMIT (built from an abstract) but I'm guessing you would add a method called "getWAMIT". I worry that packing everything into one class will create something monolithic and hard to maintain.

Another important difference is that I am proposing to separate the data from the methods, (that's what the DTOs are), which would be represented by separate boxes in the diagram. You don't have to do this, you can just pass the objects around, but you won't get a truly modular architecture. Maybe it doesn't matter?

I also don't understand the purpose of the RM3 and wavebot classes here - what does the gear ratio have to do with the geometry, for instance? What is really being stored in here? Why do you need them? How will they be reused? It feels like a conflict, because you want some sort of superclass that does everything:

myDevice = WecOptTool.getDevice('RM3','CC','parametric')

but you also want a configurable system. If I'm a developer do I have to worry about creating a class for my device that can do all these things before I can use it, or can I bring a geometry file and use the tool immediately?

[H0R5E]

God, I'd kill for a whiteboard and a day of your time. This is so difficult to do in such a static environment.

So, I wonder if #66 is the difference between us here? How does the user bring arbitrary kinematics to the problem? When is what we have already done reusable and when does some recoding need to be done?

This is really where my domain knowledge lets me down, I'm afraid. I don't really understand what kind of data this is generating and because of that I can't picture the steps required to include it or how it changes the operation of the system.

[rgcoe]

yes, a whiteboard would be so nice - funny how the pandemic means that the Atlantic Ocean is really the determining factor in not being able to meet in person

the specific kinematics of each device is probably a key factor...

let's spend some time on our call tomorrow to work through this together

[H0R5E]

Regarding the issue with the run method in the polymorphic classes, we could implement something like this:

 [D, F, M] = param_values * geometry + solver + controller;

As before, the objects are defined outside of the loop. It's a little complicated to do, but this would be hidden in the abstract so that the concrete classes didn't have to worry about it.

[rgcoe]

I'm going to close this for now, as this effort has been taken over by #140.