c++ - C++ での dgemm へのアクセス

Question

C++ で行列クラスを実装する必要があり、演算の 1 つは dgemm による行列乗算でなければなりません。私の教授はクラスで C の例を作成しましたが、何らかの理由で C++ で動作させることができません。これは私のヘッダーファイルmatrix.hです:

#include <iostream>
#include <stdlib.h>

extern "C" {
#include "blas.h"
}

[blah blah blah, matrix class here; overloaded * operator will do the matrix      
multiplication]

matrix operator* (const matrix &other){

matrix AxB(Nrows, other.Ncolumns, "(" + name + "*" + "other.name" + ")");

char TRANSA = 'N';
char TRANSB = 'N';
int M = Nrows;
int N = other.Ncolumns;
int K = Ncolumns;
double alpha = 1.;
double beta = 0.;

dgemm_ (&TRANSA,
        &TRANSB,
        &M,
        &N,
        &K,
        &alpha,
        data,
        &M,
        other.data,
        &K,
        &beta,
        AxB.data,
        &M);

return AxB;

[blah blah blah, overloaded = operator here; i'm positive this is not the problem 
 since it works for matrix addition]

主な機能：

#include "matrix.h"

main(){

// entries of A and B are randomized between 0 and 1
matrix A(5,5);
matrix B(5,5);

matrix C = A*B;

}

blas.h ヘッダーファイルは次のようになります。

void dgemm_ (char *TRANSA,
             char *TRANSB,
             int *M,
             int *N,
             int *K,
             double *ALPHA,
             double **A,
             int *LDA,
             double **B,
             int *LDB,
             double *BETA,
             double **C,
             int *LDC);

ここで概説されているサブルーチンを使用します: http://www.netlib.org/blas/dgemm.f

基本的に、C 実装のクラスで行った方法では、(double*) calloc(rows*columns, sizeof(double)) を使用して行列を 1 つの長い配列として構築しました。

私のC ++実装は次のことを行います：

double **data;

data = new double[rows];

for(int i=1; i<=rows; ++i){
    data[i-1] = new double[columns];
}

そして、data[i][j] を使用してインデックスを作成できます。しかし、dgemm サブルーチンは double *A、double *B などを受け取るはずですが、私の行列は double **A などです。これを修正するにはどうすればよいですか?

これは、valgrind から取得したエラーメッセージです。

==18845== Invalid write of size 8
==18845==    at 0x4165C1C: ATL_dgezero (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x4149DA7: ATL_dNCmmJIK (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x415BE8E: ATL_dgemm (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x4182F54: ATL_dptgemm_nt (in /usr/lib/atlas-base/atlas   
/libblas.so.3gf.0)
==18845==    by 0x4183140: ATL_dptgemm (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x40818F6: atl_f77wrap_dgemm_ (in /usr/lib/atlas-base/atlas 
/libblas.so.3gf.0)
==18845==    by 0x4E4E0004: ???
==18845==  Address 0x478d680 is 8 bytes after a block of size 16 alloc'd
==18845==    at 0x402B454: operator new[](unsigned int) (in /usr/lib/valgrind   
/vgpreload_memcheck-x86-linux.so)
==18845==    by 0x8049045: dmatrix::dmatrix(int, int, std::string) (dmatrix.hpp:77)
==18845==    by 0x8049532: dmatrix::operator*(dmatrix const&) (dmatrix.hpp:218)
==18845==    by 0x8048DCB: main (main.cpp:17)
==18845== 
==18845== Invalid write of size 8
==18845==    at 0x4165C1F: ATL_dgezero (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x4149DA7: ATL_dNCmmJIK (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x415BE8E: ATL_dgemm (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x4182F54: ATL_dptgemm_nt (in /usr/lib/atlas-base/atlas 
/libblas.so.3gf.0)
==18845==    by 0x4183140: ATL_dptgemm (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x40818F6: atl_f77wrap_dgemm_ (in /usr/lib/atlas-base/atlas   
/libblas.so.3gf.0)
==18845==    by 0x4E4E0004: ???
==18845==  Address 0x478d678 is 0 bytes after a block of size 16 alloc'd
==18845==    at 0x402B454: operator new[](unsigned int) (in /usr/lib/valgrind 
/vgpreload_memcheck-x86-linux.so)
==18845==    by 0x8049045: dmatrix::dmatrix(int, int, std::string) (dmatrix.hpp:77)
==18845==    by 0x8049532: dmatrix::operator*(dmatrix const&) (dmatrix.hpp:218)
==18845==    by 0x8048DCB: main (main.cpp:17)
==18845== 
==18845== Invalid write of size 8
==18845==    at 0x4165C22: ATL_dgezero (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x4149DA7: ATL_dNCmmJIK (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x415BE8E: ATL_dgemm (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x4182F54: ATL_dptgemm_nt (in /usr/lib/atlas-base/atlas    
/libblas.so.3gf.0)
==18845==    by 0x4183140: ATL_dptgemm (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x40818F6: atl_f77wrap_dgemm_ (in /usr/lib/atlas-base/atlas 
/libblas.so.3gf.0)
==18845==    by 0x4E4E0004: ???
==18845==  Address 0x478d690 is not stack'd, malloc'd or (recently) free'd
==18845== 
==18845== Invalid write of size 8 
==18845==    at 0x4165C25: ATL_dgezero (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x4149DA7: ATL_dNCmmJIK (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x415BE8E: ATL_dgemm (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x4182F54: ATL_dptgemm_nt (in /usr/lib/atlas-base/atlas 
/libblas.so.3gf.0)
==18845==    by 0x4183140: ATL_dptgemm (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x40818F6: atl_f77wrap_dgemm_ (in /usr/lib/atlas-base/atlas 
/libblas.so.3gf.0)
==18845==    by 0x4E4E0004: ???
==18845==  Address 0x478d688 is not stack'd, malloc'd or (recently) free'd
==18845== 
==18845== Invalid read of size 8
==18845==    at 0x4143D60: ATL_dJIK0x0x0NN0x0x0_aX_bX (in /usr/lib/atlas-base/atlas
/libblas.so.3gf.0)
==18845==    by 0x4149A9D: ATL_dNCmmJIK (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x415BE8E: ATL_dgemm (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x4182F54: ATL_dptgemm_nt (in /usr/lib/atlas-base/atlas 
/libblas.so.3gf.0)
==18845==    by 0x4183140: ATL_dptgemm (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x40818F6: atl_f77wrap_dgemm_ (in /usr/lib/atlas-base/atlas 
/libblas.so.3gf.0)
==18845==    by 0x4E4E0004: ???
==18845==  Address 0x478cf80 is not stack'd, malloc'd or (recently) free'd
==18845== 
==18845== Invalid read of size 8
==18845==    at 0x4143D64: ATL_dJIK0x0x0NN0x0x0_aX_bX (in /usr/lib/atlas-base/atlas 
/libblas.so.3gf.0)
==18845==    by 0x4149A9D: ATL_dNCmmJIK (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x415BE8E: ATL_dgemm (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x4182F54: ATL_dptgemm_nt (in /usr/lib/atlas-base/atlas 
/libblas.so.3gf.0)
==18845==    by 0x4183140: ATL_dptgemm (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x40818F6: atl_f77wrap_dgemm_ (in /usr/lib/atlas-base/atlas    
/libblas.so.3gf.0)
==18845==    by 0x4E4E0004: ???
==18845==  Address 0x478d130 is 0 bytes after a block of size 16 alloc'd
==18845==    at 0x402B454: operator new[](unsigned int) (in /usr/lib/valgrind 
/vgpreload_memcheck-x86-linux.so)
==18845==    by 0x8049045: dmatrix::dmatrix(int, int, std::string) (dmatrix.hpp:77)
==18845==    by 0x8048D49: main (main.cpp:6)
==18845== 
==18845== Invalid read of size 8 
==18845==    at 0x4143D4C: ATL_dJIK0x0x0NN0x0x0_aX_bX (in /usr/lib/atlas-base/atlas 
/libblas.so.3gf.0)
==18845==    by 0x4149A9D: ATL_dNCmmJIK (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x415BE8E: ATL_dgemm (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x4182F54: ATL_dptgemm_nt (in /usr/lib/atlas-base/atlas 
/libblas.so.3gf.0)
==18845==    by 0x4183140: ATL_dptgemm (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x40818F6: atl_f77wrap_dgemm_ (in /usr/lib/atlas-base/atlas 
/libblas.so.3gf.0)
==18845==    by 0x4E4E0004: ???
==18845==  Address 0x478d678 is 0 bytes after a block of size 16 alloc'd
==18845==    at 0x402B454: operator new[](unsigned int) (in /usr/lib/valgrind 
/vgpreload_memcheck-x86-linux.so)
==18845==    by 0x8049045: dmatrix::dmatrix(int, int, std::string) (dmatrix.hpp:77)
==18845==    by 0x8049532: dmatrix::operator*(dmatrix const&) (dmatrix.hpp:218)
==18845==    by 0x8048DCB: main (main.cpp:17)
==18845== 
==18845== Invalid write of size 8
==18845==    at 0x4143D93: ATL_dJIK0x0x0NN0x0x0_aX_bX (in /usr/lib/atlas-base/atlas 
/libblas.so.3gf.0)
==18845==    by 0x4149A9D: ATL_dNCmmJIK (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x415BE8E: ATL_dgemm (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x4182F54: ATL_dptgemm_nt (in /usr/lib/atlas-base/atlas 
/libblas.so.3gf.0)
==18845==    by 0x4183140: ATL_dptgemm (in /usr/lib/atlas-base/atlas/libblas.so.3gf.0)
==18845==    by 0x40818F6: atl_f77wrap_dgemm_ (in /usr/lib/atlas-base/atlas   
/libblas.so.3gf.0)
==18845==    by 0x4E4E0004: ???
==18845==  Address 0x478d678 is 0 bytes after a block of size 16 alloc'd
==18845==    at 0x402B454: operator new[](unsigned int) (in /usr/lib/valgrind 
/vgpreload_memcheck-x86-linux.so)
==18845==    by 0x8049045: dmatrix::dmatrix(int, int, std::string) (dmatrix.hpp:77)
==18845==    by 0x8049532: dmatrix::operator*(dmatrix const&) (dmatrix.hpp:218)
==18845==    by 0x8048DCB: main (main.cpp:17)
==18845== 
==18845== Invalid read of size 8
==18845==    at 0x8048DEA: main (main.cpp:19)
==18845==  Address 0x0 is not stack'd, malloc'd or (recently) free'd
==18845== 
==18845== 
==18845== Process terminating with default action of signal 11 (SIGSEGV)
==18845==  Access not within mapped region at address 0x0
==18845==    at 0x8048DEA: main (main.cpp:19)
==18845==  If you believe this happened as a result of a stack
==18845==  overflow in your program's main thread (unlikely but
==18845==  possible), you can try to increase the size of the
==18845==  main thread stack using the --main-stacksize= flag.
==18845==  The main thread stack size used in this run was 8388608.
==18845== 
==18845== HEAP SUMMARY:
==18845==     in use at exit: 3,581 bytes in 39 blocks
==18845==   total heap usage: 49 allocs, 10 frees, 7,760 bytes allocated
==18845== 
==18845== LEAK SUMMARY:
==18845==    definitely lost: 128 bytes in 4 blocks
==18845==    indirectly lost: 0 bytes in 0 blocks
==18845==      possibly lost: 88 bytes in 4 blocks
==18845==    still reachable: 3,365 bytes in 31 blocks
==18845==         suppressed: 0 bytes in 0 blocks
==18845== Rerun with --leak-check=full to see details of leaked memory
==18845== 
==18845== For counts of detected and suppressed errors, rerun with: -v
==18845== ERROR SUMMARY: 111 errors from 9 contexts (suppressed: 0 from 0)
Segmentation fault (core dumped)

dgemm行列乗算を実装しようとするまで、エラーもリークも発生しないため、すべての問題がdgemm実装にあると確信しています

score 1 · Accepted Answer

C バージョンのメモリレイアウトは、C++ バージョンのメモリレイアウトとは異なります。BLAS は C バージョンで使用される種類のレイアウトを想定しているため、C++ バージョンは機能しません。

したがって、C++ バージョンでも大きな 1 次元配列を割り当てる必要があります。2D 配列のインデックスを取得するために operator() をオーバーロードできます。または、製品コードの場合は、 Eigenなどのライブラリを使用します。

score 0 · Accepted Answer

申し訳ありませんが、応答中にEnterキーを押すと、応答が送信されることに気づきませんでした。

インデックスを変更しました：

double *data;

data = new double[rows*columns];

次に、インデックスを作成するたびに、行メジャーのインデックスを作成するためにdata [（i-1）* columns +（j-1）]にインデックスを付けます。コードは間違いなく機能し、0エラーと0メモリリークでコンパイルされます！唯一のことは、私が今完全にぎこちない結果を得るということです。

これは、私が行メジャーを実行しているのに対し、Fortranは列メジャーの順序付けを実行しているという事実と関係がありますか？これをどのように調整しますか？

c++ - C++ での dgemm へのアクセス

2 に答える 2

Related

Reference