More volume projection fixes

[anilyil]

Purpose

The volume projection routine changes the inverse evaluation into an optimization problem and utilizes the second derivative information to find the u,v,w parameters that correspond to given x,y,z coordinates. However, this second order Newton iteration does not always work; it looks like the final Jacobian ends up producing update vectors that make no sense. I think this might be caused by the second order derivative information.

So in this PR, I removed the contributions of the second order derivatives from the Newton solver, at least temporarily. With this change, the solver works much more robustly for me and @sseraj.

There can be 2 ways I think the second order derivatives are messing up the update:

1. The resulting optimization problem is non-convex, and therefore exactly evaluating the second order derivatives without any positive definiteness check causes the optimization step to fail.
2. There is a bug in the second derivative code. I documented a similar problem in the curve derivatives in First derivatives of curves give inaccurate results under complex step. #52

Ideally, if the reason is the first one listed above, we either need to improve additional checks in the optimization algorithm, or convert the entire thing back to a regular Newton solver. Currently, it looks like it's solving an optimization problem but only using part of the derivative information. Effectively, this is the Newton's method with a weird nonlinear preconditioning step. It does work, but it is not as easy to understand as applying the Newton's method directly.

Expected time until merged

The fix is fairly quick. I created an issue to track this. I will keep looking into the problem, until then, merging this is fine for me.

I will dig a bit deeper into this and try to figure out the issue, and then we can test these terms again.

Type of change

• Bugfix (non-breaking change which fixes an issue)
• New feature (non-breaking change which adds functionality)
• Breaking change (non-backwards-compatible fix or feature)
• Code style update (formatting, renaming)
• Refactoring (no functional changes, no API changes)
• Documentation update
• Maintenance update
• Other (please describe)

Testing

Checklist

• I have run flake8 and black to make sure the Python code adheres to PEP-8 and is consistently formatted
• I have formatted the Fortran code with fprettify or C/C++ code with clang-format as applicable
• I have run unit and regression tests which pass locally with my changes
• I have added new tests that prove my fix is effective or that my feature works
• I have added necessary documentation

[anilyil]

I am tagging @sseraj and @marcomangano because they helped me review a similar PR previously.

[anilyil]

These ADflow tests failed with this change:

test_adjoint.py:TestAdjoint_1_euler_scalar_JST_tut_wing.test_adjoint
test_adjoint.py:TestAdjoint_3_rans_tut_wing.test_adjoint
test_adjoint.py:TestAdjoint_4_Rotating_wing.test_adjoint
test_functionals.py:TestFunctionals_6_rans_tut_wing.test_dot_products

It looked like they were all small tolerance issues, I will look into these now.

[codecov]

Codecov Report

Merging #51 (af8aa87) into main (f4fbb14) will not change coverage.
The diff coverage is n/a.

@@           Coverage Diff           @@
##             main      #51   +/-   ##
=======================================
  Coverage   52.53%   52.53%           
=======================================
  Files           5        5           
  Lines        1616     1616           
=======================================
  Hits          849      849           
  Misses        767      767           

📣 We’re building smart automated test selection to slash your CI/CD build times. Learn more

[anilyil]

Here's the test output on these tests:

AssertionError: Exception raised on rank 0: 
Not equal to tolerance rtol=1e-10, atol=1e-10
Failed value for: Eval Functions Sens:: mdo_tutorial_cd: mdo_tutorial_cl: mdo_tutorial_cmz: mdo_tutorial_drag: P_mdo_tutorial: T_mdo_tutorial: alpha_mdo_tutorial: beta_mdo_tutorial: mach_mdo_tutorial: shape
Mismatched elements: 1 / 72 (1.39%)
Max absolute difference: 7.91987986e-09
Max relative difference: 1.44174571e-10
 x: array([[ 5.773086e+02, -2.913565e+02, -1.843664e+04, -1.400464e+04,
         1.320355e+03,  2.717153e+03, -8.866237e+03, -6.977819e+03,
         4.419584e+02, -1.564830e+02,  2.531239e+03,  7.378316e+03,...
 y: array([[ 5.773086e+02, -2.913565e+02, -1.843664e+04, -1.400464e+04,
         1.320355e+03,  2.717153e+03, -8.866237e+03, -6.977819e+03,
         4.419584e+02, -1.564830e+02,  2.531239e+03,  7.378316e+03,...



The following tests failed:
test_adjoint.py:TestAdjoint_1_euler_scalar_JST_tut_wing.test_adjoint

AssertionError: Exception raised on rank 0: 
Not equal to tolerance rtol=1e-10, atol=1e-10
Failed value for: Eval Functions Sens:: mdo_tutorial_cavitation: mdo_tutorial_cd: mdo_tutorial_cl: mdo_tutorial_cmz: mdo_tutorial_drag: mdo_tutorial_fx: mdo_tutorial_lift: P_mdo_tutorial: T_mdo_tutorial: alpha_mdo_tutorial: beta_mdo_tutorial: mach_mdo_tutorial: shape
Mismatched elements: 1 / 72 (1.39%)
Max absolute difference: 2.59868102e-07
Max relative difference: 1.29197405e-10
 x: array([[   -284.190245,   -5123.884097, -169289.381538, -137183.891827,
           4674.276904,    9501.87888 ,  -88389.956082,  -73500.97798 ,
          14305.727385,   26034.378961,   37310.054559,   72237.967364,...
 y: array([[   -284.190245,   -5123.884097, -169289.381538, -137183.891827,
           4674.276904,    9501.87888 ,  -88389.956083,  -73500.97798 ,
          14305.727385,   26034.378961,   37310.054559,   72237.967364,...



The following tests failed:
test_adjoint.py:TestAdjoint_3_rans_tut_wing.test_adjoint

AssertionError: Exception raised on rank 0: 
Not equal to tolerance rtol=1e-10, atol=1e-10
Failed value for: Eval Functions Sens:: mdo_tutorial_fy: mdo_tutorial_my: T_mdo_tutorial: alpha_mdo_tutorial: beta_mdo_tutorial: mach_mdo_tutorial: shape
Mismatched elements: 1 / 72 (1.39%)
Max absolute difference: 1.64891389e-08
Max relative difference: 3.30757316e-10
 x: array([[ -648.344403,   733.552076, -1291.819375,   802.485278,
          173.808769,  1152.27247 , -1766.586713,  -914.863547,
        -1991.160695, -4917.817028, -5058.395574, -3378.472428,...
 y: array([[ -648.344403,   733.552076, -1291.819375,   802.485278,
          173.808769,  1152.27247 , -1766.586713,  -914.863547,
        -1991.160695, -4917.817028, -5058.395574, -3378.472428,...



The following tests failed:
test_adjoint.py:TestAdjoint_4_Rotating_wing.test_adjoint

AssertionError: Exception raised on rank 0: 
Not equal to tolerance rtol=1e-10, atol=1e-10
Failed value for: Dot product test for xV -> F
Mismatched elements: 1 / 1 (100%)
Max absolute difference: 0.00019228
Max relative difference: 1.06597651e-10
 x: array(1803756.473321)
 y: array(1803756.473129)



The following tests failed:
test_functionals.py:TestFunctionals_6_rans_tut_wing.test_dot_products

[sseraj]

The changes are good with me. Once the ADflow tests are retrained, we can merge this and the PR on ADflow at the same time.

[marcomangano]

This also looks good to me, but we need one of the maintainers to approve for this to be merged anyway