acb_mat_eig_simple_rump segmentation faults when R is NULL but not L

[Joel-Dahne]

It can be reproduced with the following code

#include "acb.h"
#include "acb_mat.h"

int
main()
{
  slong prec = 64;
  mag_t tol;
  acb_ptr E, E_approx;
  acb_mat_t L, R, A, L_approx, R_approx;

  mag_init(tol);

  E = _acb_vec_init(3);
  E_approx = _acb_vec_init(3);

  acb_mat_init(L, 3, 3);
  acb_mat_init(R, 3, 3);
  acb_mat_init(A, 3, 3);
  acb_mat_init(L_approx, 3, 3);
  acb_mat_init(R_approx, 3, 3);

  // Setup matrix
  acb_set_d_d(acb_mat_entry(A, 0, 0), 0.6873474041954415, 0.8982563031334123);
  acb_set_d_d(acb_mat_entry(A, 0, 1), 0.7282180564881044, 0.3029712969740874);
  acb_set_d_d(acb_mat_entry(A, 0, 2), 0.07360652513458521, 0.8585014523679579);
  acb_set_d_d(acb_mat_entry(A, 1, 0), 0.000835810121029068, 0.7583002736998279);
  acb_set_d_d(acb_mat_entry(A, 1, 1), 0.9256166870757694, 0.8854763478184455);
  acb_set_d_d(acb_mat_entry(A, 1, 2), 0.5363310989411239, 0.3031103325817668);
  acb_set_d_d(acb_mat_entry(A, 2, 0), 0.07387174694790022, 0.2319572749472405);
  acb_set_d_d(acb_mat_entry(A, 2, 1), 0.4050436025621329, 0.5769840251057949);
  acb_set_d_d(acb_mat_entry(A, 2, 2), 0.20226010388885896, 0.5119507333628952);

  acb_mat_printd(A, 10);

  acb_mat_approx_eig_qr(E_approx, L_approx, R_approx, A, NULL, 0, prec);
  for (int i = 0; i < 3; i++) {
    acb_printd(E_approx + i, 10);
    flint_printf("\n");
  }
  acb_mat_printd(L_approx, 10);

  /* This one works with R equal to NULL */
  acb_mat_eig_simple_vdhoeven_mourrain(E, L, NULL, A, E_approx, R_approx, prec);
  flint_printf("vdhoeven_mourrain worked\n");
  /* This one works when R us not NULL*/
  acb_mat_eig_simple_rump(E, L, R, A, E_approx, R_approx, prec);
  flint_printf("rump worked with R != NULL\n");
  /* But crashes when R is NULL*/
  acb_mat_eig_simple_rump(E, L, NULL, A, E_approx, R_approx, prec);
  flint_printf("rump worked with R == NULL\n");
}

Giving

> ./eig_bug
[(0.6873474042 + 0.8982563031j)  +/-  (0, 0j), (0.7282180565 + 0.302971297j)  +/-  (0, 0j), (0.07360652513 + 0.8585014524j)  +/-  (0, 0j)]
[(0.000835810121 + 0.7583002737j)  +/-  (0, 0j), (0.9256166871 + 0.8854763478j)  +/-  (0, 0j), (0.5363310989 + 0.3031103326j)  +/-  (0, 0j)]
[(0.07387174695 + 0.2319572749j)  +/-  (0, 0j), (0.4050436026 + 0.5769840251j)  +/-  (0, 0j), (0.2022601039 + 0.5119507334j)  +/-  (0, 0j)]
(1.301533099 + 1.807658938j)  +/-  (0, 0j)
(0.6822580497 + 0.3940775247j)  +/-  (0, 0j)
(-0.1685669533 + 0.09394692117j)  +/-  (0, 0j)
[(0.4982598789 - 0.339235072j)  +/-  (0, 0j), (0.2559573928 - 0.6996723055j)  +/-  (0, 0j), (0.3050282979 - 0.2639139148j)  +/-  (0, 0j)]
[(-0.5412099397 + 0.3875106093j)  +/-  (0, 0j), (0.4472809802 - 0.05279366497j)  +/-  (0, 0j), (-0.07047395923 - 0.7024042544j)  +/-  (0, 0j)]
[(0.2395294487 + 0.06731000386j)  +/-  (0, 0j), (-0.192734642 + 0.4324182207j)  +/-  (0, 0j), (-0.1880734264 - 0.8237664918j)  +/-  (0, 0j)]
vdhoeven_mourrain worked
rump worked with R != NULL
[1]    485172 segmentation fault (core dumped)  ./eig_bug

The issue seems to be the part

if (L != NULL)
    {
        if (!result || !acb_mat_inv(L, R, prec))
            acb_mat_indeterminate(L);
    }

in eig_simple_rump.c which uses R without checking if it's NULL first.