:orphan: :py:mod:`tlm_adjoint.markers` ============================= .. py:module:: tlm_adjoint.markers Module Contents --------------- .. py:class:: ControlsMarker(M) Represents .. math:: m = m_\text{input}, where :math:`m` is the control and :math:`m_\text{input}` the input value for the control. The forward residual is defined .. math:: \mathcal{F} \left( m \right) = m - m_\text{input}. :arg M: A variable or a :class:`Sequence` of variables defining the control :math:`m`. May be static. .. !! processed by numpydoc !! .. py:method:: adjoint_jacobian_solve(adj_X, nl_deps, B) Compute an adjoint solution. :arg adj_X: Either `None`, or a variable or :class:`Sequence` of variables defining the initial guess for an iterative solve. May be modified or returned. :arg nl_deps: A :class:`Sequence` of variables defining values for non-linear dependencies. Should not be modified. :arg B: The right-hand-side. A variable or :class:`Sequence` of variables storing the value of the right-hand-side. May be modified or returned. :returns: A variable or :class:`Sequence` of variables storing the value of the adjoint solution. May return `None` to indicate a value of zero. .. !! processed by numpydoc !! .. py:class:: FunctionalMarker(J) Represents .. math:: J_\text{output} = J, where :math:`J` is the functional and :math:`J_\text{output}` is the output value for the functional. The forward residual is defined .. math:: \mathcal{F} \left( J_\text{output}, J \right) = J_\text{output} - J. :arg J: A variable defining the functional :math:`J`. .. !! processed by numpydoc !! .. py:method:: adjoint_derivative_action(nl_deps, dep_index, adj_x) Return the action of the adjoint of a derivative of the forward residual on the adjoint solution. This is the *negative* of an adjoint right-hand-side term. :arg nl_deps: A :class:`Sequence` of variables defining values for non-linear dependencies. Should not be modified. :arg dep_index: An :class:`int`. The derivative is defined by differentiation of the forward residual with respect to `self.dependencies()[dep_index]`. :arg adj_X: The adjoint solution. A variable or a :class:`Sequence` of variables. Should not be modified. :returns: The action of the adjoint of a derivative on the adjoint solution. Will be passed to :func:`.subtract_adjoint_derivative_action`, and valid types depend upon the adjoint variable type. Typically this will be a variable, or a two element :class:`tuple` `(alpha, F)`, where `alpha` is a :class:`numbers.Complex` and `F` a variable, with the value defined by the product of `alpha` and `F`. .. !! processed by numpydoc !! .. py:method:: adjoint_jacobian_solve(adj_x, nl_deps, b) Compute an adjoint solution. :arg adj_X: Either `None`, or a variable or :class:`Sequence` of variables defining the initial guess for an iterative solve. May be modified or returned. :arg nl_deps: A :class:`Sequence` of variables defining values for non-linear dependencies. Should not be modified. :arg B: The right-hand-side. A variable or :class:`Sequence` of variables storing the value of the right-hand-side. May be modified or returned. :returns: A variable or :class:`Sequence` of variables storing the value of the adjoint solution. May return `None` to indicate a value of zero. .. !! processed by numpydoc !! .. py:class:: AdjointActionMarker(J, X, adj_X) Represents .. math:: J_\text{output} = \lambda_x^* x, with forward residual .. math:: \mathcal{F} \left( J_\text{output}, x \right) = J_\text{output} - \lambda_x^* x. Note that :math:`\lambda_x` is *not* treated as a dependency. Can be used to initialize an adjoint calculation, and compute adjoint Jacobian actions, via the construction .. code-block:: python start_manager() X = forward(M) with paused_manager(): adj_X = ... J = Float(name="J") AdjointRHSMarker(J, X, adj_X).solve() stop_manager() # Compute the action of the adjoint of the Jacobian on the direction # defined by adj_X dJ = compute_gradient(J, M) :arg J: A variable defining the functional :math:`J`. :arg X: A variable or :class:`Sequence` of variables defining :math:`x`. :arg adj_X: A variable or :class:`Sequence` of variables defining :math:`\lambda_x`. .. !! processed by numpydoc !! .. py:method:: forward_solve(x, deps=None) Compute the forward solution. Can assume that the currently active :class:`.EquationManager` is paused. :arg X: A variable or a :class:`Sequence` of variables storing the solution. May define an initial guess, and should be set by this method. :arg deps: A :class:`tuple` of variables, defining values for dependencies. Only the elements corresponding to `X` may be modified. `self.dependencies()` should be used if not supplied. .. !! processed by numpydoc !! .. py:method:: adjoint_derivative_action(nl_deps, dep_index, adj_x) Return the action of the adjoint of a derivative of the forward residual on the adjoint solution. This is the *negative* of an adjoint right-hand-side term. :arg nl_deps: A :class:`Sequence` of variables defining values for non-linear dependencies. Should not be modified. :arg dep_index: An :class:`int`. The derivative is defined by differentiation of the forward residual with respect to `self.dependencies()[dep_index]`. :arg adj_X: The adjoint solution. A variable or a :class:`Sequence` of variables. Should not be modified. :returns: The action of the adjoint of a derivative on the adjoint solution. Will be passed to :func:`.subtract_adjoint_derivative_action`, and valid types depend upon the adjoint variable type. Typically this will be a variable, or a two element :class:`tuple` `(alpha, F)`, where `alpha` is a :class:`numbers.Complex` and `F` a variable, with the value defined by the product of `alpha` and `F`. .. !! processed by numpydoc !! .. py:method:: adjoint_jacobian_solve(adj_x, nl_deps, b) Compute an adjoint solution. :arg adj_X: Either `None`, or a variable or :class:`Sequence` of variables defining the initial guess for an iterative solve. May be modified or returned. :arg nl_deps: A :class:`Sequence` of variables defining values for non-linear dependencies. Should not be modified. :arg B: The right-hand-side. A variable or :class:`Sequence` of variables storing the value of the right-hand-side. May be modified or returned. :returns: A variable or :class:`Sequence` of variables storing the value of the adjoint solution. May return `None` to indicate a value of zero. .. !! processed by numpydoc !! .. py:method:: tangent_linear(tlm_map) Derive an :class:`.Equation` corresponding to a tangent-linear operation. :arg tlm_map: A :class:`.TangentLinearMap` storing values for tangent-linear variables. :returns: An :class:`.Equation`, corresponding to the tangent-linear operation. .. !! processed by numpydoc !!