[gelöst]#QNAN - Gauss Algorithmus; Kollisionsabfrage

Xethoras · Beitrag von **Xethoras** » 09.04.2009, 18:42

So.. hab mal wieder ein Problem :D
Zugegeben, ich glaube nicht, dass jemand viel Lust hat sich den Code, den ich jetzt posten werde komplett durchzulesen. Aber vllt hat ja doch jemand Langeweile oder jemandem fällt die Lösung sofort auf^^

Kurz nochmal zum Problem selbst: Ich benutze einen Gaussalgorithmus, bei dem ich "teilen durch NULL" eigentlich durch vorherige Abfragen vermeide. Er dient der Kollisionsabfrage. Diese funktioniert leider nicht, deshalb hab ich mir die betroffenen Daten im Log ausgeben lassen. Ich weiß nicht wirklich woher die ganzen NANs kommen, da ich weder mit undendlich hantiere, noch negative Zahlen "wurzele". Dass ich "durch Null teilen" versuche zu vermeiden hab ich ja auch schon erwähnt :D

Ich bedanke mich schon einmal im vorraus bei denen, die sich des Problems annehmen werden^^ (sofern es denn welche geben wird^^)

Der Gauss Algorithmus:

Code: Alles auswählen

#include "P12_Math.h"
#include "stdio.h"

void P12G_ChangeLines(int dim,float *mat,int z1,int z2)
{
	float *p=new float[dim+1];
	for(int i=0;i<dim+1;i++)
	{
		p[i]=mat[(z1-1)*(dim+1)+i];
		mat[(z1-1)*(dim+1)+i]=mat[(z2-1)*(dim+1)+i];
		mat[(z2-1)*(dim+1)+i]=p[i];
	}
	delete[] p;
}
void P12G_AddLineToLine(int dim,float *mat,int Src,int Dst)
{
	for(int i=0;i<dim+1;i++)
		mat[(Dst-1)*(dim+1)+i]+=mat[(Src-1)*(dim+1)+i];
}
void P12G_MultLine(int dim,float *mat,int z,float value)
{
	for(int i=0;i<dim+1;i++)
		mat[(z-1)*(dim+1)+i]*=value;
}
bool P12G_IfNull(float x)
{
	if(x<0.000001f&&x>-0.000001f)
		return true;
	else return false;
}
float * P12_GAUSS(int dim,float * mat)	
{
//Untere Dreiecksform
	for(int i=0;i<dim-1;i++)
	{
		for(int j=0;j<dim-1-i;j++)
		{
			if(P12G_IfNull(mat[(dim-j-1)*(dim+1)+i]))
				continue;
			if(P12G_IfNull(mat[(dim-j-2)*(dim+1)+i]))
			{
				P12G_ChangeLines(dim,mat,dim-j,dim-j-1);
				continue;
			}
			P12G_MultLine(dim,mat,dim-j,-1/mat[(dim-j-1)*(dim+1)+i]);
			P12G_MultLine(dim,mat,dim-j-1,1/mat[(dim-j-2)*(dim+1)+i]);
			P12G_AddLineToLine(dim,mat,dim-j-1,dim-j);			
		}
	}
//Obere Dreiecksform
	for(int i=0;i<dim-1;i++)
	{
		if(P12G_IfNull(mat[(dim-i-1)*(dim+1)+dim-i-1]))
				return NULL;
		P12G_MultLine(dim,mat,dim-i,1/mat[(dim-i-1)*(dim+1)+dim-i-1]);
		for(int j=0;j<dim-1-i;j++)
		{
			if(P12G_IfNull(mat[j*(dim+1)+dim-i-1]))
				continue;
			P12G_MultLine(dim,mat,j+1,-1/mat[j*(dim+1)+dim-i-1]);
			P12G_AddLineToLine(dim,mat,dim-i,j+1);
		}
	}
//Ergebnis zurückgeben...
	if(P12G_IfNull(mat[0]))
		return NULL;
	P12G_MultLine(dim,mat,1,1/mat[0]);
	static float *p=NULL;	//Ergebnisarray
	if(p) delete[] p;			//Speicherlecks bitte erst nach dem Beenden des Programms^^
	p=new float[dim];
	for(int i=0;i<dim;i++)
		p[i]=mat[(i+1)*(dim+1)-1];
	return p;
}

Der Code, wo er angewendet wird:

Code: Alles auswählen

		void * buffer;
		objects[i]->getMesh()->getMesh()->LockIndexBuffer(D3DLOCK_READONLY,&buffer);
		for(unsigned int j=0;j<objects[i]->getMesh()->getMesh()->GetNumFaces();j++)
		{
			float mat[12]=	{	reinterpret_cast<D3DXVECTOR3 *>(buffer)[j*3+1].x-reinterpret_cast<D3DXVECTOR3 *>(buffer)[j*3].x, reinterpret_cast<D3DXVECTOR3 *>(buffer)[j*3+2].x-reinterpret_cast<D3DXVECTOR3 *>(buffer)[j*3].x, oldPosX-newPosX, oldPosX-reinterpret_cast<D3DXVECTOR3 *>(buffer)[j*3].x,
								reinterpret_cast<D3DXVECTOR3 *>(buffer)[j*3+1].y-reinterpret_cast<D3DXVECTOR3 *>(buffer)[j*3].y, reinterpret_cast<D3DXVECTOR3 *>(buffer)[j*3+2].y-reinterpret_cast<D3DXVECTOR3 *>(buffer)[j*3].y, oldPosY-newPosY, oldPosY-reinterpret_cast<D3DXVECTOR3 *>(buffer)[j*3].y,
								reinterpret_cast<D3DXVECTOR3 *>(buffer)[j*3+1].z-reinterpret_cast<D3DXVECTOR3 *>(buffer)[j*3].z, reinterpret_cast<D3DXVECTOR3 *>(buffer)[j*3+2].z-reinterpret_cast<D3DXVECTOR3 *>(buffer)[j*3].z, oldPosZ-newPosZ, oldPosZ-reinterpret_cast<D3DXVECTOR3 *>(buffer)[j*3].z
							};
			if(!P12_GAUSS(3,mat)) continue;
//DEBUGMaßnahme^^
			*(iData->getLogFile())<<"\nAusgangspositionen:"<<oldPosX<<"\t"<<oldPosY<<"\t"<<oldPosZ<<"\t"<<newPosX<<"\t"<<newPosY<<"\t"<<newPosZ<<"\n";
			for(int i=0;i<12;i++)
			{
				*(iData->getLogFile())<<"\t"<<mat[i];
			}
//Bis hier
			if	(	mat[3]>=0.0f
				&&	mat[7]>=0.0f
				&&	mat[3]+mat[7]<=1.0f
				&&	mat[11]>=0.0f&&mat[11]<=1.0f
				)
			{
				state=true;
                       }
         }

So, jetzt das eigentlich deprimierende, ein Ausschnitt aus meiner Logdatei:

Code: Alles auswählen

Ausgangspositionen:-3.01092	1.10994	0.224866	-3.001	1.10979	0.226155
	1	-0	0	2.87457e-032	-0	1	-0	-6.5157e-016	1.55107e-035	-4.95424e-036	1	7478.59
Ausgangspositionen:-3.01092	1.10994	0.224866	-3.001	1.10979	0.226155
	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	-7.41455e-038	-3.06298e-037	1	303.625
Ausgangspositionen:-3.01092	1.10994	0.224866	-3.001	1.10979	0.226155
	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	0	1	0	1.31087e-034	3.00367e-037	2.26223e-037	1	303.625
Ausgangspositionen:-3.001	1.10979	0.226155	-2.99108	1.10964	0.227443
	1	-0	0	2.8818e-032	-0	1	-0	-6.53108e-016	1.55481e-035	-4.96621e-036	1	7495.65
Ausgangspositionen:-3.001	1.10979	0.226155	-2.99108	1.10964	0.227443
	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	-7.41455e-038	-3.06298e-037	1	302.625
Ausgangspositionen:-3.001	1.10979	0.226155	-2.99108	1.10964	0.227443
	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	0	1	0	1.31087e-034	3.00367e-037	2.26223e-037	1	302.625
Ausgangspositionen:-2.99108	1.10964	0.227443	-2.98117	1.10949	0.228732
	1	-0	0	2.87457e-032	-0	1	-0	-6.5157e-016	1.55107e-035	-4.95424e-036	1	7476.59
Ausgangspositionen:-2.99108	1.10964	0.227443	-2.98117	1.10949	0.228732
	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	-7.41455e-038	-3.06298e-037	1	301.625
Ausgangspositionen:-2.99108	1.10964	0.227443	-2.98117	1.10949	0.228732
	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	0	1	0	1.31087e-034	3.00367e-037	2.26223e-037	1	301.625
Ausgangspositionen:-2.98117	1.10949	0.228732	-2.97125	1.10935	0.23002
	1	-0	0	2.87698e-032	-0	1	-0	-6.52082e-016	1.55231e-035	-4.95822e-036	1	7481.6
Ausgangspositionen:-2.98117	1.10949	0.228732	-2.97125	1.10935	0.23002
	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	-7.41455e-038	-3.06298e-037	1	300.625
Ausgangspositionen:-2.98117	1.10949	0.228732	-2.97125	1.10935	0.23002
	-1.#QNAN	-1.#QNAN	-1.#QNAN	-1.#QNAN	0	1	0	1.31087e-034	3.00367e-037	2.26223e-037	1	300.625

Aramis · Beitrag von **Aramis** » 09.04.2009, 18:58

Hoi,

du kannst mit [ code=cpp ] das Syntaxhighlighting für C++-Code aktivieren. Ich denke damit dürfte es um einiges übersichtlicher werden.

Code: Alles auswählen

if(p) delete[] p;

Das ist übrigens nicht nötig. Ein delete[] auf einem Nullzeiger ist definiert.

Bist du dir sicher dass dein Epsilon beim Vergleich gegen 0 wirklich sinnvoll gewählt ist? Experimentier doch mal mit anderen Werten herum. Wenn du viel rechnest, kann die Fließkommaungenauigkeit größer sein als man denkt.

eXile · Beitrag von **eXile** » 09.04.2009, 19:30

Der Code ist mir zu wurstig, um ganz durchzusteigen. Auch hier gilt mit "cogito ergo sum": Divide et impera! Ist der Algorithmus mit als Sicherheit anzunehmender Wahrscheinlichkeit bis zur oberen Dreiecksform korrekt (d.h. mit mehreren Testdaten das Ergebnis überprüfen)?

Zur weiteren Recherche kann ich dir auch das Stichwort "Gaußverfahren mit Pivotierung" ans Herz legen. Das Gaußverfahren ist übrigens als solches numerisch instabil, und wenn es nicht unbedingt schnell gehen muss ist eine QR-Zerlegung mit Givens-Rotation das Mittel der Wahl.

Xethoras · Beitrag von **Xethoras** » 09.04.2009, 19:46

@Aramis: Erstmal vielen Dank für die Tipps, die modifizierung des Nullvergleichs hat bisher noch nicht das (von mir) gewünschte Resultat gebracht^^

@exile: Ich hab den Code durchaus getestet, ja. Allerdings mit i-welchen Zahlen, die man auch in der Schule oder so bekommt, da hat es einwandfrei funktioniert.Pivotierung, QR-Zerlegung, sowie Givens-Rotation höhre ich zugegebenermaßen das erste Mal, aber ich werde mich informieren^^ (wobei es eigentlich schon schnell gehen sollte^^)

eXile · Beitrag von **eXile** » 09.04.2009, 22:03

Xethoras hat geschrieben:@exile: Ich hab den Code durchaus getestet, ja. Allerdings mit i-welchen Zahlen, die man auch in der Schule oder so bekommt, da hat es einwandfrei funktioniert.Pivotierung, QR-Zerlegung, sowie Givens-Rotation höhre ich zugegebenermaßen das erste Mal, aber ich werde mich informieren^^ (wobei es eigentlich schon schnell gehen sollte^^)

Naja, wenn das Gleichungssystem singulär ist, dann ist es natürlich nicht durch den Gauß-Algorithmus lösbar. Aber in der Regel müsste das ja mit einem Test auf Division durch Null erkennbar sein.

Hier ne kurze Übersicht was es so gibt: http://de.wikipedia.org/wiki/Liste_nume ... ngssysteme

Aber im Grunde reicht ja ein Gauß-Verfahren schon aus. Bei Gelegenheit schaue ich mal näher in deinen Code ...

Xethoras · Beitrag von **Xethoras** » 10.04.2009, 21:31

Ihr braucht euch keine weiteren Gedanken um das Problem zu machen, es lag nicht am Gaussalogrithmus, sondern an den völlig falschen Daten, die ich übergeben habe^^ Vielen Dank für euer Bemühen mir zu helfen^^

ZFX + Developia

[gelöst]#QNAN - Gauss Algorithmus; Kollisionsabfrage

[gelöst]#QNAN - Gauss Algorithmus; Kollisionsabfrage

Re: #QNAN - Gauss Algorithmus; Kollisionsabfrage

Re: #QNAN - Gauss Algorithmus; Kollisionsabfrage

Re: #QNAN - Gauss Algorithmus; Kollisionsabfrage

Re: #QNAN - Gauss Algorithmus; Kollisionsabfrage

Re: [gelöst]#QNAN - Gauss Algorithmus; Kollisionsabfrage