Roughly defined, the GKZ (Gelfand-Kapranov-Zelevinsky) systems are classes of differential equations that can be solved in terms of generalised hypergeometric functions - for more details on the subject look e.g. these lecture notes...
A simple example of GKZ is the vieled Laplace equation:
$$\left( \frac{\partial}{\partial x_1} \frac{\partial}{\partial x_2} + \frac{\partial}{\partial x_3} \frac{\partial}{\partial x_4} \right) \Phi = 0$$
which is satisfied by any combination of:
$$\Phi [a,b;c] = \frac{x_3{}^{c-1}}{x_1{}^{a} x_2^{b}} \> {}_2 \mathfrak{F}_1 \left( a, b; c; -\frac{x_3 x_4}{x_1 x_2} \right) : c \not\in \mathbb{Z}_{-}$$
My question is how do we approach problems in which $\Phi$ is multivariable: say $\Phi = \Phi(x_0, x_1, \dots, x_n; y_{11}, y_{12}, \dots, y_{(n-1),n}, y_{nn})$ and for any $i,j : 1 \leq i \leq j \leq n$ it satisfies (among other things):
$$\left( \frac{\partial}{\partial y_{ij}} \frac{\partial}{\partial x_0} + \frac{\partial}{\partial x_i} \frac{\partial}{\partial x_j} \right) \Phi = 0 \, ?$$
In other words, do you use an ansatz or try something more systematic?
PS Such types of systems of diff eqns appear for Feynman integrals...