Entwicklung einer einfachen Matrizenklasse

[SeppJ]

@Finnegan sagte in Entwicklung einer einfachen Matrizenklasse:

Ansonsten sehe ich keinen Schaden das so zu schreiben und verstehe nicht ganz warum das überhaupt eine Diskussion Wert ist. Ich zwinge ja niemanden, das auch so zu machen.

Weil das komisch inkonsistenter Cargo-Cult ist. Wenn du irgendwie argumentierst, dass das sinnvoll ist, dann musst du std::string str = static_cast<std::string>("Hallo"); schreiben. Gerade dann, denn da macht der Cast ja tatsächlich etwas. Tust du aber nicht, weil das offensichtlicher Unsinn ist. Und genauso offensichtlicher Unsinn ist der Fall hier auch, wo du es plötzlich doch machen willst. Das heißt, du machst bloß irgendwelche Muster nach, die du mal gesehen hast, aber ohne sie verstanden zu haben.

Mal anders gesagt: Was würde passieren, wenn da ut(n, n) = static_cast<float>(1) stünde, aber val_t kein float ist? Wird da deiner Meinung nach dann ein float gespeichert?

[zeropage]

Wenn das so gemacht wird, ist es ok?

        template <typename val_t>
	Mat2D_T<val_t> UnitMat(const std::size_t size)
	{
		Mat2D_T<val_t> ut(size, size);

		for (std::size_t n = 0; n < size; ++n)
			ut(n, n) = { 1 };

		return ut;
	}
	template <typename val_t>
	Mat2D_T<val_t> Point2D(const val_t x, const val_t y)
	{
		Mat2D_T<val_t> p{ { x, y, { 1 } } };
		return p;
	}
	template <typename val_t>
	Mat2D_T<val_t> Point2D(const std::initializer_list<val_t>& list)
	{
		Mat2D_T<val_t> p{ list };
		return p;
	}
	template <typename val_t>
	Mat2D_T<val_t> Point3D(const val_t x, const val_t y, const val_t z) 
	{
		Mat2D_T<val_t> p{ { x, y, z } };
		return p;
	}
	template <typename val_t>
	Mat2D_T<val_t> Point3D(const std::initializer_list<val_t> list)
	{
		Mat2D_T<val_t> p{ list };
		return p;
	}

[SeppJ]

Lässt sich Zeile 7 überhaupt übersetzen?

[zeropage]

Huch, Schreck! Sollte es nicht?

Mat2D_T<float> ut = Mat2D::UnitMat<float>(4);
Mat2D::print(ut);

Unter Visual Studio 19 geht das bei mir. Ich kann aber gerne die geschweiften Klammern weglassen, wenn das wirklich kein Ding ist.
Oder eben

ut(n, n) = val_t{ 1 };

[SeppJ]

@zeropage sagte in Entwicklung einer einfachen Matrizenklasse:

Huch, Schreck! Sollte es nicht?

Jetzt, wo ich das ganze Programm gelesen habe: Doch, sollte gehen. Aber ist eine sehr ungewöhnliche Schreibweise.

Anders herum: Man sagt zurecht, wenn man Code schreibt, muss man von jedem einzelnen Zeichen genau wissen, wo und warum man es setzt. Kannst du erklären, wieso du die Klammern gesetzt hast?

[zeropage]

An sich kam das von val_t v = {};als 0. Also eine leere Menge.

Analog, auch wenn es mir in diesem Zusammenhang komisch vorkam, habe ich durch die letzten Beiträge dann eben eine Zahl eingesetzt.
Ich dachte, Werte wie 1, 2, 3wären int-Konstanten, die ich nicht einfach einem template zuweisen kann. Hätte ich das mal gemacht.

EDIT: Aber alles gut. Ist ja richtig, auf Unstimmigkeiten hinzuweisen.

[zeropage]

Vorneweg, ich habe meinen Kenntnisstand über Matrizen von hier
https://www.grund-wissen.de/mathematik/lineare-algebra-und-analytische-geometrie/matrizen.html
Fortgeschrittenere Sachen kenne ich noch nicht.

Ich habe einen Haufen Rechenoperationen hinzugefügt,

        Mat2D_T operator*(const Mat2D_T& mT) const
	{
		if (columns() != mT.rows())
			throw std::exception("operator*()");

		Mat2D_T nmT(rows(), mT.columns());
		val_t v = {};
		for (size_t c = 0; c < mT.columns(); ++c)
		{
			for (size_t r = 0; r < rows(); ++r)
			{
				v = (*this)(r, 0) * mT(0, c);
				for (size_t l = 1; l < columns(); ++l)
				{
					v += (*this)(r, l) * mT(l, c);
				}
				nmT(r, c) = v;
			}
			v = {};
		}
		return nmT;
	}

	/////neu/////////////////////////////////////////////

	Mat2D_T& operator+=(const val_t& v)
	{
		for (size_t r = 0; r < rows(); ++r)
			for (size_t c = 0; c < columns(); ++c)
				(*this)(r, c) = (*this)(r, c) + v;
	}

	Mat2D_T& operator-=(const val_t& v)
	{
		for (size_t r = 0; r < rows(); ++r)
			for (size_t c = 0; c < columns(); ++c)
				(*this)(r, c) = (*this)(r, c) - v;
	}

	Mat2D_T& operator*=(const val_t& v)
	{
		for (size_t r = 0; r < rows(); ++r)
			for (size_t c = 0; c < columns(); ++c)
				(*this)(r, c) = (*this)(r, c) * v;
	}

	Mat2D_T& operator/=(const val_t& v)
	{
		for (size_t r = 0; r < rows(); ++r)
			for (size_t c = 0; c < columns(); ++c)
				(*this)(r, c) = (*this)(r, c) / v;
	}

	Mat2D_T& operator+=(const Mat2D_T& mT)
	{
		if (rows() != mT.rows() && columns() != mT.columns())
			throw std::exception("operator+=()");

		for (size_t r = 0; r < rows(); ++r)
			for (size_t c = 0; c < columns(); ++c)
				(*this)(r, c) = (*this)(r, c) + mT(r, c);
	}

	Mat2D_T& operator-=(const Mat2D_T& mT)
	{
		if (rows() != mT.rows() && columns() != mT.columns())
			throw std::exception("operator-=()");

		for (size_t r = 0; r < rows(); ++r)
			for (size_t c = 0; c < columns(); ++c)
				(*this)(r, c) = (*this)(r, c) - mT(r, c);
	}

	Mat2D_T operator+(const val_t& v) const
	{
		Mat2D_T nmT(rows(), columns());

		for (size_t r = 0; r < rows(); ++r)
			for (size_t c = 0; c < columns(); ++c)
				nmT(r, c) = (*this)(r, c) + v;

		return nmT;
	}

	Mat2D_T operator-(const val_t& v) const
	{
		Mat2D_T nmT(rows(), columns());

		for (size_t r = 0; r < rows(); ++r)
			for (size_t c = 0; c < columns(); ++c)
				nmT(r, c) = (*this)(r, c) - v;

		return nmT;
	}

	Mat2D_T operator*(const val_t& v) const
	{
		Mat2D_T nmT(rows(), columns());

		for (size_t r = 0; r < rows(); ++r)
			for (size_t c = 0; c < columns(); ++c)
				nmT(r, c) = (*this)(r, c) * v;

		return nmT;
	}

	Mat2D_T operator/(const val_t& v) const
	{
		Mat2D_T nmT(rows(), columns());

		for (size_t r = 0; r < rows(); ++r)
			for (size_t c = 0; c < columns(); ++c)
				nmT(r, c) = (*this)(r, c) / v;

		return nmT;
	}

	Mat2D_T operator+(const Mat2D_T& mT) const
	{
		if (rows() != mT.rows() && columns() != mT.columns())
			throw std::exception("operator+()");

		Mat2D_T nmT(rows(), columns());

		for (size_t r = 0; r < rows(); ++r)
			for (size_t c = 0; c < columns(); ++c)
				nmT(r, c) = (*this)(r, c) + mT(r, c);

		return nmT;
	}

	Mat2D_T operator-(const Mat2D_T& mT) const
	{
		if (rows() != mT.rows() && columns() != mT.columns())
			throw std::exception("operator-()");

		Mat2D_T nmT(rows(), columns());

		for (size_t r = 0; r < rows(); ++r)
			for (size_t c = 0; c < columns(); ++c)
				nmT(r, c) = (*this)(r, c) - mT(r, c);

		return nmT;
	}

was ich mich nun frage, eine template-Matrix sollte ja alle möglichen Typen aufnehmen können, also auch zB einen std::vector()oder wieder eine Matrix. Nur rechnen kann man damit nur bedingt. Sollte man allso die Rechenoperationen mit einem eigenen template-Typ auslagern? Oder führt das zu weit?

EDIT: Und ist es sinnvoll, aus Gründen der Vollständigkeit auch alle undefinierten Fälle zu berücksichtigen wie zB

        Mat2D_T operator/(const Mat2D_T& mT) const
	{
		throw std::exception("operator/()");
	}

?

[zeropage]

PS: in den obigen frisch gezeigten Funktionen muss bei der Prüfung auf selbe Dimension ein orstatt ein andhin. Wie konnte das denn (bei mir) nur durchschlüpfen?

[*john 0]

@zeropage sagte in Entwicklung einer einfachen Matrizenklasse:

was ich mich nun frage, eine template-Matrix sollte ja alle möglichen Typen aufnehmen können, also auch zB einen std::vector()oder wieder eine Matrix. Nur rechnen kann man damit nur bedingt. Sollte man allso die Rechenoperationen mit einem eigenen template-Typ auslagern? Oder führt das zu weit?

Eine Matrix ist nur sinnvoll für einen Vektorraum definiert. Der Skalerenkörper des Vektorraums muss entweder die Körper- oder die Schiefkörperaxiome erfüllen. Die Matrizen, die Du bisher betrachtet hast, haben als Skalarenkörper, den Körper der reellen Zahlen. Ein anderer Körper ist der Körper der komplexen Zahlen. Komplexe Matrizen haben noch einige zusätzliche Eigenschaften, die reelle Matrizen nicht haben können z.B. hermitesche Matrix.

[zeropage]

Aha ok, danke. Schon etwas komplizierter.

[zeropage]

Habe noch einige Funktionen hinzugefügt:

        Mat_T transpose() const
	{
		Mat_T mT(columns(), rows());

		for (size_t r = 0; r < rows(); ++r)
			for (size_t c = 0; c < columns(); ++c)
				mT(c, r) = (*this)(r, c);

		return mT;
	}

	Mat_T& negate() 
	{
		return (*this) * -1;
	}
	
	Mat_T negate() const 
	{ 
		Mat_T<val_t> mT = (*this);
		mT = mT * -1;
		return mT;
	}

	Mat_T& invert()
	{
		if (!invertMatrix((*this)))
			throw std::exception("invert(): not possible");

		return (*this);
	}

	Mat_T invert() const
	{
		Mat_T<val_t> mT = (*this);

		if (!invertMatrix(mT))
			throw std::exception("invert(): not possible");

		return mT;
	}

	template<typename cast_t>
	Mat_T<cast_t> cast() const
	{
		Mat_T<cast_t> mT(rows(), columns());

		for (size_t r = 0; r < rows(); ++r)
			for (size_t c = 0; c < columns(); ++c)
				mT(c, r) = static_cast<cast_t>((*this)(r, c));

		return mT;
	}

        bool swapRow(Mat_T<val_t>& mT, size_t rowA, size_t rowB)
	{
		if (rowA >= mT.rows() || rowB >= mT.rows())
			return false;

		for (size_t c = 0; c < mT.columns(); ++c)
		{
			const val_t v = mT[rowA][c];
			mT[rowA][c] = mT[rowB][c];
			mT[rowB][c] = v;
		}

		return true;
	}

	/*Invertiert eine NxN Mat mit Hilfe des Gauß-Jordan-Algorithmus.
	return: false, falls die Mat nicht invertierbar ist.*/
	bool invertMatrix(Mat_T<val_t>& imT)
	{
		if (imT.rows() != imT.columns())
			return false;

		size_t N = imT.rows();
		// Eine Nx2N Mat für den Gauß-Jordan-Algorithmus aufbauen
		Mat_T<val_t> A(N, 2 * N);
		for (size_t r = 0; r < N; ++r)
		{
			for (size_t c = 0; c < N; ++c)
				A[r][c] = imT[r][c];
			for (size_t c = N; c < 2 * N; ++c)
				A[r][c] = (r == c - N) ? 1.0 : 0.0;
		}
		// Gauß-Algorithmus
		for (size_t k = 0; k < N - 1; ++k)
		{
			// Zeilen vertauschen, falls das Pivotelement eine Null ist
			if (A[k][k] == 0.0)
			{
				for (size_t r = k + 1; r < N; ++r)
				{
					if (A[r][k] != 0.0)
					{
						swapRow(A, k, r);
						break;
					}
					else if (r == N - 1)
						return false; // Es gibt kein Element != 0
				}
			}
			// Einträge unter dem Pivotelement eliminieren
			for (size_t r = k + 1; r < N; ++r)
			{
				const val_t v = A[r][k] / A[k][k];
				for (size_t c = k; c < 2 * N; ++c)
					A[r][c] -= A[k][c] * v;
			}
		}

		// Determinante der Mat berechnen
		val_t det = 1.0;
		for (size_t k = 0; k < N; ++k)
			det *= A[k][k];

		if (det == 0.0)  // Determinante ist =0 -> Mat nicht invertierbar
			return false;

		// Jordan-Teil des Algorithmus durchführen
		for (size_t k = N - 1; k > 0; --k)
		{
			for (int r = k - 1; r >= 0; --r)
			{
				const val_t v = A[r][k] / A[k][k];
				for (size_t c = k; c < 2 * N; ++c)
					A[r][c] -= A[k][c] * v;
			}
		}
		// Einträge in der linker Mat auf 1 normieren und in imT schreiben
		for (size_t r = 0; r < N; ++r)
		{
			const val_t v = A[r][r];
			for (size_t c = N; c < 2 * N; ++c)
				imT[r][c - N] = A[r][c] / v;
		}

		return true;
	}

[zeropage]

Wieso funktioniert das jetzt nicht?

        Mat_T invert() const
	{
		//Mat_T<val_t> mT = (*this);
		Mat_T<val_t> mT = copyMatrix(); //Ob nun so oder so

		if (!invertMatrix(mT))
			throw std::exception("invert(): not possible");

		return mT;
	}

        Mat_T copyMatrix() const
	{
		Mat_T<val_t> copyMat(rows(), columns());

		for (size_t r = 0; r < rows(); ++r)
			for (size_t c = 0; c < columns(); ++c)
				copyMat(r, c) = component(r, c);

		return copyMat;
	}

Fehlermeldung:

Error	C2662	'bool Mat_T<float>::invertMatrix(Mat_T<float> &)': cannot convert 'this' pointer from 'const Mat_T<float>' to 'Mat_T<float> &'	

Was will der converten? Er soll mT nehmen und gut?

[zeropage]

Die Methode, in der wiederum diese Methode aufgerufen wird darf nicht constsein. Manche Fragen scheinen sich erst zu lösen, wenn man tatsächlich auch fragt.

[Schlangenmensch]

Ist das eine freie Funktion, oder eine Member Funktion?

bool invertMatrix(Mat_T<val_t>& imT)

[zeropage]

Diese ist noch eine Memberfunktion.

[Schlangenmensch]

Warum ist die nicht const? Du kannst aus const Member Funktionen nur andere Member Funktionen aufrufen, die auch const sind.

[zeropage]

Ok, da scheine ich einen Denkfehler gehabt zu haben.

Was ist denn jetzt günstiger, nachvollziehbarer, besser? Diese, also bool swapRow() und bool invertMatrix(), als Memberfunktion lassen und const machen oder als freie Funktionen?

[Th69]

Alle Funktionen, welche nicht direkt auf private (oder protected) Klassenmember zugreifen, sollten als freie Funktionen definiert sein (am besten auch innerhalb eines zugehörigen Namensbereichs).

Aber wieso hast du überhaupt zwei Funktionen: invertMatrix() und invert()?

[zeropage]

Der Gedanke war, das ich invert() jeweils als const und als Referenz ausführen kann, ohne bool swapRow()und bool invertMatrix() groß zu ändern. Das ich sie zB problemlos als freie Funktionen in einem Namensraum auslagern kann.

[Th69]

Wäre es dann nicht besser, wenn invert die Memberfunktion ist, welche 'inplace' die Matrix invertiert und invertMatrix die freie Funktion, welche dann eine Kopie erzeugt und dann darauf invert aufruft (und die invertierte Kopie zurück gibt)?

Ich entwickle die Klassenfunktionen meistens so, daß dort die ganze Logik implementiert ist und die freien Funktionen dann nur 'convenience'-Funktionen darstellen (welche dann die Funktionalität der Memberfunktionen benutzt und vereinfachte Aufrufe darstellt).