Project Euler #101 (Optimum polynomial)

문제는 https://projecteuler.net/problem=101

처음에는 C언어 연습도 할 겸 malloc등을 이용해서 동적할당으로 배열 계산을 하려고 했다.

그러나.. u(10)이후 값은 unsigned long long 이나 unsigned __int64로 처리할 수 없었다.

determinant나 adjugate와 같은 함수들을 기껏 다 구현해놓고 데이터 크기때문에 사용할 수 없었다..ㅠ

그래서 결국엔 Matlab으로 했다..

---

우선 문제를 분석하면 아주 정통적인 행렬 문제란걸 알 수 있다.

예를들어 u(1), u(2), u(3)까지 커버되는 OP(3, n)의 식을 알고싶다면

이 만족되는 a, b, c를 구하면 된다. 구하고자하는 값이 3개이고, 방정식도 3개가 나오므로 노가다로 구할 수도 있지만

행렬을 공부해본 사람이면 역행렬을 곱하여 답을 쉽게 얻을 수 있다.

이 되므로 

이 사실을 이용하기위해..!

inverse matrix를 구하면 된다.

우선 u함수부터 먼저 만들어주고..

function y = u( x )
     y = 1 - x + x^2 - x^3 + x^4 - x^5 + x^6 - x^7 + x^8 - x^9 + x^10;
end

(a, b, c)앞에 오는 n X n 짜리의 커다란 덩어리 행렬을 만들어주는 함수를 만든다.

function y = cm( n )
    y = zeros(n, n);
    for j = 1:n;
        for i = 1:n;
            y(j, i) = j^(n-i);
        end
    end
end

그리고 (u(1), u(2), u(3), .., u(n))을 짜주는 함수도 만든다.

function y = my( n )
    y = zeros(n, 1);
    for i = 1:n;
        y(i) = u(i);
    end
end

이제 fit(first incorrect term)를 구한다. 구하는 방법은 지금의 경우는 n=3 이니까

이렇게 ((n+1)^(n-1), (n+1)^(n-2), ... , (n+1)^0)을 앞에 곱해주면 된다.

이걸 생성하는 함수는

function y = cv( n )
    y = zeros(1, n);
    for i = 1:n;
        y(1, i)=(n+1)^(n-i);
    end
end

이다. 결과적으로 fit는

function y = fit( m )
    y = cv(m)*(inv(cm(m))*my(m));
end

이다.

이제 남은건 fit(0) + ... + fit(10) 을 더하는 것이므로

r = 0;
for i = 1:10;
    r = r + fit(i);
end

로 얻어내진 답이 나오는데 숫자가 커서 바로 보이지 않는다.

이럴때

sprintf('%16.f',r)

로 자리수를 읽으면 된다.