[ceres-solver] AutoDiff

본고의 목적은ceres-solver AutoDiff의 실현을 해석하고 matlab 기호 연산과 유사한 방법이라는 것을 설명하는 것이다.
ceres-solver는 ceres::CostFunction를 계산 오차와 아크비의 구조로 사용한다.ceres: CostFunction은 순수한 허류로 사용자 코드가 이 클래스를 계승하고 순수한 허법을 실현하는 방법bool Evaluate(double const* const* parameters, double* residuals, double** jacobians);을 통해 최적화된 매개 변수 블록(parameters)을 사용하여 오차(residuals)와 아크비(jacobians)를 계산하는 방법을 제공합니다.아이디어를 신속하게 검증해야 하는 사용자에게는 아크비 계산이 번거롭다.
ceres는 두 가지 자동 계산 아크비 방법인 AutoDiff와NumericDiff를 제공했는데 사용자는 각각 ceres::AutoDiffCostFunction,ceres::NumericDiffCostFunction를 계승하여 이 두 가지 방법을 사용할 수 있다.이 두 가지 방법을 선택한 후 사용자 코드는ceres가parameters를 사용하여residuals를 계산하는 방법만 알려주고 jacobians가 어떻게 계산하는지는ceres가 스스로 방법을 찾는다.
ceres의 AutoDiff는 Dual Number를 사용하여 야크비를 계산합니다.Dual Number란 실수 하나를 자신의 (대량이라고 부르는 것을 편리하게 하기 위해) 와 소량 (e) 의 합으로 쓰고\(e^2 = 0\) (1 단계 도수를 계산할 때 이렇게 정의함) 을 정의하는 것이다.ceres가 실현한 Dual Number의 구조는ceres::Jet,Jet 구조의 대량은T a;,소량은Eigen::Matrix v;(이곳에서 소량은 하나의 Eigen:::Vector 표현은 다원 함수가 여러 변수에 대한 구도에 대한 고려이고 뒤에 설명할 것이다).
파일jet.h에는 Dual Number가 도수를 어떻게 계산하는지 설명하는 주석이 있습니다.지금 주석의 예를 따서 아래와 같다.

// To handle derivatives of functions taking multiple arguments, different
// infinitesimals are used, one for each variable to take the derivative of. For
// example, consider a scalar function of two scalar parameters x and y:
//
//   f(x, y) = x^2 + x * y
//
// Following the technique above, to compute the derivatives df/dx and df/dy for
// f(1, 3) involves doing two evaluations of f, the first time replacing x with
// x + e, the second time replacing y with y + e.
//
// For df/dx:
//
//   f(1 + e, y) = (1 + e)^2 + (1 + e) * 3
//               = 1 + 2 * e + 3 + 3 * e
//               = 4 + 5 * e
//
//               --> df/dx = 5
//
// For df/dy:
//
//   f(1, 3 + e) = 1^2 + 1 * (3 + e)
//               = 1 + 3 + e
//               = 4 + e
//
//               --> df/dy = 1
//

구함수f(x, y) = x^2 + x * y는 (1,3)에서 현재 나는 미적분의 수학 방식으로 도수를 계산한다.
위의 설명은 Dual Number 계산 함수\(f(x, y)=x^2+xy\)를 사용하여\((1, 3)\)에서\(x, y\)의 도수를 맞추는 과정을 설명합니다.x에 대한 편도는 x를 Dual Number 1 + e로 표시하고 y를 실수 3으로 표시하며 함수식 계산에 대입하여 최종적으로 얻은 e의 일차 항계수는 함수가 (1,3)에서 x에 대한 편도이다.사실 이것은 L'Hospital 법칙의 컴퓨터 구현이다.현재 미적분에서 배운 방법으로 도수를 계산한다.
\[\begin{align} {\partial f(x, y)\over\partial x} &=\text{lim}_{\Delta x\rightarrow 0} {f(x +\Delta x, y)\over\Delta x} otag\\&=\text{lim}_{\Delta x\rightarrow 0} {(x +\Delta x)^2 + (x +\Delta x) y\over\Delta x} otag\\&=\text{lim}_{\Delta x\rightarrow 0} {(x +\Delta x)^2 + (x +\Delta x) y\over\Delta x} otag\\&=\text{lim}_{\Delta x\rightarrow 0} {x^2 +\Delta x^2 + 2 x\Delta x + xy +\Delta x y\over\Delta x} otag\\&=\text{lim}_{\Delta x\rightarrow 0} {x^2 + xy + (2 x + y)\Delta x +\Delta x^2\over\Delta x} otag\\&=^{(*)}\text{lim}_{\Delta x\rightarrow 0} {(2 x + y) + 2\Delta x\over 1} otag\\&= 2x + y\\{\partial f(x, y)\over\partial x} |_{x=1, y = 3} &= 2 * 1 + 3 = 5\end{align}\]
주석: (*) 분자 분모가 각각\(\Deltax\)에 대해 도수를 구하는 L'Hospital 법칙을 한 번 사용합니다.
분석: 도수를 구할 때 분모는 일반적으로 1회항, 즉\(\Deltax\)이다.L'Hospital 법칙을 사용하면 분모에 대한 구도는\(\Deltax\) 0회 항목의 구도를 사라집니다.마침\(\Delta x\) 1번의 항목 구도 뒤에 상수가 있습니다.\(\Delta x\) 1회 이상 항목은 구도 후에\(\Delta x\)가 남아 극한 구도 후에 사라집니다.따라서 도수는 1회항에 대응하는 계수이고 프로그램 실현에서 e에 대응하는 계수이다.(하지만 여기서는\(\Delta x\) 0번 이하의 항목을 고려하지 않았습니다. 지금은 결정할 수 없습니다.)같은 이치로 2차 도수를 구하면\(\Deltax^2\)의 계수를 얻는다.
다음 주석에 이어 Eigen::Vector 의 해석이 나왔다.

// To take the gradient of f with the implementation of dual numbers ("jets") in
// this file, it is necessary to create a single jet type which has components
// for the derivative in x and y, and passing them to a templated version of f:
//
//   template
//   T f(const T &x, const T &y) {
//     return x * x + x * y;
//   }
//
//   // The "2" means there should be 2 dual number components.
//   // It computes the partial derivative at x=10, y=20.
//   Jet x(10, 0);  // Pick the 0th dual number for x.
//   Jet y(20, 1);  // Pick the 1st dual number for y.
//   Jet z = f(x, y);
//
//   LOG(INFO) << "df/dx = " << z.v[0]
//             << "df/dy = " << z.v[1];
//

모든 변수에 대한 도수를 직접 구하려면 Dual Number의 e 개수가 늘어나고 두 개의 변수가 있으면 e 2개가 필요해서 사용ceres::Jet해야 한다.실험 검증, AutoDiff를 사용할 때 사용자 코드로 이루어진 템플릿 함수operator()에서 일부러 템플릿을 특례화하는 오류, typename이 Jet으로 특례화되었는지 검사할 수 있으며, N은 최적화될 매개 변수의 개수(주의,parameters의 개수,parameter blocks의 개수가 아님).
Jet이 Dual Number로서 구도를 실현하는 것은 구체적으로 구도를 실현하는 일반적인 법칙과 일부 기본 함수의 도수 공식이다.이러한 관련 내용은 WikiPedia Differentiation rules를 참조하십시오.현재 일반 법칙과 기본 함수의 도수에 대해 각각 파일jet.h에서 찾은 예를 들 수 있다.
일반 법칙C++에서 기본적으로 연산하는 operator는 +, -, *,/네 가지뿐입니다. 이 네 가지 연산에 대응하는 Jet operator만 있으면 됩니다.파이톤에는 멱 연산자**가 있는데 아마도 파이톤이 실현되려면 이것을 고려해야 할 것 같다.
곱셈법은 수학에서 다음과 같이 표현할 수 있다.
\[\begin{align} (f(x)g(x))^{\prime} = f(x)g(x)^{\prime} + f(x)^{\prime}g(x)\end{align}\]
Jet에 해당하는 operator*는 다음과 같습니다.

template 
inline Jet operator*(const Jet& f, const Jet& g) {
  return Jet(f.a * g.a, f.a * g.v + f.v * g.a);
}

기본 함수의 도수는 정현 함수를 예로 들 수 있다.
정현 함수의 도수는 수학에서 다음과 같이 표현된다.
\[\begin{align} (\sin(x))^{\prime} =\cos(x)\end{align}\]
Jet에 해당하는 함수sin는 다음과 같습니다.

template 
inline Jet sin(const Jet& f) {
  return Jet(sin(f.a), cos(f.a) * f.v);
}

상기 두 가지 예는 코드에서 형성된 Jet의 소량은 수학 공식의 도수에 대응한다.
또한ceres에서 AutoDiff, 템플릿 함수operator()를 사용하여residuals를 계산하는 과정에서 사용된 기초 함수는ceres에서 얻어야 한다. 즉, 직접 사용할 수 없는std::sin 함수는 사용해야 한다ceres::sin, 이상sin 함수에서 사용된sin, cos 함수는 std::sin, std::cos이다.즉, stl의 템플릿은 ceres::Jet을 인스턴스화할 수 없습니다.
변수 음수 차멱의 처리는 코드operator/(T s, const Jet& g)를 참고할 수 있다. 즉,'Scalar 나누기 Jet'이다.
요약하면 ceres-solver는 ceres:::Jet을 사용하여 AutoDiff를 실현했다.구체적인 실현은ceres:::Jet의 풍부한 Operator와 정의된 일련의 기본 함수(의 도수)를 통해 이루어진다.