Efficient way to solve large Lyapunov equations and inequalities in MATLAB?

I'm currently facing a problem which requires solving large (somewhat) sparse Lyapunov equations of the form A'*X + X*A = Q, where Q is symmetric and positive-definite and A,Q are matrices of the size ca. 15000 x 15000. Unfortunately lyap(A',Q) and lyapchol(A',Q) seem (far) too slow for my equations.

Now technically I don't actually want to solve the Lyapunov equations, but I'm actually interested in finding a matrix X>0 such that A'*X + X*A <=0, but for that problem as well MATLAB's native algorithms seem to fail (for that size), that is if the proper function from the LMI toolbox doesn't just crash my MATLAB.

Now I've also tried the appropriate function from Lyapack to solve Lyapunov equations (with appropriately chosen Q) which - while no doubt efficient for certain kinds of problems - didn't provide any improvement. Does anyone have a good idea what I could do to solve this problem (in MATLAB) ?

Thanks in advance.

[Jon]

Have you already solved this problem? I encountered similar problem.