-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathindex.js
75 lines (67 loc) · 2.84 KB
/
index.js
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
var emlapack = require('emlapack');
exports.eig = (matrix) => {
//See: http://www.netlib.org/lapack/explore-html/d3/dfb/group__real_g_eeigen_ga104525b749278774f7b7f57195aa6798.html#ga104525b749278774f7b7f57195aa6798
var inputFormatErrorMessage = 'Input needs to be in a matrix form, for example: [ [ 3 ] ] or [ [3, 1], [1, 3] ]';
if (!Array.isArray(matrix)) {
throw new Error(inputFormatErrorMessage);
}
else {
matrix.forEach(row => {
if (!Array.isArray(row)) {
throw new Error(inputFormatErrorMessage);
}
});
}
var n = matrix.length,
dgeev = emlapack.cwrap('dgeev_', null, ['number', 'number', 'number', 'number', 'number', 'number', 'number', 'number', 'number', 'number', 'number', 'number', 'number', 'number']),
pjobvl = emlapack._malloc(1),
pjobvr = emlapack._malloc(1),
pn = emlapack._malloc(4),
pa = emlapack._malloc(n * n * 8),
plda = emlapack._malloc(4),
pwr = emlapack._malloc(n * 8),
pwi = emlapack._malloc(n * 8),
pvl = emlapack._malloc(n * n * 8),
pldvl = emlapack._malloc(4),
pvr = emlapack._malloc(n * n * 8),
pldvr = emlapack._malloc(4),
pworkopt = emlapack._malloc(n * n * 8),
plwork = emlapack._malloc(4),
pinfo = emlapack._malloc(4);
//Calculate both left and right eigenvectors
emlapack.setValue(pjobvl, 'V'.charCodeAt(0), 'i8');
emlapack.setValue(pjobvr, 'V'.charCodeAt(0), 'i8');
emlapack.setValue(pn, n, 'i32');
emlapack.setValue(plda, n, 'i32');
emlapack.setValue(pldvl, n, 'i32');
emlapack.setValue(pldvr, n, 'i32');
emlapack.setValue(plwork, -1, 'i32');
var a = new Float64Array(emlapack.HEAPF64.buffer, pa, n * n);
var w = new Float64Array(emlapack.HEAPF64.buffer, pwr, n);
var i = new Float64Array(emlapack.HEAPF64.buffer, pwi, n);
var rightEigenvectors = new Float64Array(emlapack.HEAPF64.buffer, pvr, n * n);
var leftEigenvectors = new Float64Array(emlapack.HEAPF64.buffer, pvl, n * n);
var flattenedMatrix = [];
matrix.forEach(row => {
row.forEach(element => {
flattenedMatrix.push(element);
});
});
a.set(flattenedMatrix);
//Somehow pvl and pvr are switched in the Lapack implementation (bug?)
dgeev(pjobvl, pjobvr, pn, pa, plda, pwr, pwi, pvr, pldvl, pvl, pldvr, pworkopt, plwork, pinfo);
var workopt = emlapack.getValue(pworkopt, 'double'),
pwork = emlapack._malloc(workopt * 8);
emlapack.setValue(plwork, workopt, 'i32');
dgeev(pjobvl, pjobvr, pn, pa, plda, pwr, pwi, pvr, pldvl, pvl, pldvr, pwork, plwork, pinfo);
return {
eigenvalues: {
real: w,
imaginary: i
},
eigenvectors: {
right: rightEigenvectors,
left: leftEigenvectors
}
};
};