Can this be solved even faster?

It is more efficient to first pick the integer partitions whose absolute values sum up to n before generating the permutations.

dimNM3[n_, m_] := Total[
   Map[
    Length@*Permutations,
    Pick[#, Abs[#].ConstantArray[1, 4], n] &[
     IntegerPartitions[m, {4}, Range[-n, n]
      ]
     ]
    ]
   ];

m = 20;
n = 40;
dimNM1[n, m] // AbsoluteTiming
dimNM2[n, m] // AbsoluteTiming
dimNM3[n, m] // AbsoluteTiming

{0.116977, 3802}

{0.995365, 3802}

{0.005579, 3802}

ClearAll[num];

num[n_, m_] /; OddQ[n + m] = 0;
num[n_, n_] := Binomial[n + 3, 3];
num[n_, m_] /; OddQ[n] := With[{z = Ceiling[m/2]}, (5*n^2 + 3)/2 + 2 z - (2 z^2)];
num[n_, m_] /; EvenQ[n] := With[{z = Ceiling[m/2]}, (5*n^2 + 4)/2 - (2 z^2)];

Testing vs fastest answer here at writing (Henrik Schumacher):

stop = 100;

res = Table[{n, m, dimNM3[n, m]}, {n, 1, stop}, {m, 1, n}]; // AbsoluteTiming//First
res2 = Table[{n, m, num[n, m]}, {n, 1, stop}, {m, 1, n}]; // AbsoluteTiming//First

res == res2

169.203

0.0219434

True

Large cases are a non-issue:

num[123423456, 123412348] // AbsoluteTiming

{0.0000247977, 30468069908023290}

Some quick timings:

Sorry for not knowing much Mathematica, but I have a Python solution you might be able to follow. I'm putting this on the community wiki for anyone who wants to translate it.

def count_solutions(Nm, Mm):
    firsthalves = dict()
    for m1 in range(-Nm,Nm+1):
        for m2 in range(-Nm,Nm+1):
            m = m1+m2
            n = abs(m1)+abs(m2)
            key = (m,n)
            if key in firsthalves:
                firsthalves[key] += 1
            else:
                firsthalves[key] = 1

    solutions = 0
    for m3 in range(-Nm,Nm+1):
        for m4 in range(-Nm,Nm+1):
            m = m3+m4
            n = abs(m3)+abs(m4)
            key = (Mm-m, Nm-n)
            if key in firsthalves:
                solutions += firsthalves[key]
    return solutions

This is a meet in the middle strategy. I enumerate all the possible $m1,m2$ combinations and record how many times each $m1+m2,|m1|+|m2|$ combination occurs in a dictionary.

Then I go through all the possible $m3,m4$ combinations and for each combination I calculate the necessary $m1+m2,|m1|+|m2|$ combination to make $Mm,Nm$, and I refer to the dictionary to find out how many $m1,m2$ combinations can make that.

The difference is that you go through the $m1,m2$ combination then the $m3,m4$ combinations, and the number of operations is roughly a square root of going through every $m1,m2,m3,m4$ combination. You should be able to solve for $Nm = 1000,Mn = 0$ in a few seconds.