0

私はかなり複雑な関数、たとえば func() を計算しようとしています - Fortran のいくつかの 2 次元配列のいくつかの加算、減算、乗算、除算、および三角関数を含みます。各 func() がその行と列の位置に対して独立しているという点で、計算は非常に並列です。各行列のサイズは数ギガバイトで、引数として約 12 個あります。

インテル® MKL 関数 (--mkl-parallel の呼び出し)、特に加算、減算、除算などの VML 関数を利用したいと考えています。

例: func(x,y,z) = x*y+cos(z*xx) ここで、x,y,z は数 GB の 2 次元配列です。

VML 関数に関してですが、より使い慣れた二項演算子を使用します。私の問題では、原則として、「+」や「*」などのすべての二項演算子を、引数を ?vadd(x,y) として受け取る二項関数に変換する必要があることがわかります。もちろん、これは非常に扱いにくく、大きな式の場合は見苦しくなります。Fortran で MKL/VML バージョンを優先的に使用するために、「+」、「-」などの二項算術演算子をオーバーロードする方法はありますか。例がいいでしょう！ありがとう！

4

2 に答える 2

1

私はこの答えが少し話題から外れていることを知っています。

すべての操作は要素単位であり、操作は単純であるfunc()ため、メモリ帯域幅が制限されたタスクになる可能性があります。この場合、VML の使用はパフォーマンスを最大化するための適切な選択ではない可能性があります。

各アレイのサイズが 10GB であると仮定すると、次のように VML を使用するには、少なくとも 9 x 10GB の読み取りと 5 x 10GB の書き込みが必要になります。

func(...) {
    tmp1=x*z
    tmp1=tmp1-x;
    tmp1=cos(tmp1);
    tmp2=x*y;
    return tmp1+tmp2;
}

ここで、すべての操作が 2 次元配列に対してすべてオーバーロードされます。

代わりに、次のアプローチの方がメモリ アクセスがはるかに少ない (3 x 10GB の読み取りと 1 x 10GB の書き込み) ため、高速になる可能性があります (疑似コード)。

$omp parallel for
for i in 1 to m
    for j in 1 to n
        result(i,j)= x(i,j)*y(i,j)+cos(z(i,j)*x(i,j)-x(i,j));
    end
end
于 2013-10-23T11:51:20.870 に答える
0

2 つのベクトルの加算を示す小さな例を作成しました。もう MKL をインストールしていないのでSAXPY、BLAS のコマンドを使用しました。原則は同じはずです。

最初に、適切な定義でモジュールを定義します。私の場合、これは実際の配列をデータ型に保存するための割り当て (array変数に直接アクセスすることもできるため、これは便利な関数にすぎません) と追加の定義です。どちらも、+演算子と=代入に対する新しいオーバーロードです。

プログラムでは、3 つのフィールドを定義します。そのうちの 2 つに乱数が割り当てられ、追加されて 3 番目のフィールドが取得されます。次に、最初の 2 つのフィールドが特殊変数に格納され、この加算の結果がこの型の 3 番目の変数に格納されます。

最後に、配列に直接アクセスして結果を比較します。カスタム データ型から同じデータ型への割り当ては既に定義されていることに注意してください (たとえばffield3 = ffield1、既に定義されています)。

私のモジュール:

MODULE fasttype

    IMPLICIT NONE

    PRIVATE

    PUBLIC :: OPERATOR(+), ASSIGNMENT(=)

    TYPE,PUBLIC :: fastreal
        REAL,DIMENSION(:),ALLOCATABLE :: array
    END TYPE

    INTERFACE OPERATOR(+)
        MODULE PROCEDURE fast_add
    END INTERFACE

    INTERFACE ASSIGNMENT(=)
        MODULE PROCEDURE fast_assign
    END INTERFACE

CONTAINS

    FUNCTION fast_add(fr1, fr2) RESULT(fr3)
        TYPE(FASTREAL), INTENT(IN) :: fr1, fr2
        TYPE(FASTREAL) :: fr3
        INTEGER :: L

        L = SIZE(fr2%array)
        fr3 = fr2       

        CALL SAXPY(L, 1., fr1%array, 1, fr3%array, 1)
    END FUNCTION

    SUBROUTINE fast_assign(fr1, r2)
        TYPE(FASTREAL), INTENT(OUT) :: fr1
        REAL, DIMENSION(:), INTENT(IN) :: r2
        INTEGER :: L

        IF (.NOT. ALLOCATED(fr1%array)) THEN
            L = SIZE(r2)
            ALLOCATE(fr1%array(L))
        END IF

        fr1%array = r2
    END SUBROUTINE
END MODULE

私のプログラム:

PROGRAM main
    USE fasttype
    IMPLICIT NONE

    REAL, DIMENSION(:), ALLOCATABLE :: field1, field2, field3
    TYPE(fastreal) :: ffield1, ffield2, ffield3

    ALLOCATE(field1(10),field2(10),field3(10))
    CALL RANDOM_NUMBER(field1)
    CALL RANDOM_NUMBER(field2)

    field3 = field1 + field2

    ffield1 = field1
    ffield2 = field2

    ffield3 = ffield1 + ffield2

    WRITE(*,*) field3 == ffield3%array
END PROGRAM
于 2013-10-23T13:26:46.113 に答える