[Take four, simplified again, and using more narrow hypothesis for security]
Improving on the line of thought in that other answer, we will craft two efficiently computable functions $H_1$ and $H_2$ each accepting two arbitrary strings as input, with output a fixed-size string, say 512 bit; and (I conjecture) indistinguishable from a random function (under the assumption that some parameters are secret) except for this superset of the desired property: $$\forall(i,j)\in\{1,2\}^2,\forall(s_1,s_2,x), H_i(s_1,H_j(s_2,x))=H_j(s_2,H_i(s_1,x))$$
We'll use as building blocks three 512-bit PRFs $P_0()$, $P_1()$ and $P_2()$ accepting an arbitrary bit string as input, and one 512-bit PRP $E()$ with $D()$ its reverse function, such that $D(E(a))=a$.
We can build the PRFs as $P(m)=\mathtt{HMAC}_\mathtt{SHA512}(k,m)$ and three arbitrary distinct $k$. And by this famous result we can build the PRP as the four-rounds Feistel cipher with round functions $\mathtt{HMAC}_\mathtt{SHA256}(k,m)$ and four arbitrary distinct $k$.
For $i\in\{1,2\}$, define $H_i(s,x)$ as $$H_i(s,x)=E\big(\text{ }D(P_i(s))+D(x)\bmod{2^{512}}\text{ }\big)\text{ when }x\text{ is a 512-bit string,}$$ $$H_i(s,x)=E\big(\text{ }D(P_i(s))+D(P_0(x))\bmod{2^{512}}\text{ }\big)\text{ otherwise.}$$
The desired property $H_i(s_1,H_j(s_2,x))=H_j(s_2,H_i(s_1,x))$ follows from commutativity and associativity of addition in $\mathbb Z_{2^{512}}$, and the $H_i$ functions appears much like random oracles, except for consequences of that property, to an adversary ignoring the PRP and the PRF $P(0)$ (but possibly knowing the other PRFs).
By querying the functions/oracles $H_i$, it seems impossible to learn something on the value internally added, or internally add 0, other than by chance. Computing $H_1(s,H_1(s,\cdots H_1(s,H_1(s,x))\cdots))$, where $x$ is a 512-bit string and there are $2^n$ iterations, has odds only $2^{n-511}$ to return to $x$; if we had used XOR instead of modular addition, we'd have $H_1(s,H_1(s,x))=x$.
Note: to an adversary knowing the PRP and PRF, the problem of finding collisions on $H_1$ unrelated to the property is tractable: it is a knapsack problem in $\mathbb Z_{2^{512}}$ with unlimited supply of random values to choose from.
Update: I realize that if $P_0$ is made public, a collision can be found without using the desired property: for any $(x,s)$ with $x$ not 512-bit, $H_j(s,x)=H_j(s,P_0(x))$. It can't be fixed by replacing $D(P_0(x))$ by $P_0(x)$, which allows finding collisions by solving a knapsack.
Open problems:
- Can we make a scheme secure yet fully public?
- Can this be meaningfully proven secure (or broken)?
- We do not need a full group operation (on the contrary, the neutral element and existence of an opposite are a risk), any commutative semi-group will do; is anything more suitable than $\mathbb Z_{2^{512}}$?
- Is $H_i(s_1,H_j(s_2,x))=H_j(s_2,H_i(s_1,x))$ a consequence of $H_1(s_1,H_2(s_2,x)) = H_2(s_2,H_1(s_1,x))$?
- If not, can we make a construction with only the originally asked property?