Die voreingestellte Schwelle für die Dimensions-Reduktion kann durch MF_SVDsetEdit verändert werden. Da MF_SVDsetEdit allerdings nicht fiber-sicher ist, sollte diese Funktion nicht verwendet werden, um verschiedene Schwellenwerte für verschiedene Aufrufe von MF_SVsolve usw. einzustellen. Anstatt wiederholt den Standard-Wert zu ändern, gebrauche man vielmehr die "wEdit"-Varianten von MF_SVsolve usw., also MF_SVsolvewEdit, MF_solveBySVDwEdit und MF_safeSolvewEdit.

Die Singulärwerte sind in W in willkürlicher Reihenfolge gespeichert. Für manche Anwendungen (manuelle Inspektion der Singulärwerte, Dimensions-Reduktion usw.) ist es wünschenswert, die Singulärwerte in absteigender Reihenfolge zu sortieren. Hierbei müssen die Matrizen MU und MV entsprechend mit um-geordnet werden. Diese Aufgabe übernimmt die Funktion MF_SVsort.

Wenn Sie mittels MF_SVdecompose berechnete Ergebnisse mit SVD-Routinen anderer Anbieter vergleichen, beachten Sie bitte folgende Punkte:

1. MF_SVdecompose nimmt selbst keine Dimensions-Reduktion (also keinen Ersatz sehr kleiner Singulärwerte durch 0) vor, sondern überlässt diesen Schritt nachfolgenden Funktionen oder ggf. manueller Inspektion durch den Anwender.
2. Singulärwerte können vertauscht werden, wenn gleichzeitig auch die entsprechenden Spalten sowohl von MU als auch von MV vertauscht werden. Gelegentlich wird diese Eigenschaft genutzt, um die Singulärwerte in eine bestimmte (meist fallende) Reihenfolge zu bringen. Da der Vorzug geordneter Singulärwerte aber durch erhöhten Rechenaufwand erkauft wird, verzichtet MF_SVdecompose selbst auf diese Ordnung. Nur, wenn wirklich benötigt, sollte daher nach MF_SVdecompose noch die eben erwähnte Funktion MF_SVsort aufgerufen werden.
3. Das Vorzeichen eines Singulärwertes kann vertauscht werden, wenn gleichzeitig die entsprechende Spalte von MV mit -1 multipliziert wird. MF_SVdecompose nutzt diese Eigenschaft, um alle nicht-verschwindenden Singulärwerte positiv zu machen.

#include <VDstd.h>
#include <MDstd.h>
#include <conio.h>
void SVDtest( void )
{
    unsigned ht=386, len=52,            // beliebige Werte
             maxhtlen, sizeA, sizeV;
    dMatrix MA, MU, MV, MOne, MA2, MUV2, MDiffA, MDiffOne;
    dVector W;
    double err;
 
    maxhtlen = (ht >= len ? ht : len);
    sizeA = ht*len;
    sizeV = len*len;
 
    MA = MD_matrix( ht, len );   // Eingabe-Matrix
    MU = MD_matrix( maxhtlen, len );
    MV = MD_matrix( len, len );
    W = VD_vector( len );
    MOne = MD_matrix( len, len );   // Einheits-Matrix zum Vergleich
    MD_equ1( MOne, len );
 
    MA2 = MD_matrix( maxhtlen, len );
      // wird MU*W*MVT erhalten zum Vergleich mit der Eingabe-Matrix
      // führende Dimension maxhtlen ist für unter-determinierte Matrizen nötig
    MUV2 = MD_matrix( len, len );
      // wird MUT*MU, MVT*MV, MV*MVT erhalten zum Vergleich mit MOne
    MDiffA = MD_matrix( ht, len );  // für Zwischenergebnisse
    MDiffOne = MD_matrix( len, len );  // für Zwischenergebnisse
 
    MD_random( MA, ht, len, 1, -1., +1. );
      // MA mit Zufallszahlen zwischen -1 und +1 füllen
      // alternativ auf irgendeinem anderen Weg eine Test-Matrix erzeugen oder einlesen
    MD_SVdecompose( MU, MV, W, MA, ht, len );
 
      /*  prüfen, ob MA = MU*W*MVT; Auswertung von rechts nach links: */
    MDdia_mulMT( MUV2, W, MV, len, len ); // MUV2 speichert W * MVT
    MD_mulM( MA2, MU, MUV2, maxhtlen, len, len );
    MD_subM( MDiffA, MA, MA2, ht, len );
      // für unter-determinierte Matrizen werden die Zeilen von ht bis maxhtlen-1 ignoriert
    err = VD_absmax( MDiffA[0], sizeA );
    printf( "SVDtest: Test für die SVD-Routine von OptiVec\n\n"
            "In den folgenden Tests in Double-Genauigkeit\n"
            "sind Fehler bis zu einigen Malen 1.e-14 in Ordnung\n"
            "bei Matrizen der Größenordnung 100x100 bis 1000x1000 Elemente:\n\n"
            "max. Fehler der Matrix-Rekonstruktion MA = MU*W*MVT: %lg", err );
 
    /* prüfe, ob MU orthogonal, also MUT*MU = (1): */
    MD_TmulM( MUV2, MU, MU, maxhtlen, len, len );
    MD_subM( MDiffOne, MUV2, MOne, len, len );
    err = VD_absmax( MDiffOne[0], sizeV );
    printf( "\nmax. Fehler von MUT*MU = (1): %lg", err );
 
    /* prüfe, ob MV orthogonal, also MVT*MV = (1): */
    MD_TmulM( MUV2, MV, MV, len, len, len );
    MD_subM( MDiffOne, MUV2, MOne, len, len );
    err = VD_absmax( MDiffOne[0], sizeV );
    printf( "\nmax. Fehler von MVT*MV = (1): %lg", err );
 
    /* prüfe, ob MV auch orthonormal, also MV*MVT = (1): */
    MD_mulMT( MUV2, MV, MV, len, len, len );
    MD_subM( MDiffOne, MUV2, MOne, len, len );
    err = VD_absmax( MDiffOne[0], sizeV );
    printf( "\nmax. Fehler von MV*MVT = (1): %lg", err );
    printf( "\n\nZum Beenden des Tests beliebige Taste drücken!" ); _getch();
 
    M_nfree( 8, MA, MU, MV, MOne, MA2, MUV2, MDiffA, MDiffOne );
    V_free( W );
}