Pretty sure I found the pattern for the maximum remainder sum for a given 'a.' (for even a, this should be a^2-2*a, for odd a, this should be a^2-a).

When I put this in sigma notation, I was able to get a closed formed solution 1/24(8A^3-6A^2-20A)+2 (the 2 I added when I saw that I was neglecting to add the first case with a=2, e=2).

I also found the multiplicative inverse of 24 modulo (10^9+7) so this computation can be performed quickly.

I put in a snippet of my code below (for A even with e=2). The rest are similar so I imagine the issue will come up below. I got the test cases 1 and 2 working, but the rest give wrong answers.

Any sort of helpful hint would be greatly appreciated.

const unsigned long long LARGE_MOD = pow(10,9)+7; 
const unsigned long long MULT_INVERSE_3 = 333333336; 
const unsigned long long MULT_INVERSE_24 = 41666667;  
const unsigned long long MULT_INVERSE_4 = 250000002; 

//computation of sum broken down into small chunks      ('formula part 1', 'formula part 2', etc.) for readability
long long max_remainder_modSquareAEven(long long a){  
long long form_p1 = ((long long) (8*pow(a,3)))%LARGE_MOD; 
long long form_p2 = ((long long) (6*pow(a,2)))%LARGE_MOD; 
long long first_diff = (form_p1-form_p2+LARGE_MOD)%LARGE_MOD; 

long long form_p3 = (20*a)%LARGE_MOD; 
long long second_diff = (first_diff-form_p3+LARGE_MOD)%LARGE_MOD; 

long long answer = (MULT_INVERSE_24 * second_diff+2)%LARGE_MOD;
return answer;

}