-
Notifications
You must be signed in to change notification settings - Fork 61
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.
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
Added a plotting utility for aeroelastic modes #315
base: develop
Are you sure you want to change the base?
Conversation
Added support for the Krylov as the eigendecomposition for the full model can be very slow. This required the linear source code to need updating to propagate the Krylov projectors to the postprocessor. Also cleaned up a lot of unused code for the Krylov, however, no numerics are changed.
Okay, so I realise a similar, but more limited tool, already exists under |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
LGTM! Thanks for cleaning up the code in the other files too
for i_eig, eig in enumerate(evals_c_ae): | ||
# we here round the eigenvalues to 5 decimal places to prevent a pair not matching to machine precision | ||
eig_pos = np.round(eig.real + 1j * np.abs(eig.imag), 5) | ||
if eig_pos not in eig_set: | ||
eig_set.add(eig_pos) | ||
first_eig_index.append(i_eig) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Style suggestion: It might be worth renaming i_eig
to eig_index
(or something to that effect). As this deals with complex numbers, i_eig
could be misconstrued as being the imaginary component.
I don't know if this affects the logic of this, but would a dictionary be a better data structure as sets aren't ordered? That way, we'd be able to guarantee the index corresponds to a given eig_pos
.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Good spot. I think I was lucky in the fact that because I'm doing nothing crazy with the set, they were staying ordered - have swapped it over to a dictionary. Thanks for the pointers!
node_disp_base = struct_evects @ evecs_ae_trunc[state_index['q'], :] | ||
|
||
# mode shapes of the structural nodes, bound gamma and wake gamma | ||
node_disp = np.einsum('ij,k->ijk', node_disp_base, phase).real # [num_dof, num_lambda, num_phase] |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Love the use of np.einsum
!
num_dof = node_disp.shape[0] | ||
num_node = int(num_dof // 6) | ||
|
||
i_pos = np.array(list(chain([6 * i, 6 * i + 1, 6 * i + 2] for i in range(num_node)))) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Is the chain()
equivalent to flattening the nested list?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Yes in this case it is. It is really best used to lazily evaluated generator expressions together, but it gives me the correct behaviour for lists so I'll use it
surf_number.append(i_wake) | ||
|
||
# create folder to save modes to | ||
folder = self.data.output_folder + '/aeroelastic_modes/' |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Does hard-coding the /
in the path mean that SHARPy doesn't work on non-unix systems? If so, it might be worth progressively moving towards using os.path.join()
method
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This is copied from the equivalent code in the modal solver, so have tried to keep things consistent. Can change this if anyone actually encounters issues from it, but if so the issue would need addressing elsewhere in the code too
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Yeah, I think it would cause it to break in other parts of the code. It was a thinking aloud thing about if it would be worth having it such that newer code uses the os.path.join()
and when we encounter the old style, to migrate it over. It boils down to the balance between consistent code style and robust coding that we want to achieve.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I have a few other branches of cleaned up code lying around, so when I merge those I can do a check and swap everything over. Although I would say it is quite low on the list of things to do at the moment
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Agreed, it is very low on the to do list
This includes a new
AeroelasticModal
post processor which plots the aeroelastic modes, found from the linearised aeroelastic system. Whilst looking at the aeroelastic eigenvalues is useful for determining stability, this post processor is useful for investigating what the modes look like, which can be hard to interpret from a vector of complex numbers!The output is similar to the
Modal
solver which plots the natural modes of the structure projected onto the aerodynamic grid as VTK files, however this here also includes the aerodynamic gamma terms for both the bound wing panels and the wake. As there is a phase component in the aeroelastic eigenvectors, I have included the ability to plot at evenly spaced points, which gives an animation of what the mode looks like. Please enjoy an example gif of the output below for the Pazy wing, however I have also validated it for the FLEXOP with more aerodynamic surfaces.To keep it compatible with the Krylov reduction (allowing for much faster analysis) I was required to make small changes to the linear system source code, only so that gain matrices are propagated to the post processor as well as information about the original linear system states. It can be very interesting to see what aerodynamic features the Krylov likes to keep, and what it discards. I increased the default Krylov order from 1 to 4, because 1 usually gives very poor results.
I also removed some unused code here - the functionality of this solver should however be as before.