Lo más difícil es conseguir un conjunto de condiciones límite coherentes, lo que requiere una combinación de conjeturas, conocimientos físicos, experiencia previa con problemas relacionados, cálculos detallados y ensayo-error. En resumen, es un poco de arte.
Sin embargo, una vez que se tiene un conjunto de condiciones de contorno (como en su caso las condiciones de contorno NHEK) las cosas son bastante sencillas.
Permítanme denotar la métrica asintótica de fondo por $g$ y las fluctuaciones dependientes del estado por $h$ de modo que cualquier métrica de la forma $g+O(h)$ está permitido por las condiciones de contorno.
Su objetivo es comprobar cuáles son las transformaciones gauge que preservan sus condiciones de contorno. En gravedad pura, deben ser algunos difeomorfismos generados por algún campo vectorial $\xi$ , de tal manera que \begin{equation} L_\xi (g+h) = O(h) \end{equation} donde $L_\xi$ es la derivada de Lie. Como se trata de una ecuación para un tensor simétrico en $D$ dimensiones que se obtienen $D(D+1)/2$ EDP lineales independientes de primer orden para el campo vectorial $\xi$ .
Si has intentado resolver la ecuación anterior, has hecho precisamente lo correcto, lo que espero que responda a parte de tu pregunta. Permíteme ahora abordar la otra parte, es decir, cómo resolver estas EDP.
En muchos ejemplos se puede resolver el campo vectorial más general $\xi$ compatible con la condición anterior simplemente adivinando un Ansatz adecuado y demostrando después que funciona.
En la mayoría de las aplicaciones se tiene alguna expansión en serie de la métrica asintótica en potencias de alguna coordenada radial $r$ (o en exponenciales de $r$ depende de su elección de calibre para la coordenada radial).
Para reducir el desorden, voy a suponer que los distintos componentes tensoriales de $g$ y $h$ se expresan en algunas series de Laurent de $r$ y que $r\to\infty$ corresponde al límite asintótico.
Entonces sólo hay que hacer el mismo tipo de serie de potencias Ansatz para el campo vectorial $\xi$ . Normalmente, todos los componentes del campo vectorial comienzan en $O(1)$ o más pequeño, pero esto no tiene por qué ser así. Si no tienes ninguna idea, entonces haz el Ansatz \begin{equation} \xi^0 = r^{n_0} (\xi^0_0 + \xi^0_1/r + ...) \end{equation} donde las funciones de coeficiente $\xi^0_i$ se permite que dependan de todas las coordenadas de la frontera a priori. Se hace un Ansatz similar para todas las demás componentes del campo vectorial.
La evaluación de las condiciones de la variación de Lie determina entonces los exponentes $n_i$ y puede poner restricciones a las funciones $\xi^i_j$ .
Véase, por ejemplo, el ejercicio (17.1) en los ejercicios de la semana 7 en mi página web de enseñanza http://quark.itp.tuwien.ac.at/~grumil/teaching.shtml para una visita guiada por el AdS estándar $_3$ ejemplo. Si nunca has hecho este tipo de cálculos, te recomiendo que empieces con ese ejemplo antes de romper el NHEK.
[Como nota al margen, nunca he intentado implementar este algoritmo en Mathematica, ya que suelo trabajar en 3 dimensiones donde un cálculo a mano es bastante rápido, pero no veo ninguna razón por la que no debería funcionar en Mathematica].