[petsc-users] BVNormColumn

[Jose E. Roman]

> El 19 oct 2016, a las 9:54, Jose E. Roman <jroman at dsic.upv.es> escribió:
> 
>> 
>> El 19 oct 2016, a las 0:26, Bikash Kanungo <bikash at umich.edu> escribió:
>> 
>> Hi Jose,
>> 
>> Thanks for the pointers. Here's what I observed on probing it further:
>> 
>> 	• The ||B - B^H|| norm was 1e-18. So I explicitly made it Hermitian by setting B = 0.5(B+B^H). However, this didn't help.
>> 	• Next, I checked for the conditioning of B by computing the ratio of the highest and lowest eigenvalues. The conditioning of the order 1e-9. 
>> 	• I monitored the imaginary the imaginary part of VecDot(y,x, dotXY) where y = B*x and noted that only when the imaginary part is more than 1e-16 in magnitude, the error of "The inner product is not well defined" is flagged. For the first few iterations of orhtogonalization (i.e., the one where orthogonization is successful), the values of VecDot(y,x, dotXY) are all found to be lower than 1e-16. I guess this small imaginary part might be the cause of the error. 
>> Let me know if there is a way to bypass the abort by changing the tolerance for imaginary part.
>> 
>> 
>> 
>> Regards,
>> Bikash
>> 
> 
> There is something wrong: the condition number is greater than 1 by definition, so it cannot be 1e-9. Anyway, maybe what happens is that your matrix has a very small norm. The SLEPc code needs a fix for the case when the norm of B or the norm of the vector x is very small. Please send the matrix to my personal email and I will make some tests.
> 
> Jose

I tested with your matrix and vector with two different machines, with different compilers, and in both cases the computation did not fail. The imaginary part is below the machine precision, as expected. I don't know why you are getting larger roundoff error. Anyway, the check that we currently have in SLEPc is too strict. You can try relaxing it, by editing function BV_SafeSqrt (in $SLEPC_DIR/include/slepc/private/bvimpl.h), for instance with this:

  if (PetscAbsReal(PetscImaginaryPart(alpha))>PETSC_MACHINE_EPSILON && PetscAbsReal(PetscImaginaryPart(alpha))/absal>100*PETSC_MACHINE_EPSILON) SETERRQ1(PetscObjectComm((PetscObject)bv),1,"The inner product is not well defined: nonzero imaginary part %g",PetscImaginaryPart(alpha));

Let us know if this works for you.
Thanks.
Jose