viennacl с комплексными числами

Question

Я использую viennacl для решения линейной системы уравнений (AX = B) с графической картой. Кроме того, код использует armadillo.

Моя система уравнений имеет комплексные числа. Поэтому возникает вопрос: могу ли я решить систему уравнений (с комплексными числами) с помощью viennacl?

Выше приведен пример рабочего кода с действительными числами.

// System headers 
#include <iostream> 

// Armadillo headers (disable BLAS and LAPACK to avoid linking issues) 
#define ARMA_DONT_USE_BLAS 
#define ARMA_DONT_USE_LAPACK 
#include <armadillo> 
#include <complex> 

#define VIENNACL_WITH_ARMADILLO 1 

// ViennaCL headers 
#include "viennacl/linalg/cg.hpp" 
#include "viennacl/linalg/bicgstab.hpp" 
#include "viennacl/linalg/gmres.hpp" 
#include "viennacl/io/matrix_market.hpp" 
#include "vector-io.hpp" 

//using namespace arma; 
using namespace viennacl::linalg; 
using namespace std; 

typedef arma::mat armat; 
typedef arma::vec arvec; 

typedef complex<double> dcmplx; 

int main(void) 
{ 

    int N = 500; 
    armat A(N,N); 
    A.randu(); 
    arvec B(N); 
    B.randu(); 
    arvec X(N); 
    arvec residual(N); 

    viennacl::matrix<double> vcl_A(N, N); 
    viennacl::vector<double> vcl_B(N); 
    viennacl::vector<double> vcl_X(N); 
    viennacl::vector<double> vcl_result(N); 

    viennacl::copy(A, vcl_A); 
    viennacl::copy(B, vcl_B); 
    viennacl::copy(X, vcl_X); 

    std::cout << "----- Running GMRES -----" << std::endl; 
    vcl_X = viennacl::linalg::solve(vcl_A, vcl_B,    viennacl::linalg::gmres_tag()); 

    viennacl::copy(vcl_A, A); 
    viennacl::copy(vcl_B, B); 
    viennacl::copy(vcl_X, X); 

    residual = A * X - B; 
    cout << "Relative residual: " << norm(residual)/norm(B) << endl; 
} 

Комплексная версия кода:

#include <iostream> 

// Armadillo headers (disable BLAS and LAPACK to avoid linking issues) 
#define ARMA_DONT_USE_BLAS 
#define ARMA_DONT_USE_LAPACK 
#include <armadillo> 
#include <complex> 

#define VIENNACL_WITH_ARMADILLO 1 

// ViennaCL headers 
#include "viennacl/linalg/cg.hpp" 
#include "viennacl/linalg/bicgstab.hpp" 
#include "viennacl/linalg/gmres.hpp" 
#include "viennacl/io/matrix_market.hpp" 
#include "vector-io.hpp" 


//using namespace arma; 
using namespace viennacl::linalg; 
using namespace std; 

typedef arma::cx_mat armat; 
typedef arma::cx_vec arvec; 

typedef complex<double> dcmplx; 

int main(void) 
{ 

int N = 500; 
armat A(N,N); 
A.randu(); 
arvec B(N); 
B.randu(); 
arvec X(N); 
arvec residual(N); 

viennacl::matrix<dcmplx> vcl_A(N, N); 
viennacl::vector<dcmplx> vcl_B(N); 
viennacl::vector<dcmplx> vcl_X(N); 
viennacl::vector<dcmplx> vcl_result(N); 

viennacl::copy(A, vcl_A); 
viennacl::copy(B, vcl_B); 
viennacl::copy(X, vcl_X); 

std::cout << "----- Running GMRES -----" << std::endl; 
vcl_X = viennacl::linalg::solve(vcl_A, vcl_B, viennacl::linalg::gmres_tag()); 

viennacl::copy(vcl_A, A); 
viennacl::copy(vcl_B, B); 
viennacl::copy(vcl_X, X); 

residual = A * X - B; 
cout << "Relative residual: " << norm(residual)/norm(B) << endl; 



std::cout << "----- Running BiCGStab -----" << std::endl; 
vcl_X = viennacl::linalg::solve(vcl_A, vcl_B, viennacl::linalg::bicgstab_tag()); 

viennacl::copy(vcl_A, A); 
viennacl::copy(vcl_B, B); 
viennacl::copy(vcl_X, X); 

residual = A * X - B; 
cout << "Relative residual: " << norm(residual)/norm(B) << endl; 

std::cout << "----- Running CG -----" << std::endl; 
vcl_X = viennacl::linalg::solve(vcl_A, vcl_B, viennacl::linalg::cg_tag()); 

viennacl::copy(vcl_A, A); 
viennacl::copy(vcl_B, B); 
viennacl::copy(vcl_X, X); 

residual = A * X - B; 
cout << "Relative residual: " << norm(residual)/norm(B) << endl; 

} 

Answer 1

ViennaCL в настоящее время не поддерживает комплексные числа. Основная техническая причина заключается в том, что OpenCL не обеспечивает поддержку комплексных чисел. Конечно, эмулировать сложную арифметику через реальную арифметику, мы не хотели идти по этому пути и (ошибочно?) Надеялись, что стандартизация для комплекса скоро появится.