Autodiff with Julia

예~~전에 짜둔 코드인데, 친구가 안 예쁘다고 해서 조금 코드를 손봤다.

손 보고 나니 전체적으로 깔끔해서 매우 마음에 들었다.

Autodiff 에 대한 전반적인 내용은 이전 글 (https://calofmijuck.tistory.com/28) 을 참고하면 될 것이다.

바로 구현으로 들어가자!

객체 정의

우선 DiffObject struct 를 정의해야 한다. Julia 문법이지만 대충 무슨 말인지 이해는 되니 너무 신경쓰지는 않아도 괜찮을 듯 하다.

struct DiffObject <: Number
    f::Number
    df::Number

    DiffObject(f::Number, df::Number) = new(f, df)
    DiffObject(f::Number) = new(f, one(f))
end

<: 는 subtype 이라는 의미이다. DiffObject 도 수(number)처럼 취급할 것이다. f, df 는 각각 함숫값과 미분계수이며, 생성자는 2개 만들어줬다. 하나는 함숫값과 미분계수를 같이 주는 경우이며, 나머지 하나는 함숫값만 주는 경우이다.

참고로 저기서 one 은 f 의 타입에 맞는 곱셈에 대한 항등원을 돌려주는 함수이다. 그러므로 f::Int64 = 1 이면 one(f) = 1 이고, f::Float64 = 1.0 이면 one(f) = 1.0 이며, 심지어 복소수의 경우에도 된다.

연산자 오버로딩

Julia 에서는 무려 기본 연산자를 오버로딩 할 수 있다! 우선 오버로딩을 위해 import 부터 한다.

import Base: +, -, *, /, ^

구현에 대한 설명은 이전 글에 있으므로 미분을 어떻게 했는지에 대한 설명은 전부 생략하겠다.

+(x::DiffObject, y::DiffObject) = DiffObject(x.f + y.f, x.df + y.df)

-(x::DiffObject, y::DiffObject) = DiffObject(x.f - y.f, x.df - y.df)

-(x::DiffObject) = DiffObject(-x.f, -x.df)

*(x::DiffObject, y::DiffObject) = DiffObject(x.f * y.f, x.df * y.f + x.f * y.df)

/(x::DiffObject, y::DiffObject) = DiffObject(x.f / y.f, (x.df * y.f - x.f * y.df) / y.f^2)

^(x::DiffObject, y::DiffObject) = DiffObject(x.f^y.f, y.f * x.f^(y.f - 1) * x.df + x.f^y.f * log(x.f) * y.df)

Julia 문법에 대해 조금 부연 설명을 하면, Julia 에서 연산자 또한 함수나 다름없다. 평상시에는 3 + 5 라고 표기하지만 (infix) 사실 계산할 때는 3, 5 를 인자로 해서 + 함수를 호출한다. 따라서 사실은 +(3, 5) 인 것이다. (prefix)

그러므로 +(x::DiffObject, y::DiffObject) 를 정의하는 것은 두 DiffObject 를 더하는 연산을 새롭게 정의해주는 것이다.

또한 +(x::DiffObject) (unary plus) 는 구현하지 않아도 되는데, unary plus 는 그냥 입력을 그대로 돌려준다.

추가로 Julia 문법은 수학의 표기법과 꽤나 유사하다. 함수 정의를 한 줄에 할 수 있으며, ^ 는 bitwise XOR 이 아닌 우리가 흔히 텍스트로 수학 얘기할 때 사용하는 거듭제곱이다.

Data Promotion

DiffObject 와 다른 타입 간의 연산도 구현해야 하지만, Julia 에서는 타입을 바꿔주면 그만이다.

import Base: convert, promote_rule

convert(::Type{DiffObject}, x::Real) = DiffObject(x, zero(x))

promote_rule(::Type{DiffObject}, ::Type{<:Number}) = DiffObject

convert(::Type{DiffObject}, x::Real) 는 x::Real 을 DiffObject 로 변환 (type casting / conversion) 하려고 할 때 어떻게 할지 정의하는 함수이다. zero 는 위의 one 과 유사하게 타입에 맞는 덧셈에 대한 항등원을 돌려준다.

promote_rule(::Type{DiffObject}, ::Type{<:Number}) 은 DiffObject 와 Number 의 subtype 을 연산하려고 할 때 각각을 어떤 타입으로 promote 해야할지 정의하는 함수이다. 예를 들어,

promote_rule(::Type{Float64}, ::Type{Float32}) = Float64

이다. (Float64 와 Float32 를 연산하면 당연히 각각을 Float64 로 변환하여 계산할 것이다.)

현재 상황에서는 DiffObject 로 promote 해주면 될 것이다.

그리고 promotion 함수는 대칭성을 갖고 있기 때문에 (symmetric) promote_rule(::Type{A}, ::Type{B}) 와 promote_rule(::Type{B}, ::Type{A}) 중 하나만 구현해도 충분하다. Python 에서 __radd__ 등을 구현해야 했는데 그런거 안해도 된다.

여기까지만 하면 아쉬우므로...

삼각, 지수, 로그함수

Python 으로는 하지 못했던 삼각, 지수, 로그함수도 미분할 수 있다! Julia 에서는 내장 함수도 오버로딩 할 수 있다!

import Base: sin, cos, tan, exp, log

sin(x::DiffObject) = DiffObject(sin(x.f), x.df * cos(x.f))

cos(x::DiffObject) = DiffObject(cos(x.f), -x.df * sin(x.f))

tan(x::DiffObject) = DiffObject(tan(x.f), x.df * sec(x.f)^2)

exp(x::DiffObject) = DiffObject(exp(x.f), x.df * exp(x.f))

log(x::DiffObject) = DiffObject(log(x.f), x.df / x.f)

미분계수가 왜 저렇게 되는지는 딱히 설명할 필요 없을 것 같다.

동작 확인

마지막으로 잘 동작함을 확인해보았다.

겁나 복잡한 함수로 확인!

WolframAlpha 에 넣어봤다.

너무 똑똑해서 소수로 안구하고 exact form 을 구해버림

정확한 형태로 안 준다. Plain text 를 복사해서 붙여넣었다.

(2^(4 - π/2 + π^2/4) e)/(π (1 + π^2/4)) - (2^(3 - π/2 + π^2/4) e π log(π/2))/(1 + π^2/4)^2 + (2^(3 - π/2 + π^2/4) e (π - 1) log(2) log(π/2))/(1 + π^2/4)

아무튼 맞음

소숫점 아래 14자리까지 일치한다.

Python 너도 이런거 할 줄 알아? ㅋㅋ

---

이상으로 autodiff 에 대한 이야기는 끝! 다 짜고 나니 오히려 python 보다 코드 자체는 훨씬 간결해 보이는 느낌이 든다. Julia 가 확실히 이런 쪽으로는 강력한 기능을 가지고 있는 것 같다.