Fortran DLL を呼び出して数値計算を行う C# アプリケーションがあります。DLL 呼び出しは、Fortran コンソール アプリから同じ計算を呼び出すよりも少なくとも 10 倍遅いことが判明したため、コードを C# に移植し始めました。移植プロセスは、コードをコピーし、構文を 1 行ずつ修正することです。したがって、C# と Fortran は非常によく似ています。Fortran の共通ブロックは、クラスのフィールドになります。いくつかのコア ルーチンを移植した後、テストを開始したところ、倍精度の C# コードは倍精度の Fortran よりも 30 倍遅い (単精度の Fortran よりも 50 倍遅い) ことがわかりました。テスト コードを 100 回ループして、C# JIT コンパイラのオーバーヘッドの影響を最小限に抑えます。
このコードでは、SQRT や LOG10 などの複雑な算術関数や複雑な関数が使用されています。その計算を処理する C# 構造体を提供しました。ここに問題があるのではないかと思いますが、VS2010 Pro にはプロファイラーがないようで、よくわかりません。
Fortran コンパイルは、Intel の最新バージョンです。C# または Fortran コードに対して特別な最適化は行っていません。リリースバージョンを使用して時間を比較しました。私は 2 コアの CPU しか持っていないので、並列化はおそらくあまり役に立ちません。
このコードを高速化する方法について、いくつかの提案を使用できます。
これが私の方法の1つです。C# は、私が書いたのではない Fortran と同じように見えます。
public void Forwx(double rshldr, double rbed, double[][] resdep, double[]toolparm)
{
var zex = new Complex[70];
var zey = new Complex[70];
var px2 = new Complex[70];
var px4 = new Complex[70];
var rlt = new double[2];
var trsmt = new double[2];
var fr = new double[2];
var dnc = -0.02;
var factr = 26332.65;
var rh2 = Math.Max(0.1, rbed);
var rh1 = Math.Max(0.1, rshldr);
const double e1 = 1.0;
const double e2 = 1.0;
const double er = 0.818730753077982;
const double re = 1.0 / er;
var ii = Complex.I;
const double pi = Math.PI;
const double eps0 = 8.854e-12;
const double amu0 = 4.0e-7 * pi;
for (var ktool = 3; ktool <= 6; ktool++)
{
if (ktool == 3) // Integrated 2MHz
{
dnc = -0.02;
rlt[0] = 0.2794;
rlt[1] = -0.2794;
trsmt[0] = 0.904875;
trsmt[1] = -0.904875;
fr[0] = 2000.0;
factr = 26332.65;
}
if (ktool == 4) // Integrated 400kHz
{
dnc = -0.02;
rlt[0] = 0.2794;
rlt[1] = -0.2794;
trsmt[0] = 0.904875;
trsmt[1] = -0.904875;
fr[0] = 400.0;
factr = 26811.866;
}
if (ktool == 5) // Option 5 20kHz
{
dnc = -0.1;
rlt[0] = 0.0;
rlt[1] = 0.0;
trsmt[0] = 5.75;
trsmt[1] = 5.75;
fr[0] = 20.0;
factr = 26811.866 * 2.516 * toolparm[1] / 0.28e8;
}
if (ktool == 6) // Option 6 50kHz
{
dnc = -0.1;
rlt[0] = 0.0;
rlt[1] = 0.0;
trsmt[0] = 5.75;
trsmt[1] = 5.75;
fr[0] = 50.0;
factr = 26811.866 * 6.291 * toolparm[2] / 0.7e8;
}
var r1 = trsmt[0] - rlt[0];
var r2 = trsmt[0] - rlt[1];
var omega = 2000.0 * pi * fr[0];
var k12 = omega*amu0*(omega*e1*eps0 + ii/rh1);
var k22 = omega*amu0*(omega*e2*eps0 + ii/rh2);
var krat = (k22 - k12)/k12;
for (var iz = 0; iz < 601; iz++)
{
var recx1 = new Complex(0.0, 0.0);
var rx1 = new Complex(0.0, 0.0);
var recy1 = new Complex(0.0, 0.0);
var ry1 = new Complex(0.0, 0.0);
var lam = new Complex(3.01517934056e-04 / (Math.Pow(er, 5) * r1));
Complex c1;
Complex c2;
for (var i = 0; i < 70; i++)
{
if (iz == 0)
{
lam = lam * re;
var lam2 = lam * lam;
var p11 = lam2 - k12;
var p1 = Complex.Sqrt(p11);
var p22 = lam2 - k22;
var p2 = Complex.Sqrt(p22);
zex[i] = Complex.Exp(dnc * p2);
zey[i] = Complex.Exp(dnc * p1);
c1 = p2 * k12;
c2 = p1 * k22;
var t3 = lam / p2;
var t2 = t3 * (c1 - c2) / (c1 + c2);
var q2 = lam * krat * (t2 + t3) / (p1 + p2);
px2[i] = (lam2 * q2 + lam * p2 * t2);
px4[i] = px2[i];
}
else
{
px2[i] = px2[i] * zex[i];
px4[i] = px4[i] * zey[i];
}
recx1 = recx1 + a1[i] * px2[i];
recy1 = recy1 + a1[i] * px4[i];
rx1 = rx1 + px2[i] * as1i[i];
ry1 = ry1 + px4[i] * as1i[i];
}
if (ktool <= 4)
{
c1 = recx1*r1;
c2 = rx1*r2;
c2 = c2 - Math.Pow(r1/r2,3)*c1;
resdep[12 - ktool][iz + 600] = c2.Re*factr;
c1 = recy1*r1;
c2 = ry1*r2;
c2 = c2 - Math.Pow(r1 / r2,3) * c1;
resdep[12 - ktool][600 - iz] = c2.Re*factr;
}
else
{
c1 = recx1*r1;
//c2 = rx1*r2;
//c2 = c2 - Math.Pow(r1 / r2,3) * c1;
resdep[ktool + 5][iz + 600] = c1.Re * factr;
c1 = recy1*r1;
//c2 = ry1*r2;
//c2 = c2 - Math.Pow(r1 / r2,3) * c1;
resdep[ktool + 5][600 - iz] = c1.Re * factr;
}
}
}
}
複雑な構造体のメソッドの一部を次に示します。
public static Complex Sqrt(double x)
{
return x >= 0 ? new Complex(Math.Sqrt(x)) : new Complex(0, Math.Sqrt(-x));
}
public static Complex Exp(Complex z)
{
return new Complex(Math.Exp(z.Re) * Math.Cos(z.Im), Math.Exp(z.Re) * Math.Sin(z.Im));
}
public static Complex Log(Complex z)
{
return new Complex(Math.Log(Abs(z)), Arg(z));
}
これは複雑な構造体の一部です。
public struct Complex
{
private readonly double _re;
private readonly double _im;
#region Properties
public double Re
{
get { return _re; }
//set { re = value; }
}
public double Im
{
get { return _im; }
//set { im = value; }
}
public static Complex I
{
get { return new Complex(0.0, 1.0); }
}
public static Complex Zero
{
get { return new Complex(0.0, 0.0); }
}
public static Complex One
{
get { return new Complex(1.0, 0.0); }
}
#endregion
#region constructors
public Complex(double x)
{
_re = x;
_im = 0.0;
}
public Complex(Complex z)
{
_re = z.Re;
_im = z.Im;
}
public Complex(double x, double y) //constructor
{
_re = x;
_im = y;
}
#endregion
}