1 
  2 /** 
  3 
  4 	CONVENTIONS:
  5 	In the documentation I refer to Integer and Doubles but this is more to indicate intent
  6 	because all numbers in Javascript, also referred as "number", are 64 bit floating point variables.
  7 **/
  8 "use strict";
  9 /**
 10 	Provides a mathematical matrix object
 11 	@class
 12 	@constructor
 13 	@param rows {Integer} number of rows in the matrix
 14 	@param columns {Integer} number of columns in the matrix
 15 	@param values {Double[][]} a multidimensional array with rows and column as specified by the other parameters
 16 **/
 17 	
 18 function Matrix(rows,columns,values) 
 19 {
 20 	this.rows = rows;
 21 	this.columns = columns;
 22 	this.values = values;
 23 
 24 }
 25 
 26 /**
 27 	Given an array of numbers creates a row Matrix object.
 28 	Essentially shortand for creating matrices of form [a,b,c,...]
 29 	@param {Double[]} the matrix entries 
 30 	@return {Matrix} a 1-by-X matrix
 31 **/
 32 Matrix.rowVector = function(data) {
 33 	return new Matrix(1,data.length, data);
 34 }
 35 
 36 /**
 37 	Given an array of numbers creates a column vector Matrix object.
 38 	Essentially shortand for creating matrices of form 
 39 	[[a]
 40 	 [b]
 41 	 [...]
 42 	 [c]]
 43 	@param {Double[]} the matrix entries 
 44 	@return {Matrix} a X-by-1 matrix
 45 **/
 46 Matrix.columnVector = function(data) {
 47 
 48 	var columns = [];
 49 	for(var i  =0; i < data.length; i++) {
 50 		columns[i] = [data[i]];
 51 	}
 52 	return new Matrix(data.length,1,columns);
 53 }
 54 
 55 /**
 56 	Adds one matrix to another.
 57 	@param that {Matrix} Matrix to be added with.
 58 	@return {Matrix} The resulting sum
 59 **/
 60 
 61 Matrix.prototype.add = function(that)
 62 {
 63 	if(this.rows == that.rows && this.columns == that.columns)
 64 	{
 65 		var newvalues = new Array(this.rows);
 66 		for(var i =0; i < this.rows; i++) {
 67 			newvalues[i] = [];
 68 			for(var j = 0; j < this.columns; j++) {
 69 				newvalues[i][j] = this.values[i][j]+that.values[i][j];
 70 			}
 71 			
 72 		}
 73 		return new Matrix(this.rows, this.columns, newvalues);
 74 	}else{
 75 		throw new Error("Matrix size mismatch; Addition is not defined.");
 76 	}
 77 }
 78 
 79 /**
 80 	Subtract one matrix from another
 81 	@param {Matrix} that The term being subtracted
 82 	@return {Matrix} The difference matrix
 83 **/
 84 Matrix.prototype.subtract = function(that)
 85 {
 86 	if(this.rows == that.rows && this.columns == that.columns)
 87 	{
 88 		var newvalues = new Array(this.rows);
 89 		for(var i =0; i < this.rows; i++) {
 90 			newvalues[i] = [];
 91 			for(var j = 0; j < this.columns; j++) {
 92 				newvalues[i][j] = this.values[i][j]-that.values[i][j];
 93 			}
 94 			
 95 		}
 96 		return new Matrix(this.rows, this.columns, newvalues);
 97 	}else{
 98 		throw new Error("Matrix size mismatch; Addition is not defined.");
 99 	}
100 }
101 
102 /**
103 	Multiply one matrix by another
104 	@param {Matrix} that The matrix being multiplied 
105 	@return {Matrix} The product matrix
106 **/
107 Matrix.prototype.multiply = function(that) 
108 {
109 	if(this.columns == that.rows)
110 	{
111 		var product = new Array(this.rows);
112 		for(var i=0; i < this.rows; i++) {
113 			product[i] = [];
114 			for(var j = 0; j < that.columns; j++) {
115 				product[i][j] =0;
116 			}
117 		}
118 		for(var k=0; k < this.columns; k++) {
119 			for(var i = 0; i < this.rows; i++) {
120 				for(var j =0; j < that.columns; j++) {
121 					product[i][j] += this.values[i][k]*that.values[k][j];
122 				}
123 			}
124 		}
125 		
126 		return new Matrix(this.rows, that.columns, product);
127 		
128 	}else{
129 		throw new Error("Matix size mismatch; Multiplication is not allowed.");
130 	}
131 }
132 
133 /**
134 	Runs a naive Gaussian elimination on a Matrix or a Matrix and a sum vector.
135 	Given the form Ax = b; A is the matrix and b is the sum vector; x is the result
136 	@param sumvector {Matrix (Nx1 sized)} An N row, 1 column vector containing the sums. 
137 	@return {Matrix (Nx1 sized)} The answer
138 	
139 	
140 **/
141 
142 Matrix.prototype.naiveGaussian = function(sumcolumnvector) {
143 	if(! (sumcolumnvector instanceof Matrix))
144 	{
145 		throw new Error("Parameter expects a Matrix");
146 	}
147 
148 	if(this.rows != this.columns) {
149 		throw new Error("Matrix incorrectly sized; Must be NxN matrix");
150 	}
151 	
152 	if(this.rows != sumcolumnvector.rows)
153 	{
154 		throw new Error("Sum vector row count must match the row count of the coefficient matrix");
155 	}
156 	
157 	var newvalues = new Array(this.rows);
158 	for(var i =0; i < this.rows; i++)
159 	{
160 		newvalues[i] = this.values[i].slice();
161 	}
162 	var newsums = sumcolumnvector.transpose().values.slice(); //Probably need to address how I want to reference things for this matrix...
163 	//console.log(newvalues+"");
164 	for(var k = 0; k < this.rows-1; k++) {
165 		for(var i = k+1; i < this.rows; i++) {
166 			var xmult = newvalues[i][k]/newvalues[k][k];
167 			//newvalues[i][k] = xmult;
168 			for(var j = k+1; j < this.rows; j++) {
169 				newvalues[i][j] = newvalues[i][j] - (xmult*newvalues[k][j]);
170 			}
171 			newsums[i] = newsums[i] - xmult*newsums[k];
172 			//if(i == 3) {console.log(newsums[i]+""); }
173 		}
174 		
175 	}
176 	//console.log(newvalues+"");
177 	//console.log(newsums+"");
178 	var solutions = new Array(this.rows);
179 	for(var i = this.rows-1; i >= 0; i--) {
180 		var sum = newsums[i];
181 		for(var j = i+1; j < this.rows; j++) {
182 			sum = sum - newvalues[i][j]*solutions[j];
183 		}
184 		solutions[i] = [sum/newvalues[i][i]];
185 	}
186 	return new Matrix(this.rows, 1, solutions);
187 	
188 	
189 }
190 
191 /**
192 	Runs a Gaussian elimination with scaled partial pivoting on a a Matrix and a sum vector.
193 	Given the form Ax = b; A is the matrix and b is the sum vector; x is the result
194 	@param sumvector {Matrix (Nx1 sized)} An N row, 1 column vector containing the sums. 
195 	@return {Matrix (Nx1 sized)} The answer
196 	
197 	P. 267
198 **/
199 
200 Matrix.prototype.scaledPartialPivotGaussian = function(sumvector) {
201 	if(! (sumvector instanceof Matrix))
202 	{
203 		throw new Error("Parameter expects a Matrix");
204 	}
205 
206 	if(this.rows != this.columns) {
207 		throw new Error("Matrix incorrectly sized; Must be NxN matrix");
208 	}
209 	
210 	if(this.rows != sumvector.rows)
211 	{
212 		throw new Error("Sum vector row count must match the row count of the coefficient matrix");
213 	}
214 	var newvalues = new Array(this.rows);
215 	for(var i =0; i < this.rows; i++)
216 	{
217 		newvalues[i] = this.values[i].slice();
218 	}
219 	var newsums = sumvector.values.slice();
220 	//console.log(newvalues+"");
221 	//Setup for scale and index vectors
222 	var index_vector = [];
223 	var scale_vector = [];
224 	for(var i = 0; i < this.rows; i++)
225 	{
226 		index_vector[i] = i;
227 		scale_vector[i] = 0;
228 		for(var j=0; j < this.rows; j++)
229 		{
230 			scale_vector[i] = Math.max(scale_vector[i], Math.abs(this.values[i][j]));
231 		}
232 	}
233 	//Begin Forward Elimination
234 	for(var k = 0; k < this.rows-1; k++) { 
235 		var max_ratio = 0;
236 		var max_index = 0;
237 		for(var i = k; i< this.rows; i++) {
238 			var ratio = Math.abs(newvalues[index_vector[i]][k]/scale_vector[index_vector[i]]);
239 			if(ratio > max_ratio) {
240 				max_index = i;
241 				max_ratio = ratio;
242 			}
243 		}
244 		var temp = index_vector[k];
245 		index_vector[k] = index_vector[max_index];
246 		index_vector[max_index] = temp;
247 		
248 		for(var i = k+1; i < this.rows; i++) {
249 			var xmult = newvalues[index_vector[i]][k]/newvalues[index_vector[k]][k];
250 			for(var j =k;j < this.columns; j++) {
251 				newvalues[index_vector[i]][j] = newvalues[index_vector[i]][j] - xmult*newvalues[index_vector[k]][j];
252 			}
253 			newsums[index_vector[i]] = newsums[index_vector[i]] - xmult*newsums[index_vector[k]];
254 		}
255 		
256 	}
257 	/*var test = new Matrix(this.rows,this.rows, newvalues);
258 	console.log("values: " +newvalues+"");
259 	console.log("sums: "+newsums+"");
260 	console.log("partial soln: "+test.naiveGaussian(new Matrix(this.rows,1,newsums)).values+"");*/
261 	
262 	var solutions = new Array(this.rows);
263 	for(var i = this.rows-1; i >= 0; i--) {
264 		var sum = newsums[index_vector[i]];
265 		for(var j = i+1; j < this.columns; j++) {
266 			sum = sum - newvalues[index_vector[i]][j]*solutions[j];
267 		}
268 		solutions[i] = [sum/newvalues[index_vector[i]][i]];
269 	}
270 	return new Matrix(this.rows, 1, solutions);
271 	
272 }
273 
274 /**
275 	An operation on graph adjacency matrix. 
276 	It checks if a graph contains a sink node. A node which
277 	has V-1 inbound edges and 0 outbound edges. 
278 	Running time is only O(n) 
279 	@return {boolean} 
280 **/
281 
282 Matrix.prototype.findSink =function() { 
283 	var k =0;
284 	for(var i=0; i < this.columns; i++) {
285 		if( this.values[k][i] == 1) {
286 			k = i;
287 		}
288 	}
289 	
290 	for(var i=0; i < this.columns; i++) {
291 		if(this.values[k][i] == 1) {
292 			return false;
293 		}
294 		if(this.values[i][k] == 0 && i != k) {
295 			return false;
296 		}
297 	}
298 	return true;
299 }
300 
301 /**
302 	Produces two matrices, L and U, which are lower triangular
303 	matrice and upper triangular matrices respectively.
304 	@return {Matrix[]} An array of Matrix objects, L = 
305 **/
306 Matrix.prototype.LUdecomposition = function() {
307 	if(this.rows != this.columns) {
308 		throw new Error("Matrix incorrectly sized; Must be NxN matrix");
309 	}
310 	
311 	var L = new Matrix(this.rows, this.columns, []);
312 	var U = new Matrix(this.rows, this.columns, []);
313 	var n = this.rows;
314 	var zeroarray = [];
315 	for(var k=0; k < n; k++) {
316 		zeroarray[k] = 0;
317 	}
318 	
319 	for(var k=0; k < n; k++) {
320 		L.values[k]= zeroarray.slice();
321 		U.values[k]= zeroarray.slice();
322 	}
323 	
324 	for(var k=0; k < n; k++) {
325 		L.values[k][k]=1;
326 		for(var j = k; j < n; j++) {
327 		
328 			var sum = 0;
329 			for(var s = 0; s < k; s++) {
330 				sum += L.values[k][s]*U.values[s][j];
331 			}
332 			U.values[k][j] = this.values[k][j] - sum;
333 		}
334 		
335 		for(var i = k+1; i < n; i++) {
336 			var sum =0;
337 			for(var s = 0; s < k; s++) {
338 				sum += L.values[i][s]*U.values[s][k];
339 			}
340 		
341 			L.values[i][k] = (this.values[i][k] - sum)/U.values[k][k];
342 		}
343 	}
344 	
345 	return [L,U];
346 }
347 
348 /**
349 	Solves a symmetric matric of the tri-diagonal form
350 	@author Sahil Diwan
351 	@param {Matrix} columnvector An N-by-1 matrix
352 	@return {Matrix} Solutions in an N-by-1 matrix
353 	@todo implement algorithm on page 282
354 **/
355 Matrix.prototype.tridiagonalSolve = function(columnvector) {
356 	
357 
358 }
359 
360 /**
361 	simple implementation of Floyd-Warshall All-Pairs Shortest Path algorithm for weighted graphs
362 	@return {Matrix} The APSP matrix for a given matrix
363 	@note Currently this will expect matrices to be of a prepared form
364 	var fwtest = [[0,Infinity,Infinity,Infinity,-1,Infinity],
365 		      [1,0,Infinity,2,Infinity,Infinity],
366 		      [Infinity,2,0,Infinity,Infinity,-8],
367 		      [-4,Infinity,Infinity,0,3,Infinity],
368 		      [Infinity,7,Infinity,Infinity,0,Infinity],
369 		      [Infinity,5,10,Infinity,Infinity,0]];
370       var mm = new Matrix(6,6,fwtest);
371 **/
372 Matrix.prototype.FloydWarshall = function() {
373 	if(this.rows != this.columns) {
374 		throw new Error("Matrix incorrectly sized; Must be NxN matrix");
375 	}
376 	var D = this.copy(); 
377 	for(var k=0; k < this.rows; k++) {
378 		for(var i =0; i < this.rows; i++) {
379 			for(var j = 0; j < this.rows; j++) {
380 				if(D.values[i][k]+D.values[k][j] < D.values[i][j]) {
381 					D.values[i][j] =D.values[i][k]+D.values[k][j]; 
382 				}
383 			}
384 		}
385 		$("body").append("D(k)=<br>"+D);
386 	}
387 	return D;
388 }
389 
390 /**
391 	Provide a copy contructor
392 	@return {Matrix} A copy of the existing matrix
393 **/
394 Matrix.prototype.copy = function() {
395 	var retmatrix = new Matrix(this.rows,this.columns,[]);
396 	for(var i =0; i < this.rows; i++) {
397 		retmatrix.values[i]=[];
398 		for(var j =0; j < this.columns; j++) {
399 			retmatrix.values[i][j] = this.values[i][j];
400 		}
401 	}
402 	return retmatrix;
403 }
404 
405 /**
406 	Provides the Transpose matrix operation, i.e. flips the matrix
407 	@return {Matrix} The transpose of the original matrix
408 **/
409 Matrix.prototype.transpose = function() {
410 	var retmatrix = new Matrix(this.columns, this.rows,[]);
411 	for(var i =0; i < this.columns; i++) {
412 		retmatrix.values[i]=[];
413 		for(var j =0; j < this.rows; j++) {
414 			retmatrix.values[i][j] = this.values[j][i];
415 		}
416 	}
417 	return retmatrix;
418 }
419 //I might rename this toHTML and create a different to string method
420 Matrix.prototype.toHTML = function() {
421 	var retstring = "";
422 	for(var i = 0; i < this.rows; i++) {
423 		retstring += "| ";
424 		for(var k = 0; k < this.columns; k++) {
425 			retstring += this.values[i][k] +" ";
426 		}
427 		retstring +="|<br />";
428 	}
429 	return retstring;
430 }
431 
432 Matrix.prototype.toString = function() {
433 	var retstring = "[ ";
434 	for(var i = 0; i < this.rows; i++) {
435 		retstring += "[ ";
436 		for(var k = 0; k < this.columns; k++) {
437 			retstring += this.values[i][k] +" ";
438 		}
439 		retstring +="] ";
440 	}
441 	retstring +="]";
442 	return retstring;
443 }
444 /*
445 var testdata = [[6,-2,2,4],
446 		[12,-8,6,10],
447 		[3,-13,9,3],
448 		[-6,4,1,-18]];
449 		
450 var testsumvector = [16,26,-19,-34];
451 var mymatrix = new Matrix(4,4,testdata);
452 var sumvec = new Matrix(4,1,testsumvector);
453 
454 var testdata2 = [[3,-13,9,3],
455 		 [-6,4,1,-18],
456 		 [6,-2,2,4],
457 		 [12,-8,6,10]];
458 var sumvec2 = new Matrix(4,1,[-19,-34,16,26]);
459 var sumvec22 = new Matrix(4,1,[[-19],[-34],[16],[26]]);
460 var ppmatrix = new Matrix(4,4,testdata2);
461 
462 var left2 = new Matrix(2,2,[[1,3],[4,4]]); //answer is [[7,17],[20,28]]
463 var right2 = new Matrix(2,2,[[4,2],[1,5]]);
464 
465 var left1 =new Matrix(1,2,[[1,2]]); //anser is [[10]]
466 var right1 =new Matrix(2,1,[[4],[3]]);
467 var problem2 = new Matrix(3,3,[[1,-1,2],[-2,1,-1],[4,-1,2]]);
468 problem2.ScaledPartialPivotGaussian(new Matrix(3,1,[-2,2,1]));
469 var problem3 = new Matrix(3,3,[[2,4,-2],[1,3,4],[5,2,0]]);
470 problem3.ScaledPartialPivotGaussian(new Matrix(3,1,[6,-1,2])).values+"";
471 
472 var t = [new Matrix(2,2,[[0,0],[1,0]])];
473 t[1] = new Matrix(2,2,[[0,1],[1,0]]);
474 t[2] = new Matrix(3,3,[[0,1,0],[0,0,0],[0,1,0]]);
475 t[3] = new Matrix(3,3,[[0,1,0],[0,0,0],[0,1,1]]);
476 t[4] = new Matrix(3,3,[[0,1,0],[0,0,0],[0,0,1]]);
477 t[5] = new Matrix(3,3,[[0,0,0],[0,0,0],[0,1,1]]);
478 t[6] = new Matrix(4,4,[[0,0,1,0],[0,0,1,0],[0,0,0,0],[0,0,1,0]]);
479 t[7] = new Matrix(4,4,[[0,1,1,1],[1,0,1,1],[0,0,0,0],[1,1,1,0]]);
480 t[8] = new Matrix(4,4,[[0,1,1,1],[1,0,1,1],[0,0,0,1],[1,1,1,0]]);
481 t[9] = new Matrix(4,4,[[0,1,1,1],[1,0,1,1],[0,0,0,0],[1,1,0,0]]);
482 */
483