Slow scipy.sparse.linalg.spilu for very large sparse matrix

I have a large square scipy.sparse.csc_matrix A of size 1_874_602x1_874_602 with 1_623_358 non-zero entries. I am trying to solve a linear system Ax = b using an incomplete LU factorisation.

import scipy.sparse.linalg as spla

N = A.shape[0]
A_iLU = spla.spilu(A) # This line doesn't terminate
M = spla.LinearOperator((N, N), A_iLU.solve)

h, exitcode = spla.gmres(A, b, M=M)

Unfortunately, the call spla.spilu does not end in reasonable time. Is there any obvious optimisation that I am missing? If not, is there another method that I should try?

Thank you for your help!