Updates to the Elastic mesh solver and AD-wrapped code - #760
Conversation
…on FEA solver in SU2.
…ored in geometry.
rsanfer
left a comment
There was a problem hiding this comment.
Hi @pcarruscag
thanks for the comments. Indeed, we will temporarily have 3 different mesh deformation routines, but as you say this introduces little duplication and I hope we can put some PR soon to remove a lot more lines that the ones added here.
The FSI deformation should be relatively easy to replace, as all functionalities are already here (only some work is needed on the adjoint, which goes in line with the discussion we had with @oleburghardt.
The original grid deformation requires, as you mention, to implement the symmetry boundary condition, and also to adapt the other routines to work with absolute (rather than incremental) deformations. It shouldn't be too hard, but it requires to go over a lot of test cases as the residuals will change (the resulting deformed mesh is now different). This of course will take time... But also, I think @talbring mentioned a lot of functionality is duplicated in the regression tests (and some other is not covered), so that should be part of a bigger discussion.
Cheers,
Ruben
| * \param[in] iVertex - Vertex identifier. | ||
| * \return Vector of velocities at time n. | ||
| */ | ||
| vector<passivedouble> GetFEA_Velocity_n(unsigned short iMarker, unsigned short iVertex); |
There was a problem hiding this comment.
Yep, I'll add some documentation soon, but I'd rather do it with a tutorial. I think I should be able to do it by the time we do the release.
| void Set_YoungModulus(unsigned short i_DV, su2double val_Young); | ||
|
|
||
| void SetElement_Properties(CElement *element_container, CConfig *config); | ||
| virtual void SetElement_Properties(CElement *element_container, CConfig *config); |
There was a problem hiding this comment.
Nice catch. I think I changed it at some point to overload the SetElement_Properties, but then I didn't need it anymore and I forgot to change it back. Solved.
|
|
||
| su2double WallDistance; /*!< \brief Store the wall distance in reference coordinates. */ | ||
|
|
||
| su2double *Mesh_Coord; /*!< \brief Store the reference coordinates of the mesh. */ |
There was a problem hiding this comment.
Nope, the Mesh_Coord variable stores the original coordinates of the mesh. The "Solution" vector is the displacement, so that the current coordinate of any mesh point is Mesh_Coord + Solution. At some point I had implemented it as "Displacement", but I changed it back for the exact same reason you mention here.
| delete [] SlidingStateNodes; | ||
| } | ||
|
|
||
|
|
There was a problem hiding this comment.
I just did, although that meant I needed to move the nMarker and nVertex declarations also there (no big deal). That's sorted as well then.
Proposed Changes
This is the initial PR for a long overdue refactorization of the mesh deformation routines. The mesh deformation acquires here a full "solver" status, and uses the FEA Elasticity classes as its basis. The pseudo-linear elastic solver now always refers to the initial, undeformed geometrical discretization, therefore ensuring that the problem remains traceable to the original mesh, and that when volumes are computed to deform the mesh, the initial quality of the elements is observed.
With this framework, the wall distances only need to be computed once for distance-based stiffness, and the stiffness matrix is now symmetric after applying BCs so conjugate gradient can be used to compute the mesh deformation. The adjoint sensitivities can now be propagated down to the surface mesh points without the need to use DOT, although further work is needed to propagate them to the FFD points or other DV of choice.
I tried to keep config options simple: the new mesh solver is enabled using
DEFORM_MESH = YESand the relevant markers which will move and introduce deformations into the domain are set using
MARKER_DEFORM_MESH= ( wall )Hopefully, previous mesh deformation structures can be deprecated and removed soon, but we will need to adapt the rigid/imposed boundary deformation of other kinds of grid movement. I kindly ask for your help in there, as I'm not very familiar with some of the legacy code for this purpose.
Additionally, in this PR I have done some refactoring so the AD-based adjoint routines are accessible from the python-wrapped version of the code. The SWIG wrapper now compiles the AD binaries and interfaces directly to sensitivity information, which should make it much easier to integrate SU2 within any python-based optimization framework.
Related Work
PR Checklist
Put an X by all that apply. You can fill this out after submitting the PR. If you have any questions, don't hesitate to ask! We want to help. These are a guide for you to know what the reviewers will be looking for in your contribution.