Wie ein Lineares Gleichungssytem lösen?



  • Omti schrieb:

    So sieht das Gleichungssystem in etwa aus.
    http://img368.imageshack.us/img368/9313/matrixkubsplinewo9.th.png

    Allerdings sollte es wesentlich größer werden können. Ich hatte eigentlich vor das es mindestens eine 60x61 matrix sein sollte.

    Wenn ich den Gaus richtig verstanden habe, dann subtrahiert man immer die erste Zeile von der zweiten, danach die zweite von der dritten, usw bis die matrix in dreiecks Form ist. Muss ich dann die erste Zeile mit einer vertauschen wo halt kein 0 an erster Stelle steht?

    Danke,
    Omti

    muss es selbst implementiert sein, oder kann es schon was fertiges in form einer bibliothek sein? 🙂

    edit:

    wenns letzeres ist und C sein soll, kann ich die gnu scientific library empfehlen 🙂

    würde dann zum beispiel so aussehen:

    #include <gsl/gsl_linalg.h>
    
    double _A[60*60];
    double _Y[60];
    double a = 0.0;
    int s = 0;
    int i = 0;
    gsl_matrix_view A;
    gsl_vector_view Y;
    gsl_vector *_X = NULL;
    gsl_permutation *p = NULL;
    
    // _A und _Y initialisieren
    
    A = gsl_matrix_view_array(_A, 60, 60);
    Y = gsl_vector_view_array(_Y, 60);
    *_X = gsl_vector_alloc(60);
    *p = gsl_permutation_alloc(60);
    
    gsl_linalg_LU_decomp(&A, p, &s);
    gsl_linalg_LU_solve(&A.matrix, p, &Y.vector, _X);
    
    // mit gsl_vector_get lösungen erhalten
    for (i = 0; i < 60; i++)
    {
    	a = gsl_vector_get(_X, i);
    	// sonstwas damit machen
    }
    
    gsl_permutation_free(p);
    gsl_vector_free(_X);
    


  • Im (wirklich guten) Wikipedia-Artikel findest du auch einen passenden passenden Pseudocode.

    Wenn es auch was fertiges sein darf, kannst du boost.uBLAS benutzen.

    Omti schrieb:

    Das Problem ist, wir haben das gaussche Lösungsverfahren im LK einfach übergangen, weil wir nen Taschenrechner haben der das kann... 😞

    uff. Aus dem Grund bin ich dagegen solche Taschenrechner in Schulen einzusetzen. An deiner Stelle würde ich mir den Wikipedia-Artikel wirklich gut durchlesen.



  • Walli schrieb:

    Gerne macht man sich auch die Struktur und sonstige Eigenschaften der Matrix zu Nutze (z.B. Bandmatrix, Dreiecksmatrix, Symmetrie, Definitheit, Dünnbesetztheit ...). Kannst Du was näheres zu der Struktur sagen? Von welchen Dimensionen n reden wir überhaupt maximal?

    Bei kubischen Splines ergibt sich IIRC immer eine Tridiagonalmatrix, wenn man es geschickt anstellt, ist sie sogar symmetrisch. Das steht jedenfalls bei Sedgewick irgendwo.



  • Jan schrieb:

    Walli schrieb:

    Cramer ist im Regelfall viel zu ineffizient.

    So wie ich das sehe, geht es hier nicht vorrangig um Effizienz. Der Vorteil an Cramer ist, dass es relativ einfach umzusetzen ist.

    Naja, ich weiß nicht. Bei allen Fällen, wo ich mir den Einsatz von bikubischen Splines vorstellen kann, geht es schon irgendwie um Effizienz, auch wenn man die Koeffizienten nur einmal berechnen braucht falls die Stützpunkte konstant bleiben. Und der Gauß ist auch nicht viel schwerer zu implementieren.

    Bashar schrieb:

    Walli schrieb:

    Gerne macht man sich auch die Struktur und sonstige Eigenschaften der Matrix zu Nutze (z.B. Bandmatrix, Dreiecksmatrix, Symmetrie, Definitheit, Dünnbesetztheit ...). Kannst Du was näheres zu der Struktur sagen? Von welchen Dimensionen n reden wir überhaupt maximal?

    Bei kubischen Splines ergibt sich IIRC immer eine Tridiagonalmatrix, wenn man es geschickt anstellt, ist sie sogar symmetrisch. Das steht jedenfalls bei Sedgewick irgendwo.

    Ich weiß grad auch nicht genau, ob es unbedingt tridiagonal ist, aber zumindest Bandstruktur müsste sie haben. Die kann man ja relativ einfach in O(n) lösen.



  • Wenn es sich um einen Endomorphismus handelt kann man die Matrix ja eh immer auf die Jordansche Normalform bringen, also als Bidiagonal(abbildungs)matrix darstellen.

    Aber all diese Theorie hilft dir nichts, denn ohne Gauß-Verfahren kommst du nicht sehr weit. Aber du hast ja eh nur eine 60x60 Matrix da reicht auch einfach ein unoptimierter Gauß-Algorithmus, und im Moment bereitet der dir ja schon genug Probleme.

    Hast du mal deinen Lehrer darauf angesprochen ich bin mir sicher dass er dir offene Fragen dazu beantworten kann 🙂



  • Ich hab jetzt mein Progrämmchen so weitergeschrieben dass es die Zeilen vertauscht. Das funktioniert soweit auch ganz gut, es gibt nur noch ein Problem. Die oberste Zeile, in der Nur ein y wert und eine koeffiziente (so heißt da doch :D) stehenrutscht nach ganz unten durch.

    http://img208.imageshack.us/my.php?image=kubumsorthd8.png

    Der Code sieht so aus:

    for(int j=0;j<amt;j++)  //amt ist die Größe der Matrix
    	{
        if(kub [j][j]==0) //kub ist die Matrix
        	{
            for(int k=j;k<=amt;k++)
            	{
                if(kub [k][j]!=0)
                	{
                    for(int l=0;l<=amt;l++)
                    	{
                        tausch=0;
                        tausch=kub[j][l];
                        kub[j][l]=kub[k][l];
                        kub[k][l]=tausch;
                        }
                    }
                }
    		}
    

    Jetzt weiß ich zwar dank wikipedia in etwa wie es weitergeht... (ich habe das nach dem 10ten mal durchlesen endlich halbwegs verstanden ^^)
    Jetzt weiß ich halt leider nur nicht wie ich die Matrix richtig sortiert kriege.

    Naja... Wenn ihr mir dabei helfen könntet, sollte ich es auf die Reihe kriegen.
    Danke,
    Omti



  • So als Tipp Abstraktion erleichtert dir das Leben ich würde an deiner Stelle eine Funktion schreiben, welche zwei Zeilen vertauscht, das erhöht die Übersichtlichkeit und du verlierst nicht den Blick für das Wesentliche.



  • nunja ich habe jetzt meinen Löser fertiggestellt 😃

    for(int j=0;j<amt;j++) //amt ist die größe der matrix
    	{
         for(int k=j+1;k<amt;k++)
         	{
            double x=kub[k][j]/kub[j][j]; //kub ist die Matrix
            for(int l=0;l<=amt;l++)
            	{
                 tausch=kub[k][l]-(kub[j][l]*x);
                 if ((fabs(tausch) < 0.00001) && (l==j))
                    {
    
                     tausch = 0.0;
                     }
                 kub[k][l]=tausch;
                }
            }
    
        }
    
    for(int j=amt-1;j>=0;j--)
    	{
         for(int k=j-1;k>=0;k--)
         	{
            double y=kub[k][j]/kub[j][j];
            for(int l=0;l<=amt;l++)
            	{
                 tausch=kub[k][l]-(kub[j][l]*y);
                 if ((fabs(tausch) < 0.00001) && (l==j))
                    {
    
                     tausch = 0.0;
                     }
                 kub[k][l]=tausch;
                }
            }
    
        }
     for(int j=0;j<amt;j++)  //amt ist die Dimension der Matrix
    	{
        tausch=kub[j][j];
        kub[j][j]= kub[j][j]/kub[j][j];
        kub[j][amt]=kub[j][amt]/tausch;
        }
    

    Das ganze ist jetzt für ein bereits pivotisiertes Gleichungssystem...



  • Bei natuerlichen Randbedingungen fuer Interpolation mit kubischen Splines bekommst Du ein diagonaldominates, tridiagonales Gleichungssystem, welches keine Pivotsuche erfordert. Bei Hermite-Randbedingungen uebrigens auch.



  • Die Transformation in die Jordan-Normalform ist nicht stetig. Daher ist dieses Verfahren für numerische Anwendungen sowieso wertlos.


Anmelden zum Antworten