复数的维也纳

Posted 2023-02-21

【中文标题】复数的维也纳

【英文标题】：viennacl with complex numbers

【发布时间】：2015-09-06 21:06:26

【问题描述】：

我正在使用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;

【问题讨论】：

您是否尝试将它与std::complex 一起使用？您是否遇到任何特定的错误或问题，我们可能会为您提供帮助？正如你的问题似乎过于宽泛。

是的，我做到了。我包括上面的例子。我不知道这里是否有人可以帮助我，但无论如何我都会尝试。 ;)

您观察到的错误/问题是什么？请在您的问题中逐字记录所有内容。

【参考方案1】：

ViennaCL 目前不支持复数。主要的技术原因是 OpenCL 本身并不支持复数。虽然通过实数算术模拟复数算术当然是可能的，但我们不愿意走这条路，并且（错误地？）希望复数的标准化很快就会到来。