I have an O(N^2) solution with DP.

Let dp*[j] be the answer, if our progression ends with the j-th and i-th element (j <= i; dp** == 1). Let X* be a hash-table (or just an array, if A* are small enough, or a map<> if you don’t want a hash table, for an additional O(log N) factor to the complexity), in which you store pairs (A*-A[j],j) for j < i; in case of a tie in A*-A[j], just the larger j. This table can be constructed easily from dp* in linear time.

What happens when you want to know dp*[j]? It’s just a progression ending with the j-th element, and you append the i-th element, so dp*[j] == 1+dp[k][j]. The index k is the largest index (<= j) for which A[j]-A[k] == A*-A[j], so it’s the second part of the entry (A*-A[j],k) in X[j]. If that entry doesn’t exist, then dp*[j] == 2. Thanks to the hash-table, we can find dp*[j] in O(1), so we have O(N^2) total complexity.

I find it hard to believe that there could be an algorithm in O(N log N) (note that O-notation means “or better”), so this’ll have to do…

Code (includes excluded :D)

```
int N, ans =1;
vector<int> A(N); // the input array
vector< vector<int> > dp(N, vector<int>(N,1));
vector< unordered_map<int,int> > X(N);
for(int i =0; i < N; i++) {
for(int j =0; j < i; j++) {
if(X[j].find(A*-A[j]) == X[j].end()) {
dp*[j] =2;
continue;}
int k =X[j][A*-A[j]];
dp*[j] =max(dp*[j],dp[j][k]+1);}
for(int j =0; j < i; j++) {
X*[A*-A[j]] =j;
ans =max(ans,dp*[j]);}
}
cout << ans << "
```

";