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

Let dp[i][j] be the answer, if our progression ends with the j-th and i-th element (j <= i; dp[i][i] == 1). Let X[i] be a hash-table (or just an array, if A[i] 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[i]-A[j],j) for j < i; in case of a tie in A[i]-A[j], just the larger j. This table can be constructed easily from dp[i] in linear time.

What happens when you want to know dp[i][j]? It’s just a progression ending with the j-th element, and you append the i-th element, so dp[i][j] == 1+dp[k][j]. The index k is the largest index (<= j) for which A[j]-A[k] == A[i]-A[j], so it’s the second part of the entry (A[i]-A[j],k) in X[j]. If that entry doesn’t exist, then dp[i][j] == 2. Thanks to the hash-table, we can find dp[i][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[i]-A[j]) == X[j].end()) {
dp[i][j] =2;
continue;}
int k =X[j][A[i]-A[j]];
dp[i][j] =max(dp[i][j],dp[j][k]+1);}
for(int j =0; j < i; j++) {
X[i][A[i]-A[j]] =j;
ans =max(ans,dp[i][j]);}
}
cout << ans << "\n";
```