update state generators, optimizers, and losses

This commit is contained in:
Brian Delhaisse
2019-04-06 04:56:22 +02:00
parent c7ea4af161
commit bb2c22dd29
7 changed files with 480 additions and 234 deletions
+2
View File
@@ -97,6 +97,8 @@ def concatenate(data, axis=0, out=None):
# missing torch.dstack
# missing torch.vstack
# missing torch.hstack
# np.ndim --> torch.dim()
# np.reshape --> torch.view()
# Tests
+1
View File
@@ -1,3 +1,4 @@
# import losses
from .loss import *
from .losses import *
+4 -234
View File
@@ -1,13 +1,13 @@
#!/usr/bin/env python
"""Defines the common loss functions that are used by the learning algorithm / optimizer.
"""Defines the abstract loss class that is used by the learning algorithm / optimizer.
Losses are evaluated on model parameters, data batches, and / or storages.
TODO: loss vs cost vs reward
costs and rewards are defined possibly for each time steps.
Note that rewards / costs and losses are different concepts. Losses are minimized with respect to parameters to
optimize them, while rewards / costs depends on the state(s) and action(s) and are returned by the environment.
"""
from abc import ABCMeta, abstractmethod
from abc import ABCMeta
import operator
import copy
import collections
@@ -283,229 +283,6 @@ class Loss(object):
return loss
class FixedLoss(Loss):
r"""Fixed Loss
"""
def __init__(self, value):
super(FixedLoss, self).__init__()
self.value = value
def compute(self, batch):
return self.value * batch
class L2Loss(Loss):
r"""L2 Loss
Compute the L2 loss given by: :math:`1/2 * (y_{target} - y_{predict})^2`
"""
def __init__(self, target, approximator):
super(L2Loss, self).__init__()
self.target = target
self.approximator = approximator
def compute(self, batch):
# based on approximator check what we need
return 0.5 * (self.target(arg) - self.approximator(arg)).pow(2).mean()
class ValueLoss(Loss):
r"""L2 loss for values
"""
def __init__(self):
super(ValueLoss, self).__init__()
def compute(self, batch):
returns = batch['returns']
values = batch.current['values']
return 0.5 * (returns - values).pow(2).mean()
class PGLoss(Loss):
r"""Policy Gradient Loss
Compute the policy gradient loss which is maximized and given by:
.. math:: L^{PG} = \mathbb{E}[ \log \pi_{\theta}(a_t | s_t) \psi_t ]
where :math:`\psi_t` is the associated return estimator, which can be for instance, the total reward estimator
:math:`\psi_t = R(\tau)` (where :math:`\tau` represents the whole trajectory), the state action value estimator
:math:`\psi_t = Q(s_t, a_t)`, or the advantage estimator :math:`\psi_t = A_t = Q(s_t, a_t) - V(s_t)`. Other
estimators are also possible.
The gradient with respect to the parameters :math:`\theta` is then given by:
.. math:: g = \mathbb{E}[ \nabla_\theta \log \pi_{\theta}(a_t | s_t)
References:
[1] "Proximal Policy Optimization Algorithms", Schulman et al., 2017
[2] "High-Dimensional Continuous Control using Generalized Advantage Estimation", Schulman et al., 2016
"""
def __init__(self):
super(PGLoss, self).__init__()
def compute(self, batch):
log_curr_pi = batch.current['action_distributions']
log_curr_pi = log_curr_pi.log_probs(batch.current['actions'])
estimator = batch['estimator']
loss = torch.exp(log_curr_pi) * estimator
return -loss.mean()
def latex(self):
return "\\mathbb{E}[ r_t(\\theta) A_t ]"
class CPILoss(Loss):
r"""CPI Loss
Conservative Policy Iteration objective which is maximized and defined in [1]:
.. math:: L^{CPI}(\theta) = \mathbb{E}[ r_t(\theta) A_t ]
where the expectation is taken over a finite batch of samples, :math:`A_t` is an estimator of the advantage fct at
time step :math:`t`, :math:`r_t(\theta)` is the probability ratio given by
:math:`r_t(\theta) = \frac{ \pi_{\theta}(a_t|s_t) }{ \pi_{\theta_{old}}(a_t|s_t) }`.
References:
[1] "Approximately optimal approximate reinforcement learning", Kakade et al., 2002
[2] "Proximal Policy Optimization Algorithms", Schulman et al., 2017
"""
def __init__(self):
"""
Initialize the CPI Loss.
"""
super(CPILoss, self).__init__()
def compute(self, batch): # policy_distribution, old_policy_distribution, estimator):
# ratio = policy_distribution / old_policy_distribution
log_curr_pi = batch.current['action_distributions']
log_curr_pi = log_curr_pi.log_probs(batch.current['actions'])
log_prev_pi = batch['action_distributions']
log_prev_pi = log_prev_pi.log_probs(batch['actions'])
ratio = torch.exp(log_curr_pi - log_prev_pi)
estimator = batch['estimator']
loss = ratio * estimator
return -loss.mean()
def latex(self):
return "\\mathbb{E}[ r_t(\\theta) A_t ]"
class CLIPLoss(Loss):
r"""CLIP Loss
Loss defined in [1] which is maximized and given by:
.. math:: L^{CLIP}(\theta) = \mathbb{E}[ \min(r_t(\theta) A_t, clip(r_t(\theta), 1-\epsilon, 1+\epsilon) A_t) ]
where the expectation is taken over a finite batch of samples, :math:`A_t` is an estimator of the advantage fct at
time step :math:`t`, :math:`r_t(\theta)` is the probability ratio given by
:math:`r_t(\theta) = \frac{ \pi_{\theta}(a_t|s_t) }{ \pi_{\theta_{old}}(a_t|s_t) }`.
References:
[1] "Proximal Policy Optimization Algorithms", Schulman et al., 2017
"""
def __init__(self, clip=0.2):
"""
Initialize the loss.
Args:
epsilon (float): clip parameter
"""
super(CLIPLoss, self).__init__()
self.eps = clip
def compute(self, batch): # , policy_distribution, old_policy_distribution, estimator):
log_curr_pi = batch.current['action_distributions']
log_curr_pi = log_curr_pi.log_probs(batch.current['actions'])
log_prev_pi = batch['action_distributions']
log_prev_pi = log_prev_pi.log_probs(batch['actions'])
ratio = torch.exp(log_curr_pi - log_prev_pi)
estimator = batch['estimator']
loss = torch.min(ratio * estimator, torch.clamp(ratio, 1.0-self.eps, 1.0+self.eps) * estimator)
return -loss.mean()
def latex(self):
return "\\mathbb{E}[ \\min(r_t(\\theta) A_t, clip(r_t(\\theta), 1-\\epsilon, 1+\\epsilon) A_t) ]"
class KLPenaltyLoss(Loss):
r"""KL Penalty Loss
KL Penalty to minimize:
.. math:: L^{KL}(\theta) = \mathbb{E}[ KL( \pi_{\theta_{old}}(a_t | s_t) || \pi_{\theta}(a_t | s_t) ) ]
where :math:`KL(.||.)` is the KL-divergence between two probability distributions.
"""
def __init__(self): # p, q):
"""
Initialize the KL Penalty loss.
Args:
p (torch.distributions.Distribution): 1st distribution
q (torch.distributions.Distribution): 2nd distribution
"""
super(KLPenaltyLoss, self).__init__()
# self.p = p
# self.q = q
def compute(self, batch):
"""
Compute :math:`KL(p||q)`.
"""
curr_pi = batch.current['action_distributions']
prev_pi = batch['action_distributions']
return torch.distributions.kl.kl_divergence(prev_pi, curr_pi)
def latex(self):
return "\\mathbb{E}[ KL( \\pi_{\\theta_{old}}(a_t | s_t) || \\pi_{\\theta}(a_t | s_t) ) ]"
# class ForwardKLPenaltyLoss(KLPenaltyLoss):
# r"""Forward KL Penalty Loss"""
# pass
#
#
# class ReverseKLPenaltyLoss(KLPenaltyLoss):
# r"""Reverse KL Penalty Loss"""
# pass
class EntropyLoss(Loss):
r"""Entropy Loss
Entropy loss, which is used to ensure sufficient exploration when maximized [1,2,3]:
.. math:: L^{Entropy}(\theta) = H[ \pi_{\theta} ]
where :math:`H[.]` is the Shannon entropy of the given probability distribution.
References:
[1] "Simple Statistical Gradient-following Algorithms for Connectionist Reinforcement Learning", Williams, 1992
[2] "Asynchronous Methods for Deep Reinforcement Learning", Mnih et al., 2016
[3] "Proximal Policy Optimization Algorithms", Schulman et al., 2017
"""
def __init__(self): # approximator):
super(EntropyLoss, self).__init__()
def compute(self, batch):
distribution = batch.current['action_distributions']
entropy = distribution.entropy().mean()
return entropy
# TODO: the functions defined here should be the same for `reward.py` and `loss.py`. Thus, we need to check numpy
# TODO: or torch inside the functions.
def ceil(x):
@@ -763,10 +540,3 @@ def trunc(x):
return y
else:
return torch.trunc(x)
# Tests
if __name__ == '__main__':
# compute the losses
loss = - FixedLoss(3) # + FixedLoss(2)
print(loss(2))
+314
View File
@@ -0,0 +1,314 @@
#!/usr/bin/env python
"""Defines the common loss functions that are used by the learning algorithm / optimizer.
Losses are evaluated on model parameters, data batches / storages, or transitions tuples.
"""
import torch
from pyrobolearn.losses.loss import Loss
__author__ = "Brian Delhaisse"
__copyright__ = "Copyright 2018, PyRoboLearn"
__credits__ = ["Brian Delhaisse"]
__license__ = "MIT"
__version__ = "1.0.0"
__maintainer__ = "Brian Delhaisse"
__email__ = "briandelhaisse@gmail.com"
__status__ = "Development"
class FixedLoss(Loss):
r"""Fixed Loss
"""
def __init__(self, value):
super(FixedLoss, self).__init__()
self.value = value
def compute(self, batch):
return self.value * batch
class L2Loss(Loss):
r"""L2 Loss
Compute the L2 loss given by: :math:`1/2 * (y_{target} - y_{predict})^2`
"""
def __init__(self, target, approximator):
super(L2Loss, self).__init__()
self.target = target
self.approximator = approximator
def compute(self, batch):
# based on approximator check what we need
return 0.5 * (self.target(batch) - self.approximator(batch)).pow(2).mean()
class ValueLoss(Loss):
r"""L2 loss for values
"""
def __init__(self):
super(ValueLoss, self).__init__()
def compute(self, batch):
returns = batch['returns']
values = batch.current['values']
return 0.5 * (returns - values).pow(2).mean()
class HuberLoss(Loss):
r"""Huber Loss
"In statistics, the Huber loss is a loss function used in robust regression, that is less sensitive to outliers
in data than the squared error loss." [1]
This loss is given by [1]:
.. math:: {\mathcal L}(\delta) = \left\{ \begin{array}{lc} 1/2 a^2 & for |a| \leq \delta, \\
\delta (|a| - 1/2 \delta), & \mbox{otherwise} \left \end{array}
"This function is quadratic for small values of :math:`a`, and linear for large values, with equal values and
slopes of the different sections at the two points where :math:`|a| = \delta`" [1].
In [2], this loss is used for DQN, where :math:`\delta=1`, and :math:`a` is the temporal difference error, that is,
:math:`a = Q(s,a) - (r + \gamma \max_a Q(s',a))` where :math:`(r + \gamma \max_a Q(s',a))` is the target function.
References:
[1] Huber Loss (on Wikipedia): https://en.wikipedia.org/wiki/Huber_loss
[2] "Reinforcement Learning (DQN) Tutorial":
https://pytorch.org/tutorials/intermediate/reinforcement_q_learning.html
"""
def __init__(self, loss, delta=1.):
super(HuberLoss, self).__init__()
self.loss = loss
self.delta = delta
def compute(self, batch):
a = self.loss(batch)
if abs(a) <= self.delta:
return 0.5 * torch.pow(a, 2)
return self.delta * (torch.abs(a) - 0.5 * self.delta)
class PGLoss(Loss):
r"""Policy Gradient Loss
Compute the policy gradient loss which is maximized and given by:
.. math:: L^{PG} = \mathbb{E}[ \log \pi_{\theta}(a_t | s_t) \psi_t ]
where :math:`\psi_t` is the associated return estimator, which can be for instance, the total reward estimator
:math:`\psi_t = R(\tau)` (where :math:`\tau` represents the whole trajectory), the state action value estimator
:math:`\psi_t = Q(s_t, a_t)`, or the advantage estimator :math:`\psi_t = A_t = Q(s_t, a_t) - V(s_t)`. Other
estimators are also possible.
The gradient with respect to the parameters :math:`\theta` is then given by:
.. math:: g = \mathbb{E}[ \nabla_\theta \log \pi_{\theta}(a_t | s_t)
References:
[1] "Proximal Policy Optimization Algorithms", Schulman et al., 2017
[2] "High-Dimensional Continuous Control using Generalized Advantage Estimation", Schulman et al., 2016
"""
def __init__(self):
super(PGLoss, self).__init__()
def compute(self, batch):
log_curr_pi = batch.current['action_distributions']
log_curr_pi = log_curr_pi.log_probs(batch.current['actions'])
estimator = batch['estimator']
loss = torch.exp(log_curr_pi) * estimator
return -loss.mean()
def latex(self):
return "\\mathbb{E}[ r_t(\\theta) A_t ]"
class CPILoss(Loss):
r"""CPI Loss
Conservative Policy Iteration objective which is maximized and defined in [1]:
.. math:: L^{CPI}(\theta) = \mathbb{E}[ r_t(\theta) A_t ]
where the expectation is taken over a finite batch of samples, :math:`A_t` is an estimator of the advantage fct at
time step :math:`t`, :math:`r_t(\theta)` is the probability ratio given by
:math:`r_t(\theta) = \frac{ \pi_{\theta}(a_t|s_t) }{ \pi_{\theta_{old}}(a_t|s_t) }`.
References:
[1] "Approximately optimal approximate reinforcement learning", Kakade et al., 2002
[2] "Proximal Policy Optimization Algorithms", Schulman et al., 2017
"""
def __init__(self):
"""
Initialize the CPI Loss.
"""
super(CPILoss, self).__init__()
def compute(self, batch): # policy_distribution, old_policy_distribution, estimator):
# ratio = policy_distribution / old_policy_distribution
log_curr_pi = batch.current['action_distributions']
log_curr_pi = log_curr_pi.log_probs(batch.current['actions'])
log_prev_pi = batch['action_distributions']
log_prev_pi = log_prev_pi.log_probs(batch['actions'])
ratio = torch.exp(log_curr_pi - log_prev_pi)
estimator = batch['estimator']
loss = ratio * estimator
return -loss.mean()
def latex(self):
return "\\mathbb{E}[ r_t(\\theta) A_t ]"
class CLIPLoss(Loss):
r"""CLIP Loss
Loss defined in [1] which is maximized and given by:
.. math:: L^{CLIP}(\theta) = \mathbb{E}[ \min(r_t(\theta) A_t, clip(r_t(\theta), 1-\epsilon, 1+\epsilon) A_t) ]
where the expectation is taken over a finite batch of samples, :math:`A_t` is an estimator of the advantage fct at
time step :math:`t`, :math:`r_t(\theta)` is the probability ratio given by
:math:`r_t(\theta) = \frac{ \pi_{\theta}(a_t|s_t) }{ \pi_{\theta_{old}}(a_t|s_t) }`.
References:
[1] "Proximal Policy Optimization Algorithms", Schulman et al., 2017
"""
def __init__(self, clip=0.2):
"""
Initialize the loss.
Args:
epsilon (float): clip parameter
"""
super(CLIPLoss, self).__init__()
self.eps = clip
def compute(self, batch): # , policy_distribution, old_policy_distribution, estimator):
log_curr_pi = batch.current['action_distributions']
log_curr_pi = log_curr_pi.log_probs(batch.current['actions'])
log_prev_pi = batch['action_distributions']
log_prev_pi = log_prev_pi.log_probs(batch['actions'])
ratio = torch.exp(log_curr_pi - log_prev_pi)
estimator = batch['estimator']
loss = torch.min(ratio * estimator, torch.clamp(ratio, 1.0-self.eps, 1.0+self.eps) * estimator)
return -loss.mean()
def latex(self):
return "\\mathbb{E}[ \\min(r_t(\\theta) A_t, clip(r_t(\\theta), 1-\\epsilon, 1+\\epsilon) A_t) ]"
class KLPenaltyLoss(Loss):
r"""KL Penalty Loss
KL Penalty to minimize:
.. math:: L^{KL}(\theta) = \mathbb{E}[ KL( \pi_{\theta_{old}}(a_t | s_t) || \pi_{\theta}(a_t | s_t) ) ]
where :math:`KL(.||.)` is the KL-divergence between two probability distributions.
"""
def __init__(self): # p, q):
"""
Initialize the KL Penalty loss.
Args:
p (torch.distributions.Distribution): 1st distribution
q (torch.distributions.Distribution): 2nd distribution
"""
super(KLPenaltyLoss, self).__init__()
# self.p = p
# self.q = q
def compute(self, batch):
"""
Compute :math:`KL(p||q)`.
"""
curr_pi = batch.current['action_distributions']
prev_pi = batch['action_distributions']
return torch.distributions.kl.kl_divergence(prev_pi, curr_pi)
def latex(self):
return "\\mathbb{E}[ KL( \\pi_{\\theta_{old}}(a_t | s_t) || \\pi_{\\theta}(a_t | s_t) ) ]"
# class ForwardKLPenaltyLoss(KLPenaltyLoss):
# r"""Forward KL Penalty Loss"""
# pass
#
#
# class ReverseKLPenaltyLoss(KLPenaltyLoss):
# r"""Reverse KL Penalty Loss"""
# pass
class EntropyLoss(Loss):
r"""Entropy Loss
Entropy loss, which is used to ensure sufficient exploration when maximized [1,2,3]:
.. math:: L^{Entropy}(\theta) = H[ \pi_{\theta} ]
where :math:`H[.]` is the Shannon entropy of the given probability distribution.
References:
[1] "Simple Statistical Gradient-following Algorithms for Connectionist Reinforcement Learning", Williams, 1992
[2] "Asynchronous Methods for Deep Reinforcement Learning", Mnih et al., 2016
[3] "Proximal Policy Optimization Algorithms", Schulman et al., 2017
"""
def __init__(self): # approximator):
super(EntropyLoss, self).__init__()
def compute(self, batch):
distribution = batch.current['action_distributions']
entropy = distribution.entropy().mean()
return entropy
class MSBELoss(Loss):
r"""Mean-squared Bellman error
The mean-squared Bellman error (MSBE) computes the Bellman error, also known as the one-step temporal difference
(TD).
This is given by, in the case of the off-policy Q-learning TD algorithm:
.. math:: (r + \gamma (1-d) max_{a'} Q_{\phi}(s',a')) - Q_{\phi}(s,a),
or in the case of the on-policy Sarsa TD algorithm:
.. math:: (r + \gamma (1-d) Q_{\phi}(s',a')) - Q_{\phi}(s,a)
or in the case of the TD(0):
.. math:: (r + \gamma (1-d) V_{\phi}(s')) - V_{\phi}(s)
where (r + \gamma (1-d) f(s,a)) is called (one-step return) the target. The target could also be, instead of the
one step return, the n-step return or the lambda-return value. [2]
These losses roughly tell us how closely the value functions satisfy the Bellman equation. Thus, by trying to
minimize them, we try to enforce the Bellman equations.
References:
[1] https://spinningup.openai.com/en/latest/algorithms/ddpg.html
[2] "Reinforcement Learning: An Introduction", Sutton and Barto, 2018
"""
def __init__(self):
super(MSBELoss, self).__init__()
def compute(self, batch):
pass
# Tests
if __name__ == '__main__':
# compute the losses
loss = - FixedLoss(3) # + FixedLoss(2)
print(loss(2))
+75
View File
@@ -0,0 +1,75 @@
class TrustRegion(object):
r"""Trust Region
"Trust region is a term used in mathematical optimization to denote the subset of the region of the objective
function that is approximated using a model function (often a quadratic). If an adequate model of the objective
function is found within the trust region, then the region is expanded; conversely, if the approximation is poor,
then the region is contracted. Trust-region methods are also known as restricted-step methods.
The fit is evaluated by comparing the ratio of expected improvement from the model approximation with the actual
improvement observed in the objective function. Simple thresholding of the ratio is used as the criterion for
expansion and contraction; a model function is "trusted" only in the region where it provides a reasonable
approximation.
Trust-region methods are in some sense dual to line-search methods: trust-region methods first choose a step size
(the size of the trust region) and then a step direction, while line-search methods first choose a step direction
and then a step size." [1]
References:
[1] https://en.wikipedia.org/wiki/Trust_region
"""
pass
class LineSearch(object):
r"""Line Search
"In optimization, the line search strategy is one of two basic iterative approaches to find a local minimum
:math:`\mathbf{x}^*` of an objective function :math:`f:\mathbb{R}^{n} \to \mathbb{R}`. The other approach is trust
region.
The line search approach first finds a descent direction along which the objective function :math:`f` will be
reduced and then computes a step size that determines how far :math:`\mathbf{x}` should move along that direction.
The descent direction can be computed by various methods, such as gradient descent, Newton's method and
Quasi-Newton method. The step size can be determined either exactly or inexactly.
Here is an example gradient method that uses a line search in step 4.
1. Set iteration counter k = 0, and make an initial guess :math:`\mathbf{x}_{0}` for the minimum
2. Repeat:
3. Compute a descent direction :math:`\mathbf{p}_k`
4. Choose :math:`\alpha_k` to 'loosely' minimize :math:`h(\alpha)=f(\mathbf{x}_k + \alpha \mathbf{p}_k)` over
:math:`\alpha \in \mathbb{R}_{+}`
5. Update :math:`\mathbf{x}_{k+1} = \mathbf{x}_k + \alpha_k \mathbf{p}_k`, and :math:`k = k + 1`
6. Until :math:`|| \nabla f( \mathbf{x}_k ) || < tolerance
At the line search step (4) the algorithm might either exactly minimize :math:`h`, by solving
:math:`h'(\alpha _{k})=0`, or loosely, by asking for a sufficient decrease in :math:`h`. One example of the former
is conjugate gradient method. The latter is called inexact line search and may be performed in a number of ways,
such as a backtracking line search or using the Wolfe conditions.
Like other optimization methods, line search may be combined with simulated annealing to allow it to jump over
some local minima." [1]
References:
[1] https://en.wikipedia.org/wiki/Line_search
"""
pass
class BacktrackingLineSearch(LineSearch):
r"""Backtracking Line Search
"In (unconstrained) minimization, a backtracking line search, a search scheme based on the Armijo-Goldstein
condition, is a line search method to determine the maximum amount to move along a given search direction.
It involves starting with a relatively large estimate of the step size for movement along the search direction,
and iteratively shrinking the step size (i.e., "backtracking") until a decrease of the objective function is
observed that adequately corresponds to the decrease that is expected, based on the local gradient of the
objective function." [1]
References:
[1] https://en.wikipedia.org/wiki/Backtracking_line_search
"""
pass
+67
View File
@@ -6,9 +6,11 @@ instantiating them, here we can pass the parameters at a later stage.
References:
[1] https://pytorch.org/docs/stable/optim.html
[2] https://github.com/sbarratt/torch_cg/tree/master
"""
# Pytorch optimizers
import torch
import torch.nn as nn
import torch.optim as optim
@@ -200,3 +202,68 @@ class SGD(Optimizer):
if self.max_grad_norm is not None:
nn.utils.clip_grad_norm_(params, self.max_grad_norm)
self.optimizer.step()
class CG(Optimizer):
r"""Conjugate Gradient
"In mathematics, the conjugate gradient method is an algorithm for the numerical solution of particular systems
of linear equations, namely those whose matrix is symmetric and positive-definite. The conjugate gradient method
is often implemented as an iterative algorithm, applicable to sparse systems that are too large to be handled by
a direct implementation or other direct methods such as the Cholesky decomposition. Large sparse systems often
arise when numerically solving partial differential equations or optimization problems.
Suppose we want to solve the system of linear equations :math:`\mathbf{A} \mathbf{x} = \mathbf{b}`, for the vector
:math:`\mathbf{x}`, where the known n x n matrix :math:`\mathbf{A}` is symmetric (i.e.,
:math:`\mathbf{A}^\top = \mathbf{A}`), positive-definite (i.e. :math:`\mathbf{x}^\top \mathbf{Ax} > 0` for all
non-zero vectors :math:`x \in \mathbb{R}^n`), and real, and :math:`\mathbf{b}` is known as well. We denote the
unique solution of this system by :math:`\mathbf{x}^*`." [1]
Notably, this can be used to solve :math:`\mathbf{Hx} = \mathbf{g}` where :math:`\mathbf{H}` is the Hessian
matrix, and :math:`g` is the gradient.
References:
[1] Conjugate Gradient (Wikipedia): https://en.wikipedia.org/wiki/Conjugate_gradient_method
[2] Hessian matrix (Wikipedia): https://en.wikipedia.org/wiki/Hessian_matrix#Use_in_optimization
[3] Hessian-Vector products: https://justindomke.wordpress.com/2009/01/17/hessian-vector-products/
[4] Torch CG: https://github.com/sbarratt/torch_cg
"""
def __init__(self, threshold=1.e-8, max_iters=10, *args, **kwargs):
super(CG, self).__init__(*args, **kwargs)
self._threshold = threshold
self._max_iters = int(max_iters)
def optimize(self, A, b, x=None):
"""
Return the solution :math:`x` to the system of linear equations :math:`Ax=b`.
Notably, this can be used
Args:
A (torch.Tensor): square real symmetric and positive-definite matrix.
b (torch.Tensor): known vector.
x (torch.Tensor, None): initial solution. If None, it will be the zero vector.
Returns:
torch.Tensor: return the solution x to Ax=b.
"""
if x is None:
x = torch.zeros_like(b)
r = b - A.matmul(x)
p = r
rr_old = r.t().matmul(r)
for k in range(self._max_iters):
Ap = A.matmul(p)
alpha = rr_old / p.t().matmul(Ap)
x = x + alpha * p
r = r - alpha * Ap
rr_new = r.t().matmul(r)
if torch.sqrt(rr_new) < self._threshold:
break
p = r + rr_new / rr_old * p
rr_old = rr_new
return x
@@ -197,6 +197,23 @@ class QueueStateGenerator(StateGenerator):
"""Return the size of the queue."""
return self.qsize()
def __iter__(self):
"""Return the iterator object itself."""
self.cnt = 0
return self
def __next__(self): # only valid in Python 3
"""Return the next item in the sequence."""
if self.cnt < self.qsize():
self.cnt += 1
return self.queue[self.cnt-1]
else:
raise StopIteration
def next(self): # for Python 2
"""Return the next item in the sequence."""
return self.__next__()
class FIFOQueueStateGenerator(QueueStateGenerator):
r"""FIFO Queue Initial State Generator