我将如何在 JavaScript 中实现 MATLAB 的 "eig(A, B)" 函数
How would I implement MATLAB's "eig(A, B)" function in JavaScript
根据 MATLAB 的 documentation:
[V
,D
] = eig(A
,B
) returns diagonal matrix D
of generalized eigenvalues and full matrix V
whose columns are the corresponding right eigenvectors, so that A*V = B*V*D
.
当我阅读可用的源代码时(似乎我看过 Octave、R、Scipy 的所有实现)导致 LAPACK 的 DGGEV
例程,[=31] 中没有=].
兔子洞包括 Use Emscripten with Fortran: LAPACK binding 和学习足够的 Fortran 和线性代数来自己做。
有人知道更方便的东西吗?
我最终使用 Eigen 库和 Emscripten 取得了成功。
现在,我的测试代码被硬编码为 5x5 矩阵,但这只是模板参数的问题。
我使用行主一维数组将数据传入和传出函数。
代码如下所示:
#include <Eigen/Eigenvalues>
#typedef double ArrayMat5d[25];
#typedef double ArrayVec5d[5];
#typedef Eigen::Matrix<double, 5, 5, Eigen::RowMajor> Matrix5dR;
#typedef Eigen::Matrix<double, 5, 1> Vector5d;
extern "C" void eig(const ArrayMat5d A, const ArrayMat5d B,
ArrayMat5d V, ArrayVec5d D) {
Eigen::Map<const Matrix5dR> a(A);
Eigen::Map<const Matrix5dR> b(B);
const Eigen::GeneralizedSelfAdjointEigenSolver<Matrix5dR> solver(a, b);
Eigen::Map<Matrix5dR> v(V);
Eigen::Map<Vector5d> d(D);
v = solver.eigenvectors();
d = solver.eigenvalues();
}
我正在使用以下代码编译代码:
emcc -I /usr/include/eigen3 -O2 -o eig.js -s "DISABLE_EXCEPTION_CATCHING = 1" \
-s "NO_FILESYSTEM = 1" -s "NO_BROWSER = 1" -s "EXPORTED_FUNCTIONS = ['_eig']" \
-s "NO_EXIT_RUNTIME = 1" eig.cpp
JavaScript 方面:
// builds reference to eig function with argument type checking
var eig = Module.cwrap('eig', null, ['number', 'number', 'number', 'number']);
// sets up the two matrices
var P = new Float64Array([ 92.31360, 11.75040, -15.84640, -21.88800, -0.83200, 11.75040, 15.76960, -4.37760, -0.83200, 2.11200, -15.84640, -4.37760, 4.24960, 2.11200, -1.15200, -21.88800, -0.83200, 2.11200, 15.04000, -1.44000, -0.83200, 2.11200, -1.15200, -1.44000, -2.24000 ]);
var Q = new Float64Array([ 60.16, -2.88, 0.0, 0.0, 0.0, -2.88, 17.28, -2.88, 0.0, 0.0, 0.0, -2.88, 8.96, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0 ]);
// allocates memory for input and output matrices
var matLength = 25;
var vecLength = 5;
var matSize = matLength * P.BYTES_PER_ELEMENT;
var vecSize = vecLength * P.BYTES_PER_ELEMENT;
var Pptr = Module._malloc(matSize);
var Qptr = Module._malloc(matSize);
var Vptr = Module._malloc(matSize);
var Dptr = Module._malloc(vecSize);
// gets references to Emscripten heap
var Pheap = new Uint8Array(Module.HEAPU8.buffer, Pptr, matSize);
var Qheap = new Uint8Array(Module.HEAPU8.buffer, Qptr, matSize);
var Vheap = new Uint8Array(Module.HEAPU8.buffer, Vptr, matSize);
var Dheap = new Uint8Array(Module.HEAPU8.buffer, Dptr, vecSize);
// copies input matrices into Emscripten heap
Pheap.set(new Uint8Array(P.buffer));
Qheap.set(new Uint8Array(Q.buffer));
// calls the function (finally!)
eig(Pheap.byteOffset, Qheap.byteOffset, Vheap.byteOffset, Dheap.byteOffset);
// Gets double views into Emscripten heap containing results
var Vresult = new Float64Array(Vheap.buffer, Vheap.byteOffset, P.length);
var Dresult = new Float64Array(Dheap.buffer, Dheap.byteOffset, vecLength);
console.log(Vresult);
console.log(Dresult);
// Frees up allocated memory
Module._free(Pheap.byteOffset);
Module._free(Qheap.byteOffset);
Module._free(Vheap.byteOffset);
Module._free(Dheap.byteOffset);
整个过程运行良好。在级别 -O2
,我得到大约 800 毫秒到 运行 10000 次迭代的时间,结果与我的原始 C++ 测试代码完全匹配。 (它在 -O0
时大约慢了 10 倍。)
现在完成椭圆拟合!
我制作了一个 javascript 模块来计算特征值、特征向量或 RREF(简化行梯形)。您可以使用 npm i @ahmaddynugroho/eig
安装它
示例:
const e = require('@ahmaddynugroho/eig')
const A = [[7,3],[3,-1]]
const ans = e.eig(A)
console.log(ans)
这将记录 ans
对象如下:
{
eigval: [ 8.000000000000009, -2.0000000000000004 ],
eigvec: [
[ -0.9486832980505129, 0.31622776601684044 ],
[ 0.3162277660168379, 0.9486832980505138 ]
]
}
就是这样,你只需要传递一个矩阵作为eig()
的参数。顺便说一句,这里是 the documentation
根据 MATLAB 的 documentation:
[
V
,D
] = eig(A
,B
) returns diagonal matrixD
of generalized eigenvalues and full matrixV
whose columns are the corresponding right eigenvectors, so thatA*V = B*V*D
.
当我阅读可用的源代码时(似乎我看过 Octave、R、Scipy 的所有实现)导致 LAPACK 的 DGGEV
例程,[=31] 中没有=].
兔子洞包括 Use Emscripten with Fortran: LAPACK binding 和学习足够的 Fortran 和线性代数来自己做。
有人知道更方便的东西吗?
我最终使用 Eigen 库和 Emscripten 取得了成功。
现在,我的测试代码被硬编码为 5x5 矩阵,但这只是模板参数的问题。
我使用行主一维数组将数据传入和传出函数。
代码如下所示:
#include <Eigen/Eigenvalues>
#typedef double ArrayMat5d[25];
#typedef double ArrayVec5d[5];
#typedef Eigen::Matrix<double, 5, 5, Eigen::RowMajor> Matrix5dR;
#typedef Eigen::Matrix<double, 5, 1> Vector5d;
extern "C" void eig(const ArrayMat5d A, const ArrayMat5d B,
ArrayMat5d V, ArrayVec5d D) {
Eigen::Map<const Matrix5dR> a(A);
Eigen::Map<const Matrix5dR> b(B);
const Eigen::GeneralizedSelfAdjointEigenSolver<Matrix5dR> solver(a, b);
Eigen::Map<Matrix5dR> v(V);
Eigen::Map<Vector5d> d(D);
v = solver.eigenvectors();
d = solver.eigenvalues();
}
我正在使用以下代码编译代码:
emcc -I /usr/include/eigen3 -O2 -o eig.js -s "DISABLE_EXCEPTION_CATCHING = 1" \
-s "NO_FILESYSTEM = 1" -s "NO_BROWSER = 1" -s "EXPORTED_FUNCTIONS = ['_eig']" \
-s "NO_EXIT_RUNTIME = 1" eig.cpp
JavaScript 方面:
// builds reference to eig function with argument type checking
var eig = Module.cwrap('eig', null, ['number', 'number', 'number', 'number']);
// sets up the two matrices
var P = new Float64Array([ 92.31360, 11.75040, -15.84640, -21.88800, -0.83200, 11.75040, 15.76960, -4.37760, -0.83200, 2.11200, -15.84640, -4.37760, 4.24960, 2.11200, -1.15200, -21.88800, -0.83200, 2.11200, 15.04000, -1.44000, -0.83200, 2.11200, -1.15200, -1.44000, -2.24000 ]);
var Q = new Float64Array([ 60.16, -2.88, 0.0, 0.0, 0.0, -2.88, 17.28, -2.88, 0.0, 0.0, 0.0, -2.88, 8.96, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0 ]);
// allocates memory for input and output matrices
var matLength = 25;
var vecLength = 5;
var matSize = matLength * P.BYTES_PER_ELEMENT;
var vecSize = vecLength * P.BYTES_PER_ELEMENT;
var Pptr = Module._malloc(matSize);
var Qptr = Module._malloc(matSize);
var Vptr = Module._malloc(matSize);
var Dptr = Module._malloc(vecSize);
// gets references to Emscripten heap
var Pheap = new Uint8Array(Module.HEAPU8.buffer, Pptr, matSize);
var Qheap = new Uint8Array(Module.HEAPU8.buffer, Qptr, matSize);
var Vheap = new Uint8Array(Module.HEAPU8.buffer, Vptr, matSize);
var Dheap = new Uint8Array(Module.HEAPU8.buffer, Dptr, vecSize);
// copies input matrices into Emscripten heap
Pheap.set(new Uint8Array(P.buffer));
Qheap.set(new Uint8Array(Q.buffer));
// calls the function (finally!)
eig(Pheap.byteOffset, Qheap.byteOffset, Vheap.byteOffset, Dheap.byteOffset);
// Gets double views into Emscripten heap containing results
var Vresult = new Float64Array(Vheap.buffer, Vheap.byteOffset, P.length);
var Dresult = new Float64Array(Dheap.buffer, Dheap.byteOffset, vecLength);
console.log(Vresult);
console.log(Dresult);
// Frees up allocated memory
Module._free(Pheap.byteOffset);
Module._free(Qheap.byteOffset);
Module._free(Vheap.byteOffset);
Module._free(Dheap.byteOffset);
整个过程运行良好。在级别 -O2
,我得到大约 800 毫秒到 运行 10000 次迭代的时间,结果与我的原始 C++ 测试代码完全匹配。 (它在 -O0
时大约慢了 10 倍。)
现在完成椭圆拟合!
我制作了一个 javascript 模块来计算特征值、特征向量或 RREF(简化行梯形)。您可以使用 npm i @ahmaddynugroho/eig
示例:
const e = require('@ahmaddynugroho/eig')
const A = [[7,3],[3,-1]]
const ans = e.eig(A)
console.log(ans)
这将记录 ans
对象如下:
{
eigval: [ 8.000000000000009, -2.0000000000000004 ],
eigvec: [
[ -0.9486832980505129, 0.31622776601684044 ],
[ 0.3162277660168379, 0.9486832980505138 ]
]
}
就是这样,你只需要传递一个矩阵作为eig()
的参数。顺便说一句,这里是 the documentation