A similar problem, n-body simulation*, has n² gravitational interactions. You will similarly hit a wall if you try to do it with a dense n² matrix. However, there's a hierarchical solution that takes advantage of the sparsity and exponential decay, and can solve it in (n log n) with an imperceptibly low loss of precision.

Social interactions are sparse, and group interactions can be optimized with clustering. Fine-grained simulation of the entire society is such a massive chaotic problem with so many variables, that some loss of precision from clustering is completely insignificant compared to the inevitable simplifications you'll have to make in the design of the model itself.

* I mean the naive one with a fixed timestemp, not trying to solve chaos.