It would be fantastic if `quadeig` is supported.

Without `quadeig` is supported, the toolbox is scarcely useful for a large number of people who work on quadratic eigenvalue problems.

Thank you very much for your request. We have added "quadeig" to our TODO list.

Could you please provide a bit more details regarding your request? This would help us to add support for "quadeig" faster.

- Is there simple example of quadratic eigenvalue problem to demonstrate the need for high-precision (e.g. where double precision fails due ill-conditioning, etc.)?

- In general, what kind of applied & research problems would high-precision "quadeig" enable scientists to solve?

- Are these kind of quad-eig problems common (in research and in applications)?

- What codes researchers in the field are using now? (e.g. Francoise Tisseur's library:https://www.maths.manchester.ac.uk/~ftisseur/misc/)

This information would help us a lot. 

Thank you in advance.

Thanks for the response! And sorry for not having answered your questions sooner!

- Is there simple example of quadratic eigenvalue problem to demonstrate the need for high-precision (e.g. where double precision fails due ill-conditioning, etc.)?

Yes, smashingly. In the NLEVP library (can be downloaded from Tisseur's website using the link you posted in your last reply), there is a famous quadratic eigenvalue problem called 'damped beam' as it is highly ill-conditioned.

- In general, what kind of applied & research problems would high-precision "quadeig" enable scientists to solve?

On the one hand, there are a large number of problems like 'damped beam' from engineering and other disciplines which are of ill-condition and, therefore, could benefit from extended precision computing by obtaining much more accurate or at least meaning results. On the other hand, people need extended precision arithmetic to generate reliable results for their research (for comparison purposes).

- Are these kind of quad-eig problems common (in research and in applications)?

Yes. Otherwise, you wouldn't have seen so many citations of Tisseur's papers (and papers by many others too) on quadratic eigenvalue problems and the NLEVP library.

- What codes researchers in the field are using now? (e.g. Francoise Tisseur's library:https://www.maths.manchester.ac.uk/~ftisseur/misc/)

It is indeed the code on this webpage that are mostly used/cited in the research of quadratic eigenvalue problems.

Thank you for your response.

I have done some quick tests using ill-conditioned problems from NLEVP library. The 'polyeig' provided in our toolbox works quite well on them.

I wrote simple function 'test_polyeig' which returns absolute error of eigen-decomposition for the the problem specified (see its full code at the bottom) for 'double' and extended precision ('mp') cases:

>> test_polyeig('power_plant','double')

ans =

   1.8567e+05

>> test_polyeig('power_plant','mp')

ans = 

    4.091947679918651639316533910037398e-13


>> test_polyeig('cd_player','double')

ans =

   1.1410e+04

>> test_polyeig('cd_player','mp')

ans = 

    4.953565354331632662743256645399015e-15

>> test_polyeig('damped_beam','double')  % use 'double' precision

ans =

   34.2276

>> test_polyeig('damped_beam','mp') % use extended precision

ans = 

    7.547834125835157737664345305865886e-17

Our 'polyeig' gives much more accurate result out-of-the-box. The only possible downside compared to 'quadeig' is that 'polyeig' is (much) slower. But it is perfectly usable for research right now. 

Could you please justify why you think "...without 'quadeig' is supported, the toolbox is scarcely useful for a large number of people who work on quadratic eigenvalue problems." Quadratic (and other) eigenproblems can be handled by toolbox right now using 'polyeig'. Or am I missing something? 

The code of 'test_polyeig' below. It supports 'double' and 'mp' computations:

function abserror = test_polyeig(name,float_t,varargin)

[COEFFS,~] = nlevp(name,varargin{:});

p = length(COEFFS);
A = cell(p,1);
for k = 1:p
    A{k} = cast(COEFFS{k},float_t);
end

n = size(A{1},1);

[X,E] = polyeig(A{:});
error = zeros(n*(p-1),float_t);

for j = 1:n*(p-1)
   lambda = E(j);
   x = X(:,j);
   
   if lambda ~=Inf
      T = A{1};
      for k = 2:p
          T = T + lambda^(k-1) * A{k};
      end
      
      error(j) = norm(T*x,Inf);
   else
      error(j) = norm(A{1}*x,Inf);       
   end
end

abserror = norm(error,Inf)

P.S.

To get even higher accuracy, NLEVP code must be converted to support 'mp' computations.

Because right now, NLEVP computes matrices in 'double' precision and then we use 'mp' to compute eigenvalues.

If matrices would be computed with high precision from the onset, than overall accuracy would be higher.

Any follow up on this?