Changeset 10

trunk/dsss.conf

r9	r10
7	7	[mathematics/imaging/algorithms]
8	8	[mathematics/examples/main.d]
	9	noinstall

trunk/mathematics/numerical/Matrix.d

r9	r10
92	92	***********************************************************************/
93	93
94		~~uint getR~~owCount()
	94	final uint rowCount()
95	95	{
96	96	return _rowCount;
…	…
103	103	***********************************************************************/
104	104
105		~~uint getC~~olumnCount()
	105	final uint columnCount()
106	106	{
107	107	return _columnCount;
…	…
118	118	{
119	119	T result = 0;
120		for (uint r = 0; r < ~~getR~~owCount; r++)
121		for (uint c = 0; c < ~~getC~~olumnCount; c++)
	120	for (uint r = 0; r < rowCount; r++)
	121	for (uint c = 0; c < columnCount; c++)
122	122	result += this[r,c] * this[r,c];
123	123	return result;
…	…
133	133	{
134	134	T result = 0;
135		for (uint r = 0; r < ~~getRowCount~~; r++)
136		for (uint c = 0; c < ~~getColumnCount~~; c++)
	135	for (uint r = 0; r < rowCount(); r++)
	136	for (uint c = 0; c < columnCount(); c++)
137	137	result += this[r,c];
138	138	return result;
…	…
148	148	{
149	149	char [] result = "\n";
150		for (uint r = 0; r < ~~getR~~owCount(); r++)
	150	for (uint r = 0; r < rowCount(); r++)
151	151	{
152	152	result ~= "(";
153		for (uint c = 0; c < ~~getC~~olumnCount(); c++)
	153	for (uint c = 0; c < columnCount(); c++)
154	154	version (Tango)
155	155	{
…	…
172	172	bool opEquals(MatrixBase!(T) right)
173	173	{
174		bool result = this.~~getRowCount() == right.getRowCount() && this.getColumnCount() == right.getC~~olumnCount();
	174	bool result = this.rowCount() == right.rowCount() && this.columnCount() == right.columnCount();
175	175	if (!result) return result;
176		for (uint r = 0; r < ~~getR~~owCount(); r++)
177		for (uint c = 0; c < ~~getC~~olumnCount(); c++)
	176	for (uint r = 0; r < rowCount(); r++)
	177	for (uint c = 0; c < columnCount(); c++)
178	178	if (this[r,c] != right[r,c]) return false;
179	179	return result;
…	…
191	191	Matrix!(T) newTransposeMatrix()
192	192	{
193		Matrix!(T) result = new Matrix!(T)(~~getColumnCount(), getR~~owCount());
194		for (uint r = 0; r < this.~~getR~~owCount; r++)
195		for (uint c = 0; c < ~~getC~~olumnCount; c++)
	193	Matrix!(T) result = new Matrix!(T)(columnCount(), rowCount());
	194	for (uint r = 0; r < this.rowCount; r++)
	195	for (uint c = 0; c < columnCount; c++)
196	196	result[c,r] = this[r,c];
197	197	return result;
…	…
210	210	{
211	211
212		assert( this.~~getRowCount() == right.getR~~owCount() , " RowCount error");
213		assert( this.~~getColumnCount() == right.getC~~olumnCount() , "ColumnCount doesn't match");
214		Matrix!(T) result = new Matrix!(T)(this.~~getRowCount(), this.getC~~olumnCount);
215		for (uint r = 0; r < ~~getR~~owCount(); r++)
216		for (uint c = 0; c < ~~getC~~olumnCount(); c++)
	212	assert( this.rowCount() == right.rowCount() , " RowCount error");
	213	assert( this.columnCount() == right.columnCount() , "ColumnCount doesn't match");
	214	Matrix!(T) result = new Matrix!(T)(this.rowCount(), this.columnCount);
	215	for (uint r = 0; r < rowCount(); r++)
	216	for (uint c = 0; c < columnCount(); c++)
217	217	result[r,c] = this[r,c] + right[r,c];
218	218	return result;
…	…
235	235	Matrix!(T) opSub(MatrixBase!(T) right)
236	236	{
237		assert( this.~~getRowCount() == right.getR~~owCount() , " RowCount doesn't match");
238		assert( this.~~getColumnCount() == right.getC~~olumnCount() , "ColumnCount doesn't match");
239		Matrix!(T) result = new Matrix!(T)(this.~~getRowCount(), this.getC~~olumnCount);
240		for (uint r = 0; r < ~~getR~~owCount(); r++)
241		for (uint c = 0; c < ~~getC~~olumnCount(); c++)
	237	assert( this.rowCount() == right.rowCount() , " RowCount doesn't match");
	238	assert( this.columnCount() == right.columnCount() , "ColumnCount doesn't match");
	239	Matrix!(T) result = new Matrix!(T)(this.rowCount(), this.columnCount);
	240	for (uint r = 0; r < rowCount(); r++)
	241	for (uint c = 0; c < columnCount(); c++)
242	242	result[r,c] = this[r,c] - right[r,c];
243	243	return result;
…	…
260	260	Matrix!(T) opMul(MatrixBase!(T) right)
261	261	{
262		assert( this.~~getColumnCount() == right.getR~~owCount(), "Incorrect matrix dimensions for multiplication");
263		Matrix!(T) result = new Matrix!(T)(~~getRowCount(), getC~~olumnCount());
	262	assert( this.columnCount() == right.rowCount(), "Incorrect matrix dimensions for multiplication");
	263	Matrix!(T) result = new Matrix!(T)(rowCount(), columnCount());
264	264	T intermediate;
265		for (uint r = 0; r < this.~~getR~~owCount(); r++)
266		for (uint c = 0; c < right.~~getC~~olumnCount; c++)
	265	for (uint r = 0; r < this.rowCount(); r++)
	266	for (uint c = 0; c < right.columnCount; c++)
267	267	{
268	268	intermediate = 0;
269		for (uint i = 0; i < this.~~getC~~olumnCount(); i++)
	269	for (uint i = 0; i < this.columnCount(); i++)
270	270	intermediate += this.opIndex(r,i) * right[i,c];
271	271	result[r,c] = intermediate;
…	…
289	289	Matrix!(T) opMul(real scalar)
290	290	{
291		Matrix!(T) result = new Matrix!(T)(this.~~getRowCount(), this.getC~~olumnCount());
292		for (uint r = 0; r < ~~getR~~owCount(); r++)
293		for (uint c = 0; c < ~~getC~~olumnCount(); c++)
	291	Matrix!(T) result = new Matrix!(T)(this.rowCount(), this.columnCount());
	292	for (uint r = 0; r < rowCount(); r++)
	293	for (uint c = 0; c < columnCount(); c++)
294	294	result[r,c] = this.opIndex(r, c) * scalar;
295	295	return result;
…	…
311	311	Matrix!(T) opMul_r(real scalar)
312	312	{
313		Matrix!(T) result = new Matrix!(T)(this.~~getRowCount(), this.getC~~olumnCount());
314		for (uint r = 0; r < ~~getR~~owCount(); r++)
315		for (uint c = 0; c < ~~getC~~olumnCount(); c++)
	313	Matrix!(T) result = new Matrix!(T)(this.rowCount(), this.columnCount());
	314	for (uint r = 0; r < rowCount(); r++)
	315	for (uint c = 0; c < columnCount(); c++)
316	316	result[r,c] = this.opIndex(r, c) * scalar;
317	317	return result;
…	…
327	327	Matrix!(T) opNeg()
328	328	{
329		Matrix!(T) result = new Matrix!(T) ( ~~getRowCount(), getC~~olumnCount());
330		for (uint r = 0; r < ~~getR~~owCount(); r++)
331		for (uint c = 0; c < ~~getC~~olumnCount; c++)
	329	Matrix!(T) result = new Matrix!(T) ( rowCount(), columnCount());
	330	for (uint r = 0; r < rowCount(); r++)
	331	for (uint c = 0; c < columnCount; c++)
332	332	result[r,c] = -this[r,c];
333	333	return result;
…	…
344	344	T [] getColumn(uint column)
345	345	{
346		assert( column >= 0 && column < ~~getC~~olumnCount(), "column index out of bound");
347		T[] result = new T[~~getR~~owCount()];
348		for (uint r = 0; r < ~~getR~~owCount; r++)
	346	assert( column >= 0 && column < columnCount(), "column index out of bound");
	347	T[] result = new T[rowCount()];
	348	for (uint r = 0; r < rowCount; r++)
349	349	result[r] = this[r, column];
350	350	return result;
…	…
360	360	T [] getRow(uint row)
361	361	{
362		T [] result = new T[~~getC~~olumnCount()];
	362	T [] result = new T[columnCount()];
363	363	for (uint c = 0; c < _columnCount; c++)
364	364	result[c] = this.opIndex(row, c);
…	…
374	374	Matrix!(T) newClone()
375	375	{
376		Matrix!(T) result = new Matrix!(T)(~~getRowCount(), getC~~olumnCount());
377		for (uint r = 0; r < ~~getR~~owCount(); r++)
378		for (uint c = 0; c < ~~getC~~olumnCount(); c++)
	376	Matrix!(T) result = new Matrix!(T)(rowCount(), columnCount());
	377	for (uint r = 0; r < rowCount(); r++)
	378	for (uint c = 0; c < columnCount(); c++)
379	379	result[r,c] = this[r,c];
380	380	return result;
…	…
393	393	{
394	394
395		assert(this.~~getRowCount() == 2 && this.getC~~olumnCount() == 2, "not a 2x2 matrix");
	395	assert(this.rowCount() == 2 && this.columnCount() == 2, "not a 2x2 matrix");
396	396	VectorFixed!(T, cast(uint) 2) result;
397	397	T intermediate;
398		for (uint r = 0; r < this.~~getR~~owCount(); r++)
	398	for (uint r = 0; r < this.rowCount(); r++)
399	399	{
400	400	intermediate = 0;
401		for (uint i = 0; i < this.~~getC~~olumnCount(); i++)
	401	for (uint i = 0; i < this.columnCount(); i++)
402	402	intermediate += this.opIndex(r,i) * right[i];
403	403	result[r] = intermediate;
…	…
414	414	VectorFixed!(T, cast(uint) 3) mul3x3(ref VectorFixed!(T, cast(uint) 3) right)
415	415	{
416		assert(this.~~getRowCount() == 3 && this.getC~~olumnCount() == 3, "not a 3x3 matrix");
	416	assert(this.rowCount() == 3 && this.columnCount() == 3, "not a 3x3 matrix");
417	417	VectorFixed!(T, cast(uint) 3) result;
418	418	T intermediate;
419		for (uint r = 0; r < this.~~getR~~owCount(); r++)
	419	for (uint r = 0; r < this.rowCount(); r++)
420	420	{
421	421	intermediate = 0;
422		for (uint i = 0; i < this.~~getC~~olumnCount(); i++)
	422	for (uint i = 0; i < this.columnCount(); i++)
423	423	intermediate += this.opIndex(r,i) * right[i];
424	424	result[r] = intermediate;
…	…
436	436	VectorFixed!(T, cast(uint) 4) mul4x4(ref VectorFixed!(T, cast(uint) 4) right)
437	437	{
438		assert(this.~~getRowCount() == 4 && this.getC~~olumnCount() == 4, "not a 4x4 matrix");
	438	assert(this.rowCount() == 4 && this.columnCount() == 4, "not a 4x4 matrix");
439	439	VectorFixed!(T, cast(uint) 4) result;
440	440	T intermediate;
441		for (uint r = 0; r < this.~~getR~~owCount(); r++)
	441	for (uint r = 0; r < this.rowCount(); r++)
442	442	{
443	443	intermediate = 0;
444		for (uint i = 0; i < this.~~getC~~olumnCount(); i++)
	444	for (uint i = 0; i < this.columnCount(); i++)
445	445	intermediate += this.opIndex(r,i) * right[i];
446	446	result[r] = intermediate;
…	…
457	457	VectorFixed!(T, cast(uint) 3) mul3x4(ref VectorFixed!(T, cast(uint) 4) right)
458	458	{
459		assert(this.~~getRowCount() == 3 && this.getC~~olumnCount() == 4, "not a 4x4 matrix");
	459	assert(this.rowCount() == 3 && this.columnCount() == 4, "not a 4x4 matrix");
460	460	VectorFixed!(T, cast(uint) 3) result;
461	461	T intermediate;
462		for (uint r = 0; r < this.~~getR~~owCount(); r++)
	462	for (uint r = 0; r < this.rowCount(); r++)
463	463	{
464	464	intermediate = 0;
465		for (uint i = 0; i < this.~~getC~~olumnCount(); i++)
	465	for (uint i = 0; i < this.columnCount(); i++)
466	466	intermediate += this.opIndex(r,i) * right[i];
467	467	result[r] = intermediate;
…	…
511	511
512	512
513		class Matrix~~View~~(T):MatrixBase!(T)
	513	class Matrix(T):MatrixBase!(T)
514	514	{
	515	private:
	516	static Matrix freelist;
	517	Matrix freelistNext;
	518	static uint freelistCount = 0;
515	519	protected:
516	520	T [] _data;
…	…
520	524	uint _columnStep;
521	525	uint _indexOffset;
522		bool _is~~LapackMode~~; //Column-Major
	526	bool _isColumnMajor; //Column-Major
523	527	public:
524		/***********************************************************************
525
526		Constructor for the MatrixView
	528
	529
	530	/***********************************************************************
	531
	532	Returns a view to an existing matrix, this is part of the freelist
	533	algorithm.
	534
	535	***********************************************************************/
	536
	537	static Matrix allocateView(uint row, uint column, uint rowCount, uint columnCount, Matrix matrix)
	538	{
	539	Matrix!(T) result;
	540	if (freelist)
	541	{
	542	result = freelist;
	543	freelist = result.freelistNext;
	544	result.setView( row, column, rowCount, columnCount, matrix);
	545	freelistCount--;
	546	} else
	547	result = new Matrix!(T) (row, column, rowCount, columnCount, matrix);
	548	return result;
	549	}
	550
	551	/***********************************************************************
	552
	553	Puts the matrix onto a freelist for reuse.
	554
	555	***********************************************************************/
	556
	557	static void deallocate(Matrix m)
	558	{
	559	if (freelistCount < 1000) //Freelist Buffer size.
	560	{
	561	m.freelistNext = freelist;
	562	freelist = m;
	563	freelistCount++;
	564	}
	565	}
	566
	567	/***********************************************************************
	568
	569	Constructor for the Matrix
527	570
528	571	For normal use the array should have the dimensions
…	…
557	600	_allocatedColumnCount = allocatedColumnCount;
558	601	_indexOffset = indexOffset;
559
560		_isLapackMode = columnMajor; //Just discovered that lapackMode == column-major.
561
	602	_isColumnMajor = columnMajor; //Just discovered that lapackMode == column-major.
562	603	_data = data;
563
564		if (columnMajor)
	604	if (_isColumnMajor)
565	605	{
566	606	_columnStep = _allocatedRowCount;
…	…
571	611	_columnStep = 1;
572	612	_rowStep = _allocatedColumnCount;
573		}
574
	613	}
575	614	assert( data.length > (rowCount-1) * _rowStep + (columnCount-1) * _columnStep + _indexOffset );
576	615	}
…	…
579	618	/***********************************************************************
580	619
581		Creates a Matrix~~View that is a view of another MatrixView~~
	620	Creates a Matrix that is a view of another Matrix
582	621
583	622	Params:
…	…
586	625	rowCount= Amount of rows of the view
587	626	columnCount= Amount of columns of the view
588		matrix = The other Matrix~~View~~
	627	matrix = The other Matrix
589	628
590	629	Example:
…	…
596	635	4,8,12,16]);
597	636
598		auto view1 = new Matrix~~View~~!(float) (1,1, 3,3, matrix);
	637	auto view1 = new Matrix!(float) (1,1, 3,3, matrix);
599	638	Stdout( view1 ).newline;
600	639
601		auto view2 = new Matrix~~View~~!(float) (0,0,2,2, view1);
	640	auto view2 = new Matrix!(float) (0,0,2,2, view1);
602	641	Stdout( view2 ).newline;
603	642	---
…	…
621	660
622	661
623		this (uint row, uint column, uint rowCount, uint columnCount, Matrix~~View~~ matrix)
624		{
625		assert( row >= 0 && row < matrix.~~getRowCount() && column >= 0 && column < matrix.getC~~olumnCount(), "index out of bound");
	662	this (uint row, uint column, uint rowCount, uint columnCount, Matrix matrix)
	663	{
	664	assert( row >= 0 && row < matrix.rowCount() && column >= 0 && column < matrix.columnCount(), "index out of bound");
626	665	uint offset = row * matrix._rowStep + column * matrix._columnStep + matrix._indexOffset;
627		this( rowCount, columnCount, matrix.getAllocatedRowCount, matrix.getAllocatedColumnCount, matrix.isColumnMajor, matrix.getData, offset);
628		}
	666	this( rowCount, columnCount, matrix.allocatedRowCount, matrix.allocatedColumnCount, matrix.isColumnMajor, matrix.getData, offset);
	667	}
	668
	669
629	670
630	671
…	…
639	680	Example:
640	681	---
641		auto matrix = new Matrix~~View~~!(double) (2,2,[1,2,3,4]);
	682	auto matrix = new Matrix!(double) (2,2,[1,2,3,4]);
642	683	---
643	684
…	…
654	695
655	696
656		bool isColumnMajor()
657		{
658		return _isLapackMode;
659		}
660
661		uint indexOffset()
	697	/***********************************************************************
	698
	699	Returns whether or not the matrix is in column major format.
	700
	701	***********************************************************************/
	702
	703
	704
	705	final bool isColumnMajor()
	706	{
	707	return _isColumnMajor;
	708	}
	709
	710
	711	/***********************************************************************
	712
	713	You can set the matrix to be in column major format or not.
	714
	715	***********************************************************************/
	716
	717	final void isColumnMajor(bool columnMajor)
	718	{
	719	_isColumnMajor = columnMajor;
	720	if (_isColumnMajor)
	721	{
	722	_columnStep = _allocatedRowCount;
	723	_rowStep = 1;
	724	}
	725	else
	726	{
	727	_columnStep = 1;
	728	_rowStep = _allocatedColumnCount;
	729	}
	730	assert( _data.length > (rowCount-1) * _rowStep + (columnCount-1) * _columnStep + _indexOffset );
	731	}
	732
	733
	734	/***********************************************************************
	735
	736	Returns the offset of the first index (upper left element of the matrix)
	737	in the data array.
	738
	739	***********************************************************************/
	740
	741
	742	final uint indexOffset()
662	743	{
663	744	return _indexOffset;
…	…
682	763
683	764
684		/**
685		* Returns the allocated amount of rows.
686		*/
687
688		uint getAllocatedRowCount()
	765
	766	/***********************************************************************
	767
	768	Sets a view of Matrix!(T).
	769
	770	Params:
	771	row = top left corner
	772	column = top left corner
	773	rowCount = rowCount of the view
	774	columnCount = columnCount of the view
	775	matrix = another matrix to view.
	776
	777	***********************************************************************/
	778
	779
	780
	781	final void setView( uint row, uint column, uint rowCount, uint columnCount, Matrix matrix)
	782	{
	783	assert( row >= 0 && row < matrix.rowCount() && column >= 0 && column < matrix.columnCount(), "index out of bound");
	784	_rowCount = rowCount;
	785	_columnCount = columnCount;
	786	_allocatedRowCount = matrix.allocatedRowCount;
	787	_allocatedColumnCount = matrix.allocatedColumnCount;
	788	_indexOffset = row * matrix._rowStep + column * matrix._columnStep + matrix._indexOffset;
	789	_isColumnMajor = matrix.isColumnMajor;
	790	_data = matrix.getData;
	791	if (_isColumnMajor)
	792	{
	793	_columnStep = _allocatedRowCount;
	794	_rowStep = 1;
	795	}
	796	else
	797	{
	798	_columnStep = 1;
	799	_rowStep = _allocatedColumnCount;
	800	}
	801	assert( _data.length > (rowCount-1) * _rowStep + (columnCount-1) * _columnStep + _indexOffset );
	802	}
	803
	804
	805	/***********************************************************************
	806
	807	Returns the allocated amount of rows.
	808
	809	***********************************************************************/
	810
	811	final uint allocatedRowCount()
689	812	{
690	813	return _allocatedRowCount;
691	814	}
692	815
693		/**
694		* Returns the allocated amount of columns.
695		*/
696
697		uint getAllocatedColumnCount()
	816	/***********************************************************************
	817
	818	Returns the allocated amount of columns.
	819
	820	***********************************************************************/
	821
	822
	823	final uint allocatedColumnCount()
698	824	{
699	825	return _allocatedColumnCount;
700	826	}
701
702
703		/**
704		* Allows you to set the rowCount. It must always be smaller or equal than the allocatdRowCount.
705		*/
706
707
708		/*
709		~~void setRowCount( int rowCount )~~
710		{
711		~~assert( 0 <= rowCount && rowCount <= _allocateRowCount);~~
712		~~_rowCount = rowCount;~~
713		}
714		*/
715
716
717
718
719		/**
720		* Allows you to set the columnCount. It must always be smaller or equal than the allocatedColumnCount.
721		*/
722
723		/*
724
725		~~//Need to check with constructor if data length fit model.~~
726
727		~~void setColumnCount( int columnCount)~~
728		{
729		~~assert(0 <= columnCount && columnCount <= _allocateColumnCount);~~
730		~~_columnCount = columnCount;~~
731		}
732		*/
733	827
734	828
…	…
739	833
740	834
741		T [] getData()
	835	final T [] getData()
742	836	{
743	837	return _data;
…	…
754	848
755	849
756		~~MatrixView~~ opAddAssign( MatrixBase!(T) right)
757		{
758		assert( this.~~getRowCount()==right.getR~~owCount() , "RowCount doesn't match");
759		assert( this.~~getColumnCount() == right.getC~~olumnCount() , "ColumnCount doesn't match");
760		for (uint r = 0; r < ~~getR~~owCount(); r++)
761		for (uint c = 0; c < ~~getC~~olumnCount(); c++)
	850	final Matrix opAddAssign( MatrixBase!(T) right)
	851	{
	852	assert( this.rowCount()==right.rowCount() , "RowCount doesn't match");
	853	assert( this.columnCount() == right.columnCount() , "ColumnCount doesn't match");
	854	for (uint r = 0; r < rowCount(); r++)
	855	for (uint c = 0; c < columnCount(); c++)
762	856	this[r,c] = this[r,c] + right[r,c];
763	857	return this;
…	…
776	870
777	871
778		~~MatrixView~~ opSubAssign(MatrixBase!(T) right)
779		{
780		assert( this.~~getRowCount() == right.getR~~owCount() , " RowCount doesn't match");
781		assert( this.~~getColumnCount() == right.getC~~olumnCount() , "ColumnCount doesn't match");
782		for (uint r = 0; r < ~~getR~~owCount(); r++)
783		for (uint c = 0; c < ~~getC~~olumnCount(); c++)
	872	final Matrix opSubAssign(MatrixBase!(T) right)
	873	{
	874	assert( this.rowCount() == right.rowCount() , " RowCount doesn't match");
	875	assert( this.columnCount() == right.columnCount() , "ColumnCount doesn't match");
	876	for (uint r = 0; r < rowCount(); r++)
	877	for (uint c = 0; c < columnCount(); c++)
784	878	this[r,c] = this[r,c] - right[r,c];
785	879	return this;
…	…
792	886
793	887
794		~~MatrixView~~ opMulAssign(real scalar)
795		{
796		for (uint r = 0; r < ~~getR~~owCount(); r++)
797		for (uint c = 0; c < ~~getC~~olumnCount; c++)
	888	final Matrix opMulAssign(real scalar)
	889	{
	890	for (uint r = 0; r < rowCount(); r++)
	891	for (uint c = 0; c < columnCount; c++)
798	892	this[r,c] = this[r,c] * scalar;
799	893	return this;
…	…
804	898	/*
805	899	* Provides access to an element of the matrix.
	900	*
806	901	* Example:
807	902	* float test = matrix[row, column];
808	903	*/
809	904
810		~~override~~ T opIndex(uint row, uint column)
811		{
812		assert( row >= 0 && row < this.~~getRowCount() && column >= 0 && column < this.getC~~olumnCount(), "index out of bound");
	905	final T opIndex(uint row, uint column)
	906	{
	907	assert( row >= 0 && row < this.rowCount() && column >= 0 && column < this.columnCount(), "index out of bound");
813	908	return this._data[row_rowStep + column_columnStep + _indexOffset];
814	909	}
815	910
816	911	/**
817		* Assign a value to row and column
	912	* Assign a value to row and column.
	913	*
818	914	* Example:
819	915	* matrix[row, column] = 10;
820	916	*/
821	917
822		~~override~~ T opIndexAssign( T value, uint row, uint column)
823		{
824		assert( row >= 0 && row < this.~~getRowCount() && column >= 0 && column < this.getC~~olumnCount(), "index out of bound");
	918	final T opIndexAssign( T value, uint row, uint column)
	919	{
	920	assert( row >= 0 && row < this.rowCount() && column >= 0 && column < this.columnCount(), "index out of bound");
825	921	this._data[row_rowStep + column_columnStep+ _indexOffset] = value;
826	922	return value;
…	…
833	929	*/
834	930
835		void copyFrom(MatrixBase!(T) right)
836		{
837		assert( right.~~getRowCount() >= this.getR~~owCount(), "The size of the current matrix is not large enough");
838		assert( right.~~getColumnCount() >= this.getC~~olumnCount(), "The size of the current matrix is not large enough");
839
840		for (uint r = 0; r < this.~~getR~~owCount(); r++)
841		{
842		for (uint c = 0; c < this.~~getC~~olumnCount(); c++)
	931	final void copyFrom(MatrixBase!(T) right)
	932	{
	933	assert( right.rowCount() >= this.rowCount(), "The size of the current matrix is not large enough");
	934	assert( right.columnCount() >= this.columnCount(), "The size of the current matrix is not large enough");
	935
	936	for (uint r = 0; r < this.rowCount(); r++)
	937	{
	938	for (uint c = 0; c < this.columnCount(); c++)
843	939	this[r,c] = right[r,c];
844	940	}
…	…
846	942
847	943
848		void setIdentity()
849		{
850		for (uint r=0; r < ~~getR~~owCount; r++)
851		for (uint c= 0; c < ~~getC~~olumnCount; c++)
	944	final void setIdentity()
	945	{
	946	for (uint r=0; r < rowCount; r++)
	947	for (uint c= 0; c < columnCount; c++)
852	948	this[r,c] = (r==c?1:0);
853	949	}
…	…
858	954	*/
859	955
860		void setScalingMatrix2x2(T xScale, T yScale)
861		{
862		assert( ~~getC~~olumnCount >= 2, "Your matrix is not big enough");
863		assert( ~~getR~~owCount >= 2,"Your matrix is not big enough");
	956	final void setScalingMatrix2x2(T xScale, T yScale)
	957	{
	958	assert( columnCount >= 2, "Your matrix is not big enough");
	959	assert( rowCount >= 2,"Your matrix is not big enough");
864	960
865	961	this[0,0] = xScale; this[0,1]=0;
…	…
873	969	*/
874	970
875		void setRotationMatrix2x2(T angle)
876		{
877		assert( ~~getC~~olumnCount >= 2, "Your matrix is not big enough");
878		assert( ~~getR~~owCount >= 2,"Your matrix is not big enough");
	971	final void setRotationMatrix2x2(T angle)
	972	{
	973	assert( columnCount >= 2, "Your matrix is not big enough");
	974	assert( rowCount >= 2,"Your matrix is not big enough");
879	975
880	976	version (Tango)
…	…
892	988
893	989
894		void setLookAtMatrix4x4( ref VectorFixed!(T, cast(uint) 3) center, ref VectorFixed!(T, cast(uint) 3) z,ref VectorFixed!(T, cast(uint) 3) y )
	990	final void setLookAtMatrix4x4( ref VectorFixed!(T, cast(uint) 3) center, ref VectorFixed!(T, cast(uint) 3) z,ref VectorFixed!(T, cast(uint) 3) y )
895	991	{
896	992	z = z.normal();
…	…
900	996	y = y.normal();
901	997
902		assert( ~~getR~~owCount >= 4);
903		assert( ~~getC~~olumnCount >= 4);
	998	assert( rowCount >= 4);
	999	assert( columnCount >= 4);
904	1000
905	1001	auto matrix = new Matrix!(T)(4,4);
…	…
947	1043	*/
948	1044
949		void setRotationMatrix3x3( T angle, VectorFixed!(T, 3) normal)
950		{
951		assert( ~~getC~~olumnCount >= 3, "Your matrix is not big enough columnwise");
952		assert( ~~getR~~owCount >= 3,"Your matrix is not big enough rowwise");
	1045	final void setRotationMatrix3x3( T angle, VectorFixed!(T, 3) normal)
	1046	{
	1047	assert( columnCount >= 3, "Your matrix is not big enough columnwise");
	1048	assert( rowCount >= 3,"Your matrix is not big enough rowwise");
953	1049	version (Tango)
954	1050	{
…	…
978	1074	}
979	1075
980		void setRotationMatrix4x4( T angle, VectorFixed!(T, 3) normal)
981		{
982		assert( ~~getC~~olumnCount >= 4, "Your matrix is not big enough columnwise");
983		assert( ~~getR~~owCount >= 4,"Your matrix is not big enough rowwise");
	1076	final void setRotationMatrix4x4( T angle, VectorFixed!(T, 3) normal)
	1077	{
	1078	assert( columnCount >= 4, "Your matrix is not big enough columnwise");
	1079	assert( rowCount >= 4,"Your matrix is not big enough rowwise");
984	1080	setRotationMatrix3x3(angle, normal);
985	1081	this[0,3] = 0; //translation X
…	…
1028	1124	*/
1029	1125
1030		void setFrustumMatrix4x4( T left, T right, T bottom, T top, T near, T far)
	1126	final void setFrustumMatrix4x4( T left, T right, T bottom, T top, T near, T far)
1031	1127	{
1032		assert( ~~getC~~olumnCount >= 4, "Your matrix is not big enough columnwise");
1033		assert( ~~getR~~owCount >= 4,"Your matrix is not big enough rowwise");
	1128	assert( columnCount >= 4, "Your matrix is not big enough columnwise");
	1129	assert( rowCount >= 4,"Your matrix is not big enough rowwise");
1034	1130
1035	1131	this[0,0] =1.0/(right-left);this[0,1] =0;this[0,2] = 1.0/(right-left);this[0,3] =0;
…	…
1082	1178	*/
1083	1179
1084		void setOrthoMatrix4x4(T left, T right, T bottom, T top, T near, T far)
1085		{
1086		assert( ~~getC~~olumnCount >= 4, "Your matrix is not big enough columnwise");
1087		assert( ~~getR~~owCount >= 4,"Your matrix is not big enough rowwise");
	1180	final void setOrthoMatrix4x4(T left, T right, T bottom, T top, T near, T far)
	1181	{
	1182	assert( columnCount >= 4, "Your matrix is not big enough columnwise");
	1183	assert( rowCount >= 4,"Your matrix is not big enough rowwise");
1088	1184
1089	1185	T width = right - left;
…	…
1103	1199	*/
1104	1200
1105		void setRotationXMatrix3x3(T angle)
1106		{
1107		assert( ~~getC~~olumnCount >= 3, "Your matrix is not big enough columnwise");
1108		assert( ~~getR~~owCount >= 3,"Your matrix is not big enough rowwise");
	1201	final void setRotationXMatrix3x3(T angle)
	1202	{
	1203	assert( columnCount >= 3, "Your matrix is not big enough columnwise");
	1204	assert( rowCount >= 3,"Your matrix is not big enough rowwise");
1109	1205	T s = sin(angle);
1110	1206	T c = cos(angle);
…	…
1120	1216	*/
1121	1217
1122		void setRotationYMatrix3x3(T angle)
1123		{
1124		assert( ~~getC~~olumnCount >= 3, "Your matrix is not big enough columnwise");
1125		assert( ~~getR~~owCount >= 3,"Your matrix is not big enough rowwise");
	1218	final void setRotationYMatrix3x3(T angle)
	1219	{
	1220	assert( columnCount >= 3, "Your matrix is not big enough columnwise");
	1221	assert( rowCount >= 3,"Your matrix is not big enough rowwise");
1126	1222	T s = sin(angle);
1127	1223	T c = cos(angle);
…	…
1137	1233	*/
1138	1234
1139		void setRotationZMatrix3x3(T angle)
1140		{
1141		assert( ~~getC~~olumnCount >= 3, "Your matrix is not big enough columnwise");
1142		assert( ~~getR~~owCount >= 3,"Your matrix is not big enough rowwise");
	1235	final void setRotationZMatrix3x3(T angle)
	1236	{
	1237	assert( columnCount >= 3, "Your matrix is not big enough columnwise");
	1238	assert( rowCount >= 3,"Your matrix is not big enough rowwise");
1143	1239	T s = sin(angle);
1144	1240	T c = cos(angle);
…	…
1155	1251	*/
1156	1252
1157		void setTranslationMatrix3x3(T x, T y)
1158		{
1159		assert( ~~getC~~olumnCount >= 3, "Your matrix is not big enough columnwise");
1160		assert( ~~getR~~owCount >= 3,"Your matrix is not big enough rowwise");
	1253	final void setTranslationMatrix3x3(T x, T y)
	1254	{
	1255	assert( columnCount >= 3, "Your matrix is not big enough columnwise");
	1256	assert( rowCount >= 3,"Your matrix is not big enough rowwise");
1161	1257	this[0,0]=1; this[0,1]=0; this[0,2]=x;
1162	1258	this[1,0]=0; this[1,1]=1; this[1,2]=y;
…	…
1165	1261
1166	1262
1167		void setTranslationMatrix4x4(T x, T y, T z)
1168		{
1169		assert( ~~getC~~olumnCount >= 4);
1170		assert( ~~getR~~owCount >= 4);
	1263	final void setTranslationMatrix4x4(T x, T y, T z)
	1264	{
	1265	assert( columnCount >= 4);
	1266	assert( rowCount >= 4);
1171	1267	this[0,0]=1; this[0,1]=0; this[0,2]=0; this[0,3]= x;
1172	1268	this[1,0]=0; this[1,1]=1; this[1,2]=0; this[1,3]= y;
…	…
1180	1276	*/
1181	1277
1182		void setScalingMatrix3x3(T xScale, T yScale, T zScale)
1183		{
1184		assert( ~~getC~~olumnCount >= 3, "Your matrix is not big enough columnwise");
1185		assert( ~~getR~~owCount >= 3,"Your matrix is not big enough rowwise");
	1278	final void setScalingMatrix3x3(T xScale, T yScale, T zScale)
	1279	{
	1280	assert( columnCount >= 3, "Your matrix is not big enough columnwise");
	1281	assert( rowCount >= 3,"Your matrix is not big enough rowwise");
1186	1282	this[0,0]=xScale; this[0,1]=0; this[0,2]=0;
1187	1283	this[1,0]=0; this[1,1]=yScale; this[1,2]=0;
…	…
1199	1295	*******************************************************************************/
1200	1296
1201		class Matrix(T):MatrixView!(T)
	1297
	1298	/*
	1299
	1300	class Matrix(T):Matrix!(T)
1202	1301	{
1203	1302	public:
1204
1205		/***********************************************************************
1206
1207		~~Constructor for the column-major rectangular matrix.~~
1208
1209		~~Example:~~
1210		~~---~~
1211		~~auto mat1 = new Matrix!(float)(2,2, [1,2,3,4]);~~
1212		~~---~~
1213
1214		***********************************************************************/
1215
1216	1303	this(uint rowCount, uint columnCount, T [] data)
1217	1304	{
…	…
1225	1312	}
1226	1313	}
	1314
	1315	*/
1227	1316
1228	1317
…	…
1246	1335	super(ROWCOUNT, COLUMNCOUNT);
1247	1336	}
1248
1249		override Matrix!(T) opMul(MatrixBase!(T) right)
1250		{
1251		return super.opMul(right);
1252		}
1253
1254		/**
1255		* Multiplication of a matrix and a vector.
1256		*/
1257
1258		VectorFixed!(T, ROWCOUNT) opMul( VectorFixed!(T, COLUMNCOUNT) right)
	1337
	1338	/***********************************************************************
	1339
	1340	Multiplication of a matrix with a vector, returns a vector.
	1341
	1342	***********************************************************************/
	1343
	1344	final VectorFixed!(T, ROWCOUNT) opMul( VectorFixed!(T, COLUMNCOUNT) right)
1259	1345	{
1260	1346	VectorFixed!(T, ROWCOUNT) result;
1261	1347	T intermediate;
1262		for (uint r = 0; r < this.~~getR~~owCount(); r++)
	1348	for (uint r = 0; r < this.rowCount(); r++)
1263	1349	{
1264	1350	intermediate = 0;
1265		for (uint ~~i = 0; i < this.getColumnCount(); i~~++)
1266		intermediate += ~~_data[i * ROWCOUNT+r] * right[i~~];
	1351	for (uint c = 0; c < this.columnCount(); c++)
	1352	intermediate += this[r,c] * right[c];
1267	1353	result[r] = intermediate;
1268	1354	}
1269	1355	return result;
1270		}
1271
1272		/**
1273		* Returns the transpose matrix.
1274		*/
1275
1276		MatrixSized!(T, COLUMNCOUNT, ROWCOUNT) newTransposeMatrix()
	1356
	1357	}
	1358
	1359	/***********************************************************************
	1360
	1361	Returns a new Transposed Matrix.
	1362
	1363	***********************************************************************/
	1364
	1365	final MatrixSized!(T, COLUMNCOUNT, ROWCOUNT) newTransposeMatrix()
1277	1366	{
1278	1367	MatrixSized!(T, COLUMNCOUNT, ROWCOUNT) result = new MatrixSized!(T, COLUMNCOUNT, ROWCOUNT);
1279		for (uint r = 0; r < this.~~getR~~owCount; r++)
1280		for (uint c = 0; c < ~~getC~~olumnCount; c++)
	1368	for (uint r = 0; r < this.rowCount; r++)
	1369	for (uint c = 0; c < columnCount; c++)
1281	1370	result[c,r] = this[r,c];
1282	1371	return result;
1283	1372	}
1284	1373
1285		~~override T opIndex(uint row, uint column)~~
1286		{
1287		~~assert( row >= 0 && row < this.getRowCount() && column >= 0 && column < this.getColumnCount(), "index out of bound");~~
1288		~~//return this._data[row_rowStep + column_columnStep];~~
1289		~~assert(_rowStep == 1);~~
1290		~~return this._data[row+column*_columnStep];~~
1291		}
1292
1293		~~override T opIndexAssign( T value, uint row, uint column)~~
1294		{
1295		~~assert( row >= 0 && row < this.getRowCount() && column >= 0 && column < this.getColumnCount(), "index out of bound");~~
1296		~~assert(_rowStep == 1);~~
1297		~~this._data[row + column*_columnStep] = value;~~
1298		~~return value;~~
1299		}
1300	1374
1301	1375	static if (ROWCOUNT == COLUMNCOUNT)
1302	1376	{
1303		MatrixSized opMul(MatrixSized right)
1304		{
1305		assert( this.getColumnCount() == right.getRowCount(), "Incorrect matrix dimensions for multiplication");
	1377
	1378
	1379	/***********************************************************************
	1380
	1381	Multiplication of two squared matrices of the same dimensions, where
	1382	rowCount == columnCount.
	1383
	1384	***********************************************************************/
	1385
	1386	final MatrixSized opMul(MatrixSized right)
	1387	{
	1388	assert( this.columnCount() == right.rowCount(), "Incorrect matrix dimensions for multiplication");
1306	1389	MatrixSized result = new MatrixSized;
1307	1390	T intermediate;
1308		for (uint r = 0; r < this.~~getR~~owCount(); r++)
1309		for (uint c = 0; c < right.~~getC~~olumnCount; c++)
	1391	for (uint r = 0; r < this.rowCount(); r++)
	1392	for (uint c = 0; c < right.columnCount; c++)
1310	1393	{
1311	1394	intermediate = 0;
1312		for (uint i = 0; i < this.~~getC~~olumnCount(); i++)
	1395	for (uint i = 0; i < this.columnCount(); i++)
1313	1396	intermediate += this.opIndex(r,i) * right[i,c];
1314	1397	result[r,c] = intermediate;
…	…
1374	1457	{
1375	1458
1376		assert(this.~~getRowCount == right.getRowCount() && this.getColumnCount() == right.getC~~olumnCount(),
	1459	assert(this.rowCount == right.rowCount() && this.columnCount() == right.columnCount(),
1377	1460	"dimensions do not match");
1378		MatrixDiagonal result = new MatrixDiagonal(~~getRowCount, getC~~olumnCount());
	1461	MatrixDiagonal result = new MatrixDiagonal(rowCount, columnCount());
1379	1462
1380	1463	for (int i = 0; i < this.getData().length; i++)
…	…
1385	1468	override T opIndex(uint row, uint column)
1386	1469	{
1387		assert( row >= 0 && row < ~~getRowCount() && column >= 0 && column < getC~~olumnCount() , "index out of bound");
	1470	assert( row >= 0 && row < rowCount() && column >= 0 && column < columnCount() , "index out of bound");
1388	1471	if (row == column) return _data[row];
1389	1472	return 0;
…	…
1392	1475	override T opIndexAssign( T value, uint row, uint column)
1393	1476	{
1394		assert( row >= 0 && row < ~~getRowCount() && column >= 0 && column < getC~~olumnCount() , "index out of bound");
	1477	assert( row >= 0 && row < rowCount() && column >= 0 && column < columnCount() , "index out of bound");
1395	1478	_data[row] = value;
1396	1479	return value;
…	…
1404	1487	uint getMinimumDimension()
1405	1488	{
1406		return ~~getRowCount() < getColumnCount()? getRowCount(): getC~~olumnCount();
	1489	return rowCount() < columnCount()? rowCount(): columnCount();
1407	1490	}
1408	1491
…	…
1473	1556	assert(row >= 0, "row");
1474	1557	assert(column >= 0, "column");
1475		assert(row < ~~getR~~owCount, "row");
1476		assert(column < ~~getC~~olumnCount, "column");
	1558	assert(row < rowCount, "row");
	1559	assert(column < columnCount, "column");
1477	1560	if (row <= column)
1478	1561	return _data[row + column*(column+cast(uint)1)/cast(uint)2];
…	…
1493	1576	assert(row >= 0, "row");
1494	1577	assert(column >= 0, "column");
1495		assert(row < ~~getR~~owCount, "row");
1496		assert(column < ~~getC~~olumnCount, "column");
	1578	assert(row < rowCount, "row");
	1579	assert(column < columnCount, "column");
1497	1580
1498	1581	assert(row <= column, "row > column");
…	…
1568	1651	assert(row >= 0, "row");
1569	1652	assert(column >= 0, "column");
1570		assert(row < ~~getR~~owCount, "row");
1571		assert(column < ~~getC~~olumnCount, "column");
	1653	assert(row < rowCount, "row");
	1654	assert(column < columnCount, "column");
1572	1655	if (row >= column)
1573		return _data[row + ~~getR~~owCount * column - column * (column + 1) / 2];
	1656	return _data[row + rowCount * column - column * (column + 1) / 2];
1574	1657	else
1575	1658	return 0;
…	…
1588	1671	assert(row >= 0, "row");
1589	1672	assert(column >= 0, "column");
1590		assert(row < ~~getR~~owCount, "row");
1591		assert(column < ~~getC~~olumnCount, "column");
	1673	assert(row < rowCount, "row");
	1674	assert(column < columnCount, "column");
1592	1675
1593	1676	assert(row >= column, "column > row");
1594	1677	//return _data[row + (2size-column)(column -1) /2 ] = value;
1595		_data[ row + ~~getR~~owCount * column - column*(column + 1) / 2] = value;
	1678	_data[ row + rowCount * column - column*(column + 1) / 2] = value;
1596	1679	return value;
1597	1680	}
…	…
1674	1757	assert( trianglMul1 == new Matrix!(float) (2,2,[1,2,0,3]) * new Matrix!(float) (2,2,[1,2,0,3]));
1675	1758
1676		auto matView1 = new Matrix~~View~~!(float)(2,2);
	1759	auto matView1 = new Matrix!(float)(2,2);
1677	1760	matView1.copyFrom( mat1 );
1678	1761	assert( matView1 == mat1 );

trunk/mathematics/numerical/Vector.d

r8	r10
42	42	//Integer values check out Point(T) or PointFixed(T)
43	43
44
45	44	/**
46	45	* A class for representing a vector. This provides a means to inherit from the class.
…	…
53	52	private:
54	53	T [] _data;
55		uint _stepSize;
56		uint _startIndex;
57		uint _count;
	54	uint _stepSize; //The step to the next element
	55	uint _startIndex; //The index this vector starts with. //startindex can be removed.
	56	uint _count; //Number of elements in this vector.
58	57	public:
	58
59	59
60	60	/**
…	…
106	106	}
107	107
108		Vector opAdd(Vector right)
	108	/**
	109	Returns the pointer to the first element of the vector.
	110	*/
	111
	112	final T * iterStart()
	113	{
	114	return &_data[_startIndex];
	115	}
	116
	117	/**
	118	Returns the pointer to the last element of the vector.
	119	*/
	120
	121	final T * iterStop()
	122	{
	123	return &_data[_startIndex + _count - 1];
	124	}
	125
	126	final Vector opAdd(Vector right)
109	127	{
110	128	assert(_count == right._count);
111	129	Vector result = new Vector(_count);
	130	T * resultStop = result.iterStop;
	131	for (T * resultIter = result.iterStart, dataIter = this.iterStart, rightIter = right.iterStart; resultIter <= resultStop; resultIter++, dataIter++, rightIter++)
	132	resultIter = dataIter + *rightIter;
	133	return result;
	134	}
	135
	136	/*
	137	final Vector opAdd(Vector right)
	138	{
	139	assert(_count == right._count);
	140	Vector result = new Vector(_count);
112	141	for (int x = 0; x < _count; x++)
113	142	result[x] = this[x] + right[x];
114	143	return result;
115	144	}
116
117
118		Vector opAddAssign(Vector right)
	145	*/
	146
	147
	148	final Vector opAddAssign(Vector right)
119	149	{
120	150	assert(_count == right._count );
…	…
125	155
126	156
127		Vector opSub(Vector right)
	157	final Vector opSub(Vector right)
128	158	{
129	159	assert(_count == right._count);
…	…
135	165
136	166
137		Vector opSubAssign(Vector right)
	167	final Vector opSubAssign(Vector right)
138	168	{
139	169	assert(_count == right._count);
…	…
150	180
151	181
152		Vector opMul( Vector right )
	182	final Vector opMul( Vector right )
153	183	{
154	184	assert(getCount == right.getCount );
…	…
170	200
171	201
172		Vector opMulAssign( Vector right)
	202	final Vector opMulAssign( Vector right)
173	203	{
174	204	assert( getCount == right.getCount );
…	…
180	210	}
181	211
182		Vector opMul(real scale)
	212	final Vector opMul(real scale)
183	213	{
184	214	Vector result = this.clone();
…	…
187	217	}
188	218
189		Vector opMulAssign( real scale )
	219	final Vector opMulAssign( real scale )
190	220	{
191	221	for (int x = 0; x < _count; x++)
…	…
194	224	}
195	225
196		Vector opMul_r(real scale)
	226	final Vector opMul_r(real scale)
197	227	{
198	228	Vector result = this.clone();
…	…
205	235	*/
206	236
207		real norm()
	237	final real norm()
208	238	{
209	239	real result = 0;
…	…
218	248	*/
219	249
220		Vector normal()
	250	final Vector normal()
221	251	{
222	252	Vector result = this.clone();
…	…
230	260	*/
231	261
232		void normalize()
	262	final void normalize()
233	263	{
234	264	real normValue = cast(real)1.0/norm;
…	…
237	267	}
238	268
239		Vector opNeg()
	269	final Vector opNeg()
240	270	{
241	271	Vector result = this.clone();
…	…
246	276	}
247	277
248		T opIndex(int i)
	278	final T opIndex(int i)
249	279	{
250	280	assert(0 <= i && i < _count, "class Vector: Index out of bound");
…	…
253	283	}
254	284
255		T opIndexAssign( T value, int i)
	285	final T opIndexAssign( T value, int i)
256	286	{
257	287	assert(0 <= i && i < _count, "class Vector: Index out of bound");
…	…
261	291	}
262	292
263		int opEquals(Vector right)
	293	final int opEquals(Vector right)
264	294	{
265	295	assert(_count == right._count);
…	…
271	301	}
272	302
273		int opCmp(Vector right )
	303	final int opCmp(Vector right )
274	304	{
275	305	assert(_count == right._count, "class Vector: _count are not the same");
…	…
287	317	*/
288	318
289		int getCount()
	319	final int getCount()
290	320	{
291	321	return _count;
292	322	}
293	323
294		char [] toString()
	324	final char [] toString()
295	325	{
296	326	char[] result = "";
…	…
311	341	*/
312	342
313		T scalarMul( Vector right )
	343	final T scalarMul( Vector right )
314	344	{
315	345	assert( _count == right._count);
…	…
324	354	*/
325	355
326		Vector vectorMul( Vector right)
	356	final Vector vectorMul( Vector right)
327	357	{
328	358	assert(this._count == 3 && right._count == 3);
…	…
341	371
342	372
343		Vector clone() //not an exact clone no more
	373	final Vector clone() //not an exact clone no more
344	374	{
345	375	Vector result = new Vector(_count);
…	…
362	392
363	393
364		void mulAssign( MatrixBase!(T) matrix, Vector vec)
	394	final void mulAssign( MatrixBase!(T) matrix, Vector vec)
365	395	{
366	396	// MxN Nx1 = Mx1
367		assert( matrix.~~getC~~olumnCount() == vec.getCount(), "matrix size and right vector size do not match");
368		assert(_data.length >= matrix.~~getR~~owCount(), "length of the vector data is not high enough");
369		_count = matrix.~~getR~~owCount();
	397	assert( matrix.columnCount() == vec.getCount(), "matrix size and right vector size do not match");
	398	assert(_data.length >= matrix.rowCount(), "length of the vector data is not high enough");
	399	_count = matrix.rowCount();
370	400	T intermediate;
371		for (int r = 0; r < matrix.~~getR~~owCount(); r++)
	401	for (int r = 0; r < matrix.rowCount(); r++)
372	402	{
373	403	intermediate = 0;
374		for (int c = 0; c < matrix.~~getC~~olumnCount(); c++)
	404	for (int c = 0; c < matrix.columnCount(); c++)
375	405	intermediate += matrix[r,c] * vec[c];
376	406	this[r] = intermediate;
377	407	}
	408	}
	409
	410	//TODO:
	411
	412	final void mulAssign(Matrix!(T) matrix, Vector vec)
	413	{
	414
378	415	}
379	416	}

trunk/mathematics/numerical/VectorFixed.d

r8	r10
34	34
35	35
	36	char [] assertIterRange(char [] iterator, char [] list)
	37	{
	38	return "assert(" ~ iterator ~ ">=&" ~ list ~ "[0] && "
	39	~ iterator ~ "<=&" ~ list ~ "[$-1], \"Iterator out of range\");";
	40	}
	41
	42
	43
	44
	45
36	46
37	47	/**
…	…
47	57	Adds two homogeneous vectors together.
48	58	*/
49
50	59	VectorFixed homAdd(ref VectorFixed right)
51	60	{
…	…
196	205
197	206	//Start normal routines
	207
	208
	209
	210	//TODO: VectorFixed opAdd( ref const VectorFixed right )
	211
	212	VectorFixed opAdd( ref VectorFixed right )
	213	{
	214	VectorFixed result;
	215	T * resultStop = &result.data[$-1];
	216	for (T * resultIter = &result.data[0], rightIter = &right.data[0], dataIter = &data[0]; resultIter <= resultStop; resultIter++, rightIter++, dataIter++)
	217	{
	218	mixin(assertIterRange("resultIter", "result.data"));
	219	mixin(assertIterRange("rightIter", "right.data"));
	220	mixin(assertIterRange("dataIter", "data"));
	221	resultIter = dataIter + *rightIter;
	222	}
	223	return result;
	224	}
	225
	226
198	227
199	228	/**
…	…
201	230	*/
202	231
	232
	233	/*
203	234	VectorFixed opAdd(ref VectorFixed right)
204	235	{
…	…
208	239	return result;
209	240	}
	241	*/
210	242
211	243	/**

trunk/mathematics/numerical/algorithms/LuDecomposition.d

r9	r10
53	53	int luDecomposition(T)(Matrix!(T) aMatrix, Matrix!(T) bMatrix, out int [] ipivot )
54	54	{
55		assert( aMatrix.~~getRowCount() == aMatrix.getC~~olumnCount(), "not a square matrix");
56		assert( aMatrix.~~getRowCount() == bMatrix.getR~~owCount(), "first matrix must be nxn second nxnhrs");
	55	assert( aMatrix.rowCount() == aMatrix.columnCount(), "not a square matrix");
	56	assert( aMatrix.rowCount() == bMatrix.rowCount(), "first matrix must be nxn second nxnhrs");
57	57	//Check if matrix is symmetric.
58	58
59		int n = aMatrix.~~getR~~owCount();
60		int nrhs = bMatrix.~~getC~~olumnCount();
	59	int n = aMatrix.rowCount();
	60	int nrhs = bMatrix.columnCount();
61	61	T [] a = aMatrix.getData(); //NxN matrix
62	62	int lda = aMatrix.getAllocatedRowCount();
…	…
92	92
93	93	unittest {
	94	Stdout("LU Decomposition").newline;
94	95	Matrix!(double) aMatrix = new Matrix!(double)(2,2,[1,2,3,4]);
95	96	Matrix!(double) bMatrix = new Matrix!(double)(2,2,[1,0,0,1]);
	97	Stdout("Input:").newline;
	98	Stdout( aMatrix ).newline;
	99	Stdout( bMatrix ).newline;
	100	Stdout("Output:").newline;
96	101	int [] rowPivot;
97		luDecomposition!(double) (aMatrix, bMatrix, rowPivot);
	102	int result = luDecomposition!(double) (aMatrix, bMatrix, rowPivot);
	103	if (result != 0) Stdout("An error has occured").newline;
	104	Stdout("resultCode:")(result).newline;
98	105	Stdout( aMatrix ).newline;
99	106	Stdout( bMatrix ).newline;

trunk/mathematics/numerical/algorithms/SingularValueDecomposition.d

r8	r10
63	63	void singularValueDecomposition(T)(Matrix!(T) matrix,out Matrix!(T) uMatrix,out MatrixDiagonal!(T) sMatrix, out Matrix!(T) vTMatrix)
64	64	{
65		uMatrix = new Matrix!(T) (matrix.~~getRowCount(), matrix.getR~~owCount());
66		sMatrix = new MatrixDiagonal!(T) (matrix.~~getRowCount(), matrix.getC~~olumnCount());
67		vTMatrix = new Matrix!(T) (matrix.~~getColumnCount(), matrix.getC~~olumnCount());
	65	uMatrix = new Matrix!(T) (matrix.rowCount(), matrix.rowCount());
	66	sMatrix = new MatrixDiagonal!(T) (matrix.rowCount(), matrix.columnCount());
	67	vTMatrix = new Matrix!(T) (matrix.columnCount(), matrix.columnCount());
68	68
69	69	T bestLWork;
70	70	char jobu='A';
71	71	char jobvt='A';
72		int m = matrix.~~getR~~owCount();
73		int n = matrix.~~getC~~olumnCount();
	72	int m = matrix.rowCount();
	73	int n = matrix.columnCount();
74	74	T [] a = matrix.getData(); //is this the correct address?
75	75	int lda = matrix.getAllocatedRowCount(); //LDA needs to figure out more about this.
76	76	T [] s= sMatrix.getData(); //min(m,n)
77	77	T [] u = uMatrix.getData(); //LDU, M for 'A'
78		int ldu = uMatrix.~~getR~~owCount(); //leading dimensions of u
	78	int ldu = uMatrix.rowCount(); //leading dimensions of u
79	79	T [] vt = vTMatrix.getData(); //N x N for 'A'
80		int ldvt = vTMatrix.~~getR~~owCount(); //N for 'A'
	80	int ldvt = vTMatrix.rowCount(); //N for 'A'
81	81	int lwork = -1; //dimensions of array work //Two Calls Find Optimal LWORK
82	82	int info; //0 success, <0 -i argument had illegal value

mathematics

Navigation

Legend:

trunk/dsss.conf

trunk/mathematics/numerical/Matrix.d

trunk/mathematics/numerical/Vector.d

trunk/mathematics/numerical/VectorFixed.d

trunk/mathematics/numerical/algorithms/LuDecomposition.d

trunk/mathematics/numerical/algorithms/SingularValueDecomposition.d

Download in other formats: