Skip to content

correct negative damping in BEM results#134

Merged
rgcoe merged 35 commits intoSNL-WaterPower:masterfrom
rgcoe:checkBemIrregFreqs
Sep 15, 2020
Merged

correct negative damping in BEM results#134
rgcoe merged 35 commits intoSNL-WaterPower:masterfrom
rgcoe:checkBemIrregFreqs

Conversation

@rgcoe
Copy link
Member

@rgcoe rgcoe commented May 1, 2020

Description

Nemoh is well known to have issues with irregular frequencies (see discussion in #35). This can create cases where there is negative radiation damping at certain frequencies (which is not physical for diagonal elements). If the sum of the radiation damping and friction are less than zero, this can create stability issues and is believed to be the source of the convergence problems noted in #4.

This PR adds a function (checkNemoh) to look for negative values in the diagonal radiation damping terms in the hydro structure produced by getNemoh. checkNemoh sets values less than 0 to 0 and returns a tracking of where negative values were found.

Screen Shot 2020-05-01 at 11 20 46 AM

Before
----------
    0.0523    2.7783   15.8544   23.5220   20.8416   16.0492   11.6885    8.3137    1.3735    3.6513    2.6734
    0.0523    2.7783   15.8543   23.5220   20.8416   16.0492   11.6885    8.3137    1.3735    3.6513    2.6734
   58.1025  103.7792   82.5054   46.8996   21.5337    9.7066    0.7951    0.5628    0.0681   -0.0089   -0.0419
    0.8039   32.5602  121.1623  102.2461   45.1301   13.4254    1.9090    1.3394  -14.7975    0.6399    0.5099
    0.8039   32.5602  121.1622  102.2458   45.1301   13.4250    1.9093    1.3394  -14.7976    0.6399    0.5099
   -0.0000   -0.0000    0.0001   -0.0005    0.0009   -0.0006   -0.0001    0.0002   -0.0002    0.0000    0.0003
    0.0264   -0.0025    0.0000   -0.0000   -0.0000   -0.0000    0.0000    0.0000   -0.0000   -0.0000   -0.0000
    0.0264   -0.0025    0.0000   -0.0000   -0.0000   -0.0000    0.0000    0.0000   -0.0000   -0.0000   -0.0000
    2.0205    0.4124    0.0026   -0.0013   -0.0017   -0.0053   -0.0014   -0.0001    0.0000   -0.0000   -0.0000
   51.0431   -3.9977    0.0049   -0.0001   -0.0001   -0.0002    0.0000    0.0001   -0.0005   -0.0000   -0.0000
   51.0555   -3.9983    0.0049   -0.0001   -0.0001   -0.0002    0.0000    0.0001   -0.0005   -0.0000   -0.0000
    0.0000   -0.0000    0.0000   -0.0000   -0.0000   -0.0000   -0.0000   -0.0000   -0.0000    0.0000    0.0000


After
---------
    0.0523    2.7783   15.8544   23.5220   20.8416   16.0492   11.6885    8.3137    1.3735    3.6513    2.6734
    0.0523    2.7783   15.8543   23.5220   20.8416   16.0492   11.6885    8.3137    1.3735    3.6513    2.6734
   58.1025  103.7792   82.5054   46.8996   21.5337    9.7066    0.7951    0.5628    0.0681         0         0
    0.8039   32.5602  121.1623  102.2461   45.1301   13.4254    1.9090    1.3394         0    0.6399    0.5099
    0.8039   32.5602  121.1622  102.2458   45.1301   13.4250    1.9093    1.3394         0    0.6399    0.5099
         0         0    0.0001         0    0.0009         0         0    0.0002         0    0.0000    0.0003
    0.0264         0    0.0000         0         0         0    0.0000    0.0000         0         0         0
    0.0264         0    0.0000         0         0         0    0.0000    0.0000         0         0         0
    2.0205    0.4124    0.0026         0         0         0         0         0    0.0000         0         0
   51.0431         0    0.0049         0         0         0    0.0000    0.0001         0         0         0
   51.0555         0    0.0049         0         0         0    0.0000    0.0001         0         0         0
    0.0000         0    0.0000         0         0         0         0         0         0    0.0000    0.0000

Fixes #4

test_results.pdf

Checklist:

  • Added the results of running the test suite (in pdf form)
  • All new files contain the GPL header
  • If example.m has been modified, the content / line numbers in
    docs\user\example.rst are still valid or have been fixed

Ryan Coe added 2 commits May 1, 2020 11:08
currently checks for negative vals in diagonal FRFs
may need to further consider whether/how to pass a warning up to higher
level functions and/or set debug/verbosity level across WecOptTool
@rgcoe rgcoe added bug Something isn't working enhancement New feature or request labels May 1, 2020
'flagged for negative values'])
end

function testNegativeRadDampingKnownBad(testCase)
Copy link
Member Author

Choose a reason for hiding this comment

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

One thing that I have not yet rigorously checked is whether this PR resolves the convergence issues originally noted in #4. My own naivety is the limiting factor. We need to write a unit test case (or some snippet that we just call to test this PR ad-hoc) that does the following:

  1. Tries to find the power for known nonconvergent cases mentioned in noncoverged PS cases? #4 (RM3_getPow(S, 'PS', 'parametric', [14.76, 15.71, 2.58, 37.27]))
  2. Check that the PS problem does NOT converge when checkNemoh is not used (note that the application of checkNemoh is currently hardcoded in getNemoh)
  3. Check that the PS problem does converged when checkNemoh is used.

@H0R5E - I would do this, but I don't quite know how yet with the toolbox/OOP structure, I can certainly finish it through if you'll get me started or point to a relevant example.

Copy link
Member

Choose a reason for hiding this comment

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

I think the best way to deal with this is to break up getNemoh, so that the call to WEC_Sim.Read_NEMOH is in a different function (getHydro?), and remove the call to checkNemoh in there also. Then the assumption is that getNemoh (or whatever is a sensible name for the rest of the function) will take a checked (or not) WEC-SIM hydro object. Then you can test with or without checking.

I can do this for you, if you want, but I don't have push access to your fork and I don't think you clicked the allow edits from maintainers box, so you would need to grant me rights on your fork.

@rgcoe rgcoe requested a review from H0R5E May 1, 2020 17:35
@rgcoe rgcoe mentioned this pull request May 4, 2020
@H0R5E
Copy link
Member

H0R5E commented May 5, 2020

I'm going to switch this to a draft while we're updating it.

@H0R5E H0R5E marked this pull request as draft May 5, 2020 13:19
H0R5E and others added 4 commits May 5, 2020 15:24
This is done through the parametric geometry option, which takes
an option called fixNegative, that defaults to true.
…to checkBemIrregFreqs

# Conflicts:
#	toolbox/+WecOptLib/+nemoh/getNemoh.m
@H0R5E
Copy link
Member

H0R5E commented May 5, 2020

OK so this does not appear to fix #4. I added two tests to +WecOptLib\+tests\+unit\RM3DeviceModelTest.m called testPSNoConverge and testPSConverge where one uses the damping fix and the other doesn't and testPSConverge errors. Interestingly the first test in that file also fails.

Back to you!

@H0R5E
Copy link
Member

H0R5E commented May 5, 2020

Hang on, that's probably my fault! I didn't test the exit flag properly. One sec.

@H0R5E
Copy link
Member

H0R5E commented May 5, 2020

OK, so after I check the value of the exitflag (doh), it seems like your expected failure case doesn't actually fail. Have I done something else wrong? Maybe it's sensitive to S?

@rgcoe
Copy link
Member Author

rgcoe commented May 5, 2020

OK, so after I check the value of the exitflag (doh), it seems like your expected failure case doesn't actually fail. Have I done something else wrong? Maybe it's sensitive to S?

Yes, a sensitivity to S (specifically to the frequencies at which we run NEMOH) is definitely likely. The frequencies that would be evaluated in a case has been updated by 1e2866a (#129) and probably a few other times since that specific problem case was noted in #4...

My thinking is we should complete this PR, but keep an eye on this (your addition of the error check on the exitflag with 65164f6 should allow us to do so).

@H0R5E
Copy link
Member

H0R5E commented May 5, 2020

Seems OK to me. Comment out the failing test so we can reuse later if needed?

Ryan Coe added 2 commits May 5, 2020 10:06
our test case parameters for this
do not currently create a non-
convergent solution, so we must
wait to find a set of parameters
that result in such an error
verifyError(testCase, test, "WecOptTool:getPSPhasePower:NoConvergence");
end

function testPSConverge(testCase)
Copy link
Member Author

@rgcoe rgcoe May 5, 2020

Choose a reason for hiding this comment

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

Should we be somehow performing an operation on testCase in testPSConverge (somehow explicitly pass to the testCase whether an error has occurred)?

Copy link
Member

Choose a reason for hiding this comment

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

Well, the test passes if the case doesn't error, so you don't strictly have to check something. Checking its numerical value (again what it is producing now) would probably be useful so that we will see if it changes for some reason. I was just lazy putting it in so I could show you what was happening.

Copy link
Member

Choose a reason for hiding this comment

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

You might consider replacing the case in RM3DeviceModelParametricCCTest with this one, so that we are not doing more NEMOH runs.

Copy link
Member Author

Choose a reason for hiding this comment

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

My little brain cannot quite digest this... I guess my main hang up is that I understand how to reuse data within a single file (e.g., RM3DeviceModelParametricCCTest.m), but not how to share between files.

I also feel like the +tests directory could use some clean up. I don't understand where I should look between all of the files and how tests should be broken up between files. Is there a suggested structure for this?

Copy link
Member

Choose a reason for hiding this comment

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

You can't share between files, I don't think, so anything you want to share should be in a single file.

The +unit package is getting a bit unruly for sure. We could add another package level below (it won't upset the testing system) for when there are several tests for one source file. Perhaps something else for our coding policy.

I'll fix these things, so you can see what I mean.

@H0R5E
Copy link
Member

H0R5E commented May 6, 2020

So, I pushed that update. It's interesting because on the second run of the PS controller with the stored NEMOH result the output is about 3% less (I added a relative tolerance of 5% to make it pass). It's likely this is a bug somewhere.

@H0R5E
Copy link
Member

H0R5E commented May 6, 2020

Further to that this number it's producing: 6.556208812689474e+09 is too big right? So, I'm not sure this is working, but we are not capturing why it's not working.

@rgcoe
Copy link
Member Author

rgcoe commented May 6, 2020

Further to that this number it's producing: 6.556208812689474e+09 is too big right? So, I'm not sure this is working, but we are not capturing why it's not working.

this sounds like an unconverged case... can you point me to a specific test case and/or MWE?

@H0R5E
Copy link
Member

H0R5E commented May 6, 2020

It's the same test as before, but it's now in toolbox/+WecOptLib/+tests/+unit/+models/RM3DeviceModelParametricPSTest.m.

It's being run in the getPower method and then the value of that run is being checked in test_runParametric. I just copied the value it was producing to expSol, but I didn't really think about it.

Note, this test is being run with the negative damping fix and it's not throwing an error in toolbox/+WecOptLib/+models/+RM3/getPSPhasePower.m (line 66) so there must be something we are not checking.

@H0R5E
Copy link
Member

H0R5E commented May 13, 2020

Having played around with the optimiser for a bit, I think it's a red herring. I think it's worth looking at the start of the chain and seeing if the NEMOH results are sound. @ryancoe do you have any advice on how best to check the quality of the NEMOH model?

@H0R5E H0R5E marked this pull request as ready for review September 14, 2020 13:52
@H0R5E
Copy link
Member

H0R5E commented Sep 14, 2020

OK, here is a list of the pertinent changes in this PR:

  • Added the WecOptTool.Hydrodynamics class to take standardised output of solvers and do stuff. Currently it fixes any negative damping found on the diagonal of the B matrix.
  • Added a warning to the RM3 model if the maximum absolute combined damping (i.e. B39) exceeds the maximum absolute damping of the individual bodies.
  • Added plotMesh which can plot the mesh structures generated by the AxiMesh class. An example of its usage has been added to the RM3/optimization.m example.

The only query I have from this (which is unrelated to the damping issue) is the relationship between wave amplitude and spectral density, which is used in a few places in the code. Is there some reference for this somewhere? I was wondering if it is missing a factor of pi?

@H0R5E H0R5E removed their request for review September 14, 2020 14:00
@H0R5E
Copy link
Member

H0R5E commented Sep 14, 2020

@ryancoe can you review this, please?

Comment on lines 84 to 87
% Calculate wave amplitude
waveAmp = SS.getAmpSpectrum(w, interpMethod);
% TODO: Why times root(2) here?
waveAmpSS = sqrt(2 * SS.dw * SS.S);
waveAmp = interp1(SS.w, waveAmpSS, w, interpMethod, 'extrap');

Copy link
Member

@H0R5E H0R5E Sep 14, 2020

Choose a reason for hiding this comment

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

@ryancoe, I couldn't figure out this relationship between wave amplitude and spectral density. What is the root(2) doing in there and shouldn't there be a factor of pi, also?

Copy link
Member Author

Choose a reason for hiding this comment

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

I believe it looks correct to me. You would need a 2pi if the frequency were in Hz.

image
https://ocw.mit.edu/courses/mechanical-engineering/2-019-design-of-ocean-systems-spring-2011/lecture-notes/MIT2_019S11_OWE.pdf

Copy link
Member

@H0R5E H0R5E Sep 14, 2020

Choose a reason for hiding this comment

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

Doesn't that indicate that A is sqrt(2 * S * dw) / 2? Unless amplitude is supposed to say wave height?

Copy link
Member

Choose a reason for hiding this comment

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

Hmm, I think that half is incorrect in those note you posted. If you look at https://www.usna.edu/NAOE/_files/documents/Courses/EN455/EN455_Chapter3.pdf page 40, then it has the formula used in the code, and the derivation there makes sense to me.

Copy link
Member Author

Choose a reason for hiding this comment

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

Yes, agreed. I was going to send you a photo from a textbook, but then decided to try to respect copyright rules and find something publicly available. Even MIT is wrong sometimes :)

@H0R5E
Copy link
Member

H0R5E commented Sep 14, 2020

Test results:

test_results.pdf

@H0R5E H0R5E mentioned this pull request Sep 14, 2020
1 task
@rgcoe
Copy link
Member Author

rgcoe commented Sep 15, 2020

Tests all passing on macOS: test_results.pdf

@rgcoe
Copy link
Member Author

rgcoe commented Sep 15, 2020

Tests still passing after my attempts to break it: test_results.pdf

@rgcoe rgcoe merged commit 00c6603 into SNL-WaterPower:master Sep 15, 2020
@H0R5E
Copy link
Member

H0R5E commented Sep 16, 2020

I think we need to put a little more effort into preparing the commit messages for the merged pull requests. I don't think you would have a clue what was added from the default commit message that was generated from this one.

@rgcoe
Copy link
Member Author

rgcoe commented Sep 16, 2020

Agreed, I only recently understood how this works with squash and merge

@H0R5E
Copy link
Member

H0R5E commented Sep 16, 2020

Added rule 10 to the DevOps Rule Book

@rgcoe rgcoe deleted the checkBemIrregFreqs branch September 25, 2020 16:36
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

bug Something isn't working enhancement New feature or request

Projects

None yet

Development

Successfully merging this pull request may close these issues.

noncoverged PS cases?

2 participants