This is really simple. Check out the idea below.

for j in s:
    c = ord(j) - ord('a')
    for i in range(26):
        c4[c][i][i][c] += c3[c][i][i]
        c3[i][c][c] += c2[i][c]
        c2[i][c] += c1[i]
    c1[c] += 1

Note that I omitted the modulo operations. All you need is some matrices. The result is the sum of the elements in c4.

True, in Python this solution times out, but it's nonsense. Python has its limitations and we all know that. Use PyPy when needed. Also, the C++ version of the code above gets in with no problem.

Remember: use the best tool for the problem, don't force yourself to only use Python or only use C++ or some other language.

In case someone is interested, here was my approach(longest test case completed in 1.66sec out of 4seconds for pypy3).

I start with this algorithm:
http://stackoverflow.com/questions/6877249/find-the-number-of-occurrences-of-a-subsequence-in-a-string

However, instead of having an inner loop where you iterate over each char in the sequence, i unrolled it and tracked the current row(4 units long) for every possible inner letter pair(aa, bb, cc, ect.). So, I'd calculate all possible a??a on the first pass, all b??b on the second pass, ect.. Then, I summed the answers for every pass for my answer.

M8... if you are getting error for test cases 30-37, even after using unsigned long long everywhere, start taking modulo one trillion and 7 for your array numbers. It really helps with the overflow and all that.

Here is my solution (passes everything) using that popular DP stackoverflow link. All this optimization making me crazy.

CPP:

#include <iostream>
using namespace std;
typedef unsigned long long ull;
#define FOR(i,n) for(ull (i)=0;(i)<(n);((i)++))

int main() {
    ull iArr[26][26][4]={}, iLen, iCum;
    char iString[1000001]; cin >> iString;
    iLen=iCum=0;
    while(iString[++iLen]){}
    FOR(m,iLen){
        char i = iString[m]-'a';
        FOR(j,26){
            iArr[i][j][3]+=iArr[i][j][2]%1000000007;
            iArr[j][i][2]+=iArr[j][i][1]%1000000007;
            iArr[j][i][1]+=iArr[j][i][0]%1000000007;
            iArr[i][j][0]++;
        }
    }

    FOR(i,26){
        FOR(j,26){
            iCum+=iArr[i][j][3];
        }
    }
    cout << iCum%1000000007;
    return 0;
}

Plain Python3 solution without pypy3.

Passes all the tests without timeout.

No manual unrolling the loop required.

def shortPalindrome(s):
    mod = 10**9 + 7
    C1 = [0] * 26
    C2 = [0] * 26 * 26
    C3 = [0] * 26 * 26
    count = 0
    r26 = list(range(26))
    for c in s:
        k = ord(c) - 97
        p = 26 * k - 1
        q = k - 26
        for i in r26:
            q += 26
            p += 1
            count += C3[q]
            C3[p] += C2[p]
            C2[p] += C1[i]
        C1[k] += 1

    return count%mod
   

I agree that this is not a medium problem. After much messing around, I finally managed to get in running on Python with no timeouts. Even with an O(n) time and O(1) space algorithm, Python can be limiting just based on lookup times and aritmetic operation time. If you are close but can't quite make it without timeouts, try PyPy (I was able to use the code I wrote for Python 3 without alteration). Also, in a similar vein list calls seem to be much more efficient than dictionary calls in PyPy (unlike Python where I found them to be similar) which is important for an algorithm where elements are accessed so often.