Hi Francois,
I could reproduce the issue with BiCGStab. The reason for the divergence
is the search direction vector becoming degenerate, which happens
because of the simple diagonal structure of the matrix. BiCGStab is
known for this possibility of failure, yet it rarely happens in
practical settings. (I agree that an explicit exception warning about
'inf' would be good.)
One way to get correct results is to switch to the CG method, as your
diagonal matrix is clearly symmetric positive definite.
Best regards,
Karli
On 08/04/2015 12:15 PM, François Fayard wrote:
> #include <iostream>
>
> #define VIENNACL_WITH_OPENMP
> #include <viennacl/scalar.hpp>
> #include <viennacl/vector.hpp>
> #include <viennacl/matrix.hpp>
> #include <viennacl/linalg/bicgstab.hpp>
>
> int main(int argc, const char *argv[]) {
> typedef double ScalarType;
> auto side = std::vector<std::size_t>{1, 3, 5, 10, 20, 30, 100};
>
> for (auto n : side) {
> std::cout << "n: " << n << std::endl;
>
> auto vcl_matrix = viennacl::matrix<ScalarType>(n, n);
> for (auto k = std::size_t{0}; k < n; ++k) {
> vcl_matrix(k, k) = 2.0;
> }
>
> auto vcl_rhs = viennacl::vector<ScalarType>(n);
> for (std::size_t k = 0; k < n; ++k) {
> vcl_rhs[k] = 1.0;
> }
>
> auto vcl_result =
> viennacl::vector<ScalarType>(viennacl::linalg::solve(
> vcl_matrix, vcl_rhs, viennacl::linalg::bicgstab_tag()));
>
> std::cout << 0 << " " << vcl_result[0] << std::endl;
> std::cout << n - 1 << " " << vcl_result[n - 1] << std::endl;
> std::cout << std::endl;
> }
>
> return 0;
> }
------------------------------------------------------------------------------
_______________________________________________
ViennaCL-support mailing list
[email protected]
https://lists.sourceforge.net/lists/listinfo/viennacl-support
