Currently we have:

In [1]: GF(3).unify(ZZ)
Out[1]: ZZ

Note firstly that SymPy's GF is incorrectly named (gh-9544) so really what is called GF(n) here is what is usually called Z/(n) or Z_n (the integers mod n) rather than a finite field (F_p or F_q). I will call refer to Z_n below but what I mean is the domain currently referred to as GF(p) (or FF(p)) in SymPy.

I can sort of see the argument that ZZ is a "superset" of ZZ_n and that usually unify chooses a superset. Really though neither ZZ nor ZZ_n is a superset of each other and ZZ_n does not fit in the scheme of sets that are otherwise subsets of the complex numbers (ZZ, QQ, RR, CC etc).

Neither ZZ_n.unify(ZZ) -> ZZ nor ZZ_n.unify(ZZ) -> ZZ_n will always make sense but I would say that if there is an objective criterion that should be used then it is that if we have A.unify(B) -> C then there should be natural homomorphisms from A -> C and B -> C so that arithmetic in either of A or B maps homorphically to arithmetic in C. Whenever A is a subset of B there will be such a homorphism. In the case of ZZ_n and ZZ there is a homomorphism from ZZ -> ZZ_n but not vice versa so I think that the only sensible result here is ZZ_n.unify(ZZ) -> ZZ_n.

We also have this:

In [2]: GF(3).unify(GF(2))
Out[2]: GF(3)

In this case there are no nontrivial homomorphisms in any direction because of the different characteristic so this should just be an error.

See also #19070 (comment) and gh-20327