Let \(N \coloneqq |S|\), and \(\sigma\) be the number of characters that may appear. In our problem, \(\sigma = 26\).

The string obtained by concatenating \(S_i, S_j\), and \(S_k\) in this order is a palindrome if and only if \(S_i = S_k\).

After rephrasing like this, there are multiple possible approaches. Here, we introduce two of them.

(1) Focus on \(j\)

Let us find the count for each fixed \(j\). The number of \((i, k)\) is the number of pairs with \(S_i = S_k\), \(1 \leq i < j\), and \(j < k \leq N\).

When we fix \(c \coloneqq S_i\) (note that we do not fix the index, but the character, which is one of only \(26\) kinds), what we want is the number of \(1 \leq i < j\) with \(S_i = c\), and \(j < k \leq N\) with \(S_k = c\).

In the same manner as cumulative sums, we may count the occurrences of \(c\) in the \(1, 2, \ldots , i\)-th characters for each \(i\) and \(c\) firsthand, so that the sought counts can be obtained in a constant time; this way, it can be processed in a total of \(O(\sigma N)\) time.

One can also use the fact that the sought counts, the numbers of \(1 \leq i < j\) with \(S_i = c\) and \(j < k \leq N\) with \(S_k = c\) is almost the same for \(j\) and \((j + 1)\). By applying appropriate differential updates, the complexity can be optimized from \(O(\sigma N)\) to \(O(N)\) time.

(2) Focus on \(i\) and \(k\)

If \(S_i = S_k\), any \(j\) with \(i < j < k\) is applicable; there are \(k - i - 1\) choices. Thus, the answer is the sum of \(k - i - 1\) over all pairs \((i, k)\) with \(S_i = S_k\).

If we define \(X_c\) as the sequence of the indices \(i\) with \(S_i = c\), the answer is the sum of \(\displaystyle\sum_{s = 1}^{|X_c|}\displaystyle\sum_{t = s + 1}^{|X_c|} (X_{c, t} - X_{c, s} - 1)\) over \(c =\) A, B, \(\ldots\), Z.

This can be evaluated in \(O(|X_c|)\) time by considering the contribution of each \(X_{c, s}\), or calculating the cumulative sums \(\displaystyle\sum_{s = 1}^{s'}X_{c,s}\). Thus, by applying this method for each \(c =\) A, B, \(\ldots\), Z, the answer can be found in a total of \(O(N)\) time.

Sample code (1)

#include <bits/stdc++.h>
using namespace std;

#define rep(i, n) for (ll i = 0; i < (n); i++)
using ll = long long;

int main() {
	string s;
	cin >> s;
	int n = s.size();
	vector<vector<int>> sum(26, vector<int>(n + 1));
	rep(i, n) {
		rep(j, 26) {
			sum[j][i + 1] = sum[j][i];
		}
		sum[s[i] - 'A'][i + 1]++;
	}
	ll ans = 0;
	for (int i = 1; i < n - 1; i++) {
		rep(j, 26) {
			ll l = sum[j][i];
			ll r = sum[j][n] - sum[j][i + 1];
			ans += l * r;
		}
	}
	cout << ans << '\n';
}

Sample code (2)

#include <bits/stdc++.h>
using namespace std;

#define rep(i, n) for (ll i = 0; i < (n); i++)
using ll = long long;

int main() {
	string s;
	cin >> s;
	vector<ll> cnt(26), sum(26);
	ll ans = 0;
	int n = s.size();
	rep(i, n) {
		int v = s[i] - 'A';
		ans += (i - 1) * cnt[v] - sum[v];
		cnt[v]++;
		sum[v] += i;
	}
	cout << ans << '\n';
}