CODEWITHSUNDEEP
openmdao
BeeperTeeper
BeeperTeeper

Reputation: 161

Distinguishing between incorrect, small partials and correct partials not matching due to finite difference step size

The code in question:

class SecondPanelCoordinateInPanelReferenceFrame(om.ExplicitComponent):

    def initialize(self):
        self.options.declare('n', types=int)

    def setup(self):
        n = self.options['n']

        self.add_input('panel_angles', shape=(n,))
        self.add_input('second_panel_coordinate_relative_positions', shape=(2, n))

        self.add_output('second_panel_coordinate_transformed_positions', shape=(2, n))

        self.declare_partials('second_panel_coordinate_transformed_positions', 'panel_angles',
                              rows=np.arange(n),
                              cols=np.arange(n))
        self.declare_partials('second_panel_coordinate_transformed_positions', 'second_panel_coordinate_relative_positions',
                              rows=np.tile(np.arange(n), 2),
                              cols=np.arange(2*n))

    def compute(self, inputs, outputs, discrete_inputs=None, discrete_outputs=None):
        self.second_panel_coordinate_relative_positions = inputs['second_panel_coordinate_relative_positions']
        self.panel_angles = inputs['panel_angles']

        self.x2r = self.second_panel_coordinate_relative_positions[0, :]  # x2r = x coordinate number 2, relative
        self.y2r = self.second_panel_coordinate_relative_positions[1, :]

        self.cos_angles = np.cos(self.panel_angles)
        self.sin_angles = np.sin(self.panel_angles)

        second_panel_coordinate_transformed_positions = np.zeros_like(self.second_panel_coordinate_relative_positions)
        second_panel_coordinate_transformed_positions[0, :] = self.x2r*self.cos_angles + self.y2r*self.sin_angles

        outputs['second_panel_coordinate_transformed_positions'] = second_panel_coordinate_transformed_positions

    def compute_partials(self, inputs, partials, discrete_inputs=None):

        partials['second_panel_coordinate_transformed_positions',
                 'panel_angles'] = self.y2r*self.cos_angles - self.x2r*self.sin_angles

        partials['second_panel_coordinate_transformed_positions',
                 'second_panel_coordinate_relative_positions'] = np.concatenate([self.cos_angles, self.sin_angles])

I'm pretty confident that I have correctly written the partials for second_panel_coordinate_transformed_positions w.r.t panel_angles. However, when running check_partials(), it flags up a pretty large difference in terms of orders of magnitude.

  InfluenceCoefficientMatrixGroup.SecondPanelCoordinateInPanelReferenceFrame: 'second_panel_coordinate_transformed_positions' wrt 'panel_angles'
    Forward Magnitude : 8.659561e-17
         Fd Magnitude : 7.071696e-07 (fd:forward)
    Absolute Error (Jfor  - Jfd) : 7.071696e-07

    Relative Error (Jfor  - Jfd) : 1.000000e+00 *

    Raw Forward Derivative (Jfor)
[[6.123234e-17 0.000000e+00]
 [0.000000e+00 6.123234e-17]
 [0.000000e+00 0.000000e+00]
 [0.000000e+00 0.000000e+00]]

    Raw FD Derivative (Jfd)
[[-5.0004445e-07  0.0000000e+00]
 [ 0.0000000e+00 -5.0004445e-07]
 [ 0.0000000e+00  0.0000000e+00]
 [ 0.0000000e+00  0.0000000e+00]]

Both the forward and finite difference partials are very close to 0, and I'm unsure whether my partials are correct or whether the error due to finite step is showing up more strongly than I expect here. How can I work out where the problem is?

Upvotes: 1

Views: 73

Answers (1)

Rob Falck
Rob Falck

Reputation: 2704

When the absolute value of a derivative is very small, its generally fine to ignore the relative error.

In this case, since your compute function is complex-safe, we can use complex-step to get a more accurate approximation of the partials.

import openmdao.api as om
import numpy as np

class SecondPanelCoordinateInPanelReferenceFrame(om.ExplicitComponent):

    def initialize(self):
        self.options.declare('n', types=int)

    def setup(self):
        n = self.options['n']

        self.add_input('panel_angles', shape=(n,))
        self.add_input('second_panel_coordinate_relative_positions', shape=(2, n))
        self.add_output('second_panel_coordinate_transformed_positions', shape=(2, n))

        self.declare_partials('second_panel_coordinate_transformed_positions', 'panel_angles',
                              rows=np.arange(n),
                              cols=np.arange(n))
        self.declare_partials('second_panel_coordinate_transformed_positions', 'second_panel_coordinate_relative_positions',
                              rows=np.tile(np.arange(n), 2),
                              cols=np.arange(2*n))

    def compute(self, inputs, outputs, discrete_inputs=None, discrete_outputs=None):
        self.second_panel_coordinate_relative_positions = inputs['second_panel_coordinate_relative_positions']
        self.panel_angles = inputs['panel_angles']

        self.x2r = self.second_panel_coordinate_relative_positions[0, :]  # x2r = x coordinate number 2, relative
        self.y2r = self.second_panel_coordinate_relative_positions[1, :]

        self.cos_angles = np.cos(self.panel_angles)
        self.sin_angles = np.sin(self.panel_angles)

        second_panel_coordinate_transformed_positions = np.zeros_like(self.second_panel_coordinate_relative_positions)
        second_panel_coordinate_transformed_positions[0, :] = self.x2r*self.cos_angles + self.y2r*self.sin_angles

        outputs['second_panel_coordinate_transformed_positions'] = second_panel_coordinate_transformed_positions

    def compute_partials(self, inputs, partials, discrete_inputs=None):
        self.second_panel_coordinate_relative_positions = inputs['second_panel_coordinate_relative_positions']
        self.panel_angles = inputs['panel_angles']

        self.x2r = self.second_panel_coordinate_relative_positions[0, :]  # x2r = x coordinate number 2, relative
        self.y2r = self.second_panel_coordinate_relative_positions[1, :]

        self.cos_angles = np.cos(self.panel_angles)
        self.sin_angles = np.sin(self.panel_angles)

        partials['second_panel_coordinate_transformed_positions',
                 'panel_angles'] = self.y2r*self.cos_angles - self.x2r*self.sin_angles

        partials['second_panel_coordinate_transformed_positions',
                 'second_panel_coordinate_relative_positions'] = np.concatenate([self.cos_angles, self.sin_angles])


if __name__ == '__main__':
    n = 5

    p = om.Problem(model=om.Group())

    ivc = p.model.add_subsystem('ivc', om.IndepVarComp(), promotes_outputs=['*'])

    ivc.add_output('panel_angles', shape=(n,))
    ivc.add_output('second_panel_coordinate_relative_positions', shape=(2, n))

    p.model.add_subsystem('panel_comp', SecondPanelCoordinateInPanelReferenceFrame(n=n), promotes_inputs=['*'])

    p.setup(force_alloc_complex=True)

    p.set_val('panel_angles', np.random.rand(n))
    p.set_val('second_panel_coordinate_relative_positions', np.random.rand(2, n))

    with np.printoptions(linewidth=1024, edgeitems=1000):
        p.check_partials(method='cs')

Which gives us:

------------------------------------------------------------------
Component: SecondPanelCoordinateInPanelReferenceFrame 'panel_comp'
------------------------------------------------------------------
  panel_comp: 'second_panel_coordinate_transformed_positions' wrt 'panel_angles'
    Forward Magnitude : 9.957115e-01
         Fd Magnitude : 9.957115e-01 (cs:None)
    Absolute Error (Jfor  - Jfd) : 3.695246e-16

    Relative Error (Jfor  - Jfd) : 3.711162e-16

    Raw Forward Derivative (Jfor)
[[ 0.18045286  0.          0.          0.          0.        ]
 [ 0.         -0.54661688  0.          0.          0.        ]
 [ 0.          0.         -0.07489634  0.          0.        ]
 [ 0.          0.          0.         -0.15068612  0.        ]
 [ 0.          0.          0.          0.          0.79484104]
 [ 0.          0.          0.          0.          0.        ]
 [ 0.          0.          0.          0.          0.        ]
 [ 0.          0.          0.          0.          0.        ]
 [ 0.          0.          0.          0.          0.        ]
 [ 0.          0.          0.          0.          0.        ]]

    Raw FD Derivative (Jfd)
[[ 0.18045286  0.          0.          0.          0.        ]
 [ 0.         -0.54661688  0.          0.          0.        ]
 [ 0.          0.         -0.07489634  0.          0.        ]
 [ 0.          0.          0.         -0.15068612  0.        ]
 [ 0.          0.          0.          0.          0.79484104]
 [ 0.          0.          0.          0.          0.        ]
 [ 0.          0.          0.          0.          0.        ]
 [ 0.          0.          0.          0.          0.        ]
 [ 0.          0.          0.          0.          0.        ]
 [ 0.          0.          0.          0.          0.        ]]

 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  panel_comp: 'second_panel_coordinate_transformed_positions' wrt 'second_panel_coordinate_relative_positions'
    Forward Magnitude : 2.236068e+00
         Fd Magnitude : 2.236068e+00 (cs:None)
    Absolute Error (Jfor  - Jfd) : 1.922963e-16

    Relative Error (Jfor  - Jfd) : 8.599751e-17

    Raw Forward Derivative (Jfor)
[[0.86465361 0.         0.         0.         0.         0.50236852 0.         0.         0.         0.        ]
 [0.         0.54817618 0.         0.         0.         0.         0.83636288 0.         0.         0.        ]
 [0.         0.         0.8643305  0.         0.         0.         0.         0.50292423 0.         0.        ]
 [0.         0.         0.         0.82015556 0.         0.         0.         0.         0.57214059 0.        ]
 [0.         0.         0.         0.         0.98318242 0.         0.         0.         0.         0.18262619]
 [0.         0.         0.         0.         0.         0.         0.         0.         0.         0.        ]
 [0.         0.         0.         0.         0.         0.         0.         0.         0.         0.        ]
 [0.         0.         0.         0.         0.         0.         0.         0.         0.         0.        ]
 [0.         0.         0.         0.         0.         0.         0.         0.         0.         0.        ]
 [0.         0.         0.         0.         0.         0.         0.         0.         0.         0.        ]]

    Raw FD Derivative (Jfd)
[[0.86465361 0.         0.         0.         0.         0.50236852 0.         0.         0.         0.        ]
 [0.         0.54817618 0.         0.         0.         0.         0.83636288 0.         0.         0.        ]
 [0.         0.         0.8643305  0.         0.         0.         0.         0.50292423 0.         0.        ]
 [0.         0.         0.         0.82015556 0.         0.         0.         0.         0.57214059 0.        ]
 [0.         0.         0.         0.         0.98318242 0.         0.         0.         0.         0.18262619]
 [0.         0.         0.         0.         0.         0.         0.         0.         0.         0.        ]
 [0.         0.         0.         0.         0.         0.         0.         0.         0.         0.        ]
 [0.         0.         0.         0.         0.         0.         0.         0.         0.         0.        ]
 [0.         0.         0.         0.         0.         0.         0.         0.         0.         0.        ]
 [0.         0.         0.         0.         0.         0.         0.         0.         0.         0.        ]]

 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -

One note about your compute...there are some uncommon circumstances in which compute_partials is called before compute is called with the same inputs. As a result its generally not recommended that you cache your calculations.

Upvotes: 1

Related Questions