1 文本格式
using System;
namespace Legalsoft.Truffer
{
public class Multinormaldev : Ran
{
public Cholesky chol { get; set; } = null;
private int mm { get; set; }
private double[] mean { get; set; }
private double[,] xvar { get; set; }
private double[] spt { get; set; }
private double[] pt { get; set; }
public Multinormaldev(ulong j, double[] mmean, double[,] vvar) : base(j)
{
this.mm = mmean.Length;
this.mean = Globals.CopyFrom(mmean);
this.xvar = Globals.CopyFrom(vvar);
this.chol = new Cholesky(xvar);
this.spt = new double[mm];
this.pt = new double[mm];
if (xvar.GetLength(1) != mm || xvar.GetLength(0) != mm)
{
throw new Exception("bad sizes");
}
}
public double[] dev()
{
double u;
double v;
double q;
for (int i = 0; i < mm; i++)
{
do
{
u = doub();
v = 1.7156 * (doub() - 0.5);
double x = u - 0.449871;
double y = Math.Abs(v) + 0.386595;
q = Globals.SQR(x) + y * (0.19600 * y - 0.25472 * x);
} while (q > 0.27597 && (q > 0.27846 || Globals.SQR(v) > -4.0 * Math.Log(u) * Globals.SQR(u)));
spt[i] = v / u;
}
chol.elmult(spt, pt);
for (int i = 0; i < mm; i++)
{
pt[i] += mean[i];
}
return pt;
}
}
}
2 代码格式
using System;
namespace Legalsoft.Truffer
{
public class Multinormaldev : Ran
{
public Cholesky chol { get; set; } = null;
private int mm { get; set; }
private double[] mean { get; set; }
private double[,] xvar { get; set; }
private double[] spt { get; set; }
private double[] pt { get; set; }
public Multinormaldev(ulong j, double[] mmean, double[,] vvar) : base(j)
{
this.mm = mmean.Length;
this.mean = Globals.CopyFrom(mmean);
this.xvar = Globals.CopyFrom(vvar);
this.chol = new Cholesky(xvar);
this.spt = new double[mm];
this.pt = new double[mm];
if (xvar.GetLength(1) != mm || xvar.GetLength(0) != mm)
{
throw new Exception("bad sizes");
}
}
public double[] dev()
{
double u;
double v;
double q;
for (int i = 0; i < mm; i++)
{
do
{
u = doub();
v = 1.7156 * (doub() - 0.5);
double x = u - 0.449871;
double y = Math.Abs(v) + 0.386595;
q = Globals.SQR(x) + y * (0.19600 * y - 0.25472 * x);
} while (q > 0.27597 && (q > 0.27846 || Globals.SQR(v) > -4.0 * Math.Log(u) * Globals.SQR(u)));
spt[i] = v / u;
}
chol.elmult(spt, pt);
for (int i = 0; i < mm; i++)
{
pt[i] += mean[i];
}
return pt;
}
}
}