Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Default sprectral weighting #137

Closed
wants to merge 5 commits into from
Closed

Default sprectral weighting #137

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
2123ebc
set default spectral weight for example
rgcoe May 5, 2020
8b136bc
use deal instead of for loop to set default weighting
rgcoe May 5, 2020
8153a10
update line numbers to correspond w/ example.m
rgcoe May 5, 2020
cf7747c
store mu in mat file instead of adding in 'example8Spectra.m'
rgcoe May 7, 2020
705da40
check for weighting in checkSpectrum during creation of study object
rgcoe May 7, 2020
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
32 changes: 16 additions & 16 deletions docs/user/example.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -48,9 +48,9 @@ passed to the main functions of the toolbox.

.. literalinclude:: /../example.m
:language: matlab
:lines: 24-26
:lines: 21-23
:linenos:
:lineno-start: 24
:lineno-start: 21

Define a Sea-State
==================
Expand All @@ -62,9 +62,9 @@ toolbox or preset spectra from WecOptTool.

.. literalinclude:: /../example.m
:language: matlab
:lines: 27-35
:lines: 24-32
:linenos:
:lineno-start: 27
:lineno-start: 24

Spectra are formatted following the convention of the WAFO_ MATLAB toolbox, but
can be generated in via any means (e.g., from buoy measurements) as long as the
Expand Down Expand Up @@ -107,9 +107,9 @@ The desired spectrum or spectra can then be added to the study object.

.. literalinclude:: /../example.m
:language: matlab
:lines: 36-38
:lines: 33-35
:linenos:
:lineno-start: 36
:lineno-start: 33

Add a controller to the study
=============================
Expand Down Expand Up @@ -153,9 +153,9 @@ sub-package.

.. literalinclude:: /../example.m
:language: matlab
:lines: 39-42
:lines: 36-39
:linenos:
:lineno-start: 39
:lineno-start: 36

Define design variables
=======================
Expand All @@ -181,19 +181,19 @@ The initial values,``x0``, lower bounds, ``lb``, and upper bounds,

.. literalinclude:: /../example.m
:language: matlab
:lines: 43-50
:lines: 40-47
:linenos:
:lineno-start: 43
:lineno-start: 40

Alternatively, a simpler study with a single scalar design variable can be
employed. In this case, instead of scaling various dimensions of the device
individually, the entire device is scaled based on a single design variable.

.. literalinclude:: /../example.m
:language: matlab
:lines: 51-58
:lines: 48-57
:linenos:
:lineno-start: 51
:lineno-start: 48

The options for design variables are defined as classes in the
:mat:mod:`~+WecOptTool.+geom` sub-package.
Expand All @@ -216,9 +216,9 @@ MATLAB's ``fmincon`` optimization solver is used in |example.m|_.

.. literalinclude:: /../example.m
:language: matlab
:lines: 59-63
:lines: 56-60
:linenos:
:lineno-start: 59
:lineno-start: 56

Run the study and view results
==============================
Expand All @@ -229,9 +229,9 @@ follows:

.. literalinclude:: /../example.m
:language: matlab
:lines: 64-72
:lines: 61-68
:linenos:
:lineno-start: 64
:lineno-start: 61

From :mat:func:`~+WecOptTool.result`, we can see that the study has produced the
following design.
Expand Down
6 changes: 0 additions & 6 deletions example.m
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -18,9 +18,6 @@
% You should have received a copy of the GNU General Public License
% along with WecOptTool. If not, see <https://www.gnu.org/licenses/>.

%% Suppress warnings
warning('off', 'WaveSpectra:NoWeighting')

%% Create an RM3 study object
study = WecOptTool.RM3Study();

Expand Down Expand Up @@ -69,6 +66,3 @@

%% Plot the results
WecOptTool.plot(study);

%% Re-enable warnings
warning('on', 'WaveSpectra:NoWeighting')
37 changes: 0 additions & 37 deletions toolbox/+WecOptLib/+models/DeviceModelTemplate.m
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -125,43 +125,6 @@
powSSs = zeros(NSS, 1);
powPerFreqs = cell(NSS);
freqs = cell(NSS);
n_mu = 0;

% Check sea-state weights
for iSS = 1:NSS

S = SS(iSS);

if isfield(S, 'mu')
n_mu = n_mu + 1;
end

end

if NSS == 1

% Single sea-state requires no weighting
SS(1).mu = 1;

elseif n_mu == 0

% Equalise weightings for multi-sea-states if not given
for iSS = 1:NSS
SS(iSS).mu = 1;
end

warn = ['Provided wave spectra have no weightings ' ...
'(field mu). Equal weighting presumed.'];
warning('WaveSpectra:NoWeighting', warn);

elseif n_mu ~= NSS

% Don't allow partial weightings
msg = ['Weighting field mu must be set for all spectra ' ...
'or for none of them.'];
error(msg)

end

% Iterate over Sea-States
for iSS = 1:NSS % TODO - consider parfor?
Expand Down
Binary file modified toolbox/+WecOptLib/+tests/+data/8spectra.mat
View file
Open in desktop
Binary file not shown.
3 changes: 3 additions & 0 deletions toolbox/+WecOptLib/+tests/+data/example8Spectra.m
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -21,9 +21,12 @@
function SS = example8Spectra
%EXAMPLE8SPECTRA Example Bretschneider spectrum with varying HHm0s, Tps,
% Nbins, and range

% load data
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is a superfluous comment. Remove or edit the docstring, if necessary.

p = mfilename('fullpath');
[filepath, ~, ~] = fileparts(p);
dataPath = fullfile(filepath, '8spectra.mat');
example_data = load(dataPath);
SS = example_data.SS;

end
Binary file modified toolbox/+WecOptLib/+tests/+data/spectrum.mat
View file
Open in desktop
Binary file not shown.
14 changes: 14 additions & 0 deletions toolbox/+WecOptLib/+tests/+unit/checkSpectrumTest.m
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -58,3 +58,17 @@ function testCheck_Positive(testCase)
eID = 'WecOptTool:invalidSpectrum:negativeFrequencies';
verifyError(testCase,@() WecOptLib.utils.checkSpectrum(S),eID)
end

function testCheck_missingWeights(testCase)
S = WecOptLib.tests.data.example8Spectra();
S = rmfield(S,'mu');
wID = 'WaveSpectra:NoWeighting';
verifyWarning(testCase,@() WecOptLib.utils.checkSpectrum(S),wID)
end

function testCheck_WeightingsInconsistent(testCase)
S = WecOptLib.tests.data.example8Spectra();
S(1).mu = [];
eID = 'WecOptTool:invalidSpectrum:invalidWeightings';
verifyError(testCase,@() WecOptLib.utils.checkSpectrum(S),eID)
end
163 changes: 93 additions & 70 deletions toolbox/+WecOptLib/+utils/checkSpectrum.m
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,100 +1,123 @@

% Copyright 2020 National Technology & Engineering Solutions of Sandia,
% LLC (NTESS). Under the terms of Contract DE-NA0003525 with NTESS, the
% Copyright 2020 National Technology & Engineering Solutions of Sandia,
% LLC (NTESS). Under the terms of Contract DE-NA0003525 with NTESS, the
% U.S. Government retains certain rights in this software.
%
% This file is part of WecOptTool.
%
%
% WecOptTool is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
%
% WecOptTool is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
% GNU General Public License for more details.
%
%
% You should have received a copy of the GNU General Public License
% along with WecOptTool. If not, see <https://www.gnu.org/licenses/>.

function [] = checkSpectrum(S)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'm not sure if checkSpectrum is the correct name for this now if it does two jobs, checking and fixing.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

agreed, but I actually do prefer that its just one function call to do both

Copy link
Member

@H0R5E H0R5E May 15, 2020 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think I am going to push you on this, but not just because it's nasty to have one function with two purposes (ick). I'm most worried that, as we change to the new structure, there is a risk that this function ends up back inside the optimisation loop without the spectra being corrected first, rendering the original purpose of stopping the warning showing up in every iteration moot.

If we have an explicit function that fixes the mu prior to use in any optimisation that the user creates, then checkSpectrum can simply generate an error if it's not set (for more than one spectrum), which is its typical behaviour, allowing it to be safely used anywhere in the code. Does that make sense?

% checkSprectrum(S)
%
% Checks whether the input S is a valid spectrum structure (following
% WAFO).
%
% Inputs
% S spectrum structure (can be arrary) in the style of WAFO
% with the fields:
% S.w column vector of frequencies in [rad/s]
% S.S column vector of spectral density in [m^2 rad/ s]
%
% Outputs
% (none) (will throw error if not valid)
%
%
% Example
% Hm0 = 5;
% Tp = 8;
% S = bretschneider([],[Hm0,Tp]);
% WecOptLib.utils.checkSpectrum(S)
%
% -------------------------------------------------------------------------


inds = 1:length(S);

try
arrayfun(@(Spect,idx) checkFields(Spect, idx), S, inds);
arrayfun(@(Spect,idx) checkLengths(Spect, idx), S, inds);
arrayfun(@(Spect,idx) checkCol(Spect, idx), S, inds);
arrayfun(@(Spect,idx) checkPositive(Spect, idx), S, inds);
catch MEs
throw(MEs)
end

% checkSprectrum(S)
%
% Checks whether the input S is a valid spectrum structure (following
% WAFO).
%
% Inputs
% S spectrum structure (can be arrary) in the style of WAFO
% with the fields:
% S.w column vector of frequencies in [rad/s]
% S.S column vector of spectral density in [m^2 rad/ s]
%
% Outputs
% (none) (will throw error if not valid)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Would need changed if S is returned

%
%
% Example
% Hm0 = 5;
% Tp = 8;
% S = bretschneider([],[Hm0,Tp]);
% WecOptLib.utils.checkSpectrum(S)
%
% -------------------------------------------------------------------------


inds = 1:length(S);

% check if there are no weightings, if necessary populate and warn
S = checkWeightingsPresent(S);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

For this to be useful S must be returned by checkSpectrum


try
arrayfun(@(Spect,idx) checkFields(Spect, idx), S, inds);
arrayfun(@(Spect,idx) checkLengths(Spect, idx), S, inds);
arrayfun(@(Spect,idx) checkCol(Spect, idx), S, inds);
arrayfun(@(Spect,idx) checkPositive(Spect, idx), S, inds);
checkWeightingsConsistent(S)
catch MEs
throw(MEs)
end
end

function [] = checkFields(S, idx)
fns = {'S','w'};
msg = 'Spectrum #%i in array does not contain required S.S and S.w fields';
ID = 'WecOptTool:invalidSpectrum:missingFields';
try
assert(all(isfield(S,fns)),ID,msg,idx)
catch ME
throw(ME)
end
fns = {'S','w'};
msg = 'Spectrum #%i in array does not contain required S.S and S.w fields';
ID = 'WecOptTool:invalidSpectrum:missingFields';
try
assert(all(isfield(S,fns)),ID,msg,idx)
catch ME
throw(ME)
end
end

function [] = checkLengths(S, idx)
msg = 'Spectrum #%i in array S.S and S.w fields are not same length';
ID = 'WecOptTool:invalidSpectrum:mismatchedLengths';
try
assert(length(S.S) == length(S.w),ID,msg, idx)
catch ME
throw(ME)
end
msg = 'Spectrum #%i in array S.S and S.w fields are not same length';
ID = 'WecOptTool:invalidSpectrum:mismatchedLengths';
try
assert(length(S.S) == length(S.w),ID,msg, idx)
catch ME
throw(ME)
end
end

function [] = checkCol(S, idx)
msg = 'Spectrum #%i in array S.S and S.w fields are not column vectors';
ID = 'WecOptTool:invalidSpectrum:notColumnVectors';
try
assert(iscolumn(S.S) && iscolumn(S.w),ID,msg, idx)
catch ME
throw(ME)
msg = 'Spectrum #%i in array S.S and S.w fields are not column vectors';
ID = 'WecOptTool:invalidSpectrum:notColumnVectors';
try
assert(iscolumn(S.S) && iscolumn(S.w),ID,msg, idx)
catch ME
throw(ME)
end
end
end


function [] = checkPositive(S, idx)
msg = ['Frequency in Spectrum #%i contains negative values. Frequency'...
' values must be positive'];
ID = 'WecOptTool:invalidSpectrum:negativeFrequencies';
try
assert(all(S.w >=0) ,ID,msg, idx)
catch ME
throw(ME)
msg = ['Frequency in Spectrum #%i contains negative values. Frequency'...
' values must be positive'];
ID = 'WecOptTool:invalidSpectrum:negativeFrequencies';
try
assert(all(S.w >=0) ,ID,msg, idx)
catch ME
throw(ME)
end
end

function [] = checkWeightingsConsistent(S)
emsg = ['Some spectra are missing the weighting field mu'];
eID = 'WecOptTool:invalidSpectrum:invalidWeightings';
try
n_mu = sum(arrayfun(@(x) ~isempty(x.mu), S));
n_S = length(S);
assert(n_mu == n_S, eID, emsg);
catch ME
throw(ME)
end
end

function [S] = checkWeightingsPresent(S)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Consider making this a standalone function, then checkSpectrum doesn't need to return S and it's purpose remains clear. I think the name could be changed as well, to something like "fixWeightings", just to show that it changes the input.

wmsg = ['Missing spetra weighting fields, set to equal weighting'];
wID = 'WaveSpectra:NoWeighting';
if all(arrayfun(@(x) ~isfield(x,'mu'), S))
warning(wID,wmsg)
[S.mu] = deal(1);
end
end
Loading