[C] metodo di Gauss-Jordan

[kobeilprofeta]

ho fatto questa funzione per scambiare due righe tra loro:

float sca (int x, int y)
{
      for (i=1;i<=c;i++)
      {
      a[0][i]=a[x][i];
      a[x][i]=a[y][i];
      a[x][i]=a[0][i];  
      }
}

...dove a[x][y] è l'elemento di posto x della riga y.

Ovviamente non funziona. se richiamo ad esempio sca (2,4) non accade nulla.
ho pensato di dover usare i puntatori... ma non so come: ho provato qualcosa del genere ma non va:

float sca (int *x, int *y)
{
      for (i=1;i<=c;i++)
      {
      a[0][i]=a[*x][i];
      a[*x][i]=a[*y][i];
      a[*x][i]=a[0][i];  
      }
}

grazie.

[Super Squirrel]

#include<iostream>

void mostra_matrice(int v[10][10], int rig, int col)
{
    for (int i = 0; i < rig; i++)
    {
        std::cout << std::endl;
        for (int j = 0; j < col; j++)
        {
            std::cout << v[i][j] << "\t";
        }
    }
}

void inserisci_matrice(int v[10][10], int rig, int col)
{
    for (int i = 0; i < rig; i++)
    {
        for (int j = 0; j < col; j++)
        {
            std::cin >> v[i][j];
        }
    }
}

void scambia_righe(int v[10][10], int riga_1, int riga_2, int col)
{
    int p;
    for(int i = 0; i < col; i++)
    {
        p = v[riga_1][i];
        v[riga_1][i] = v[riga_2][i];
        v[riga_2][i] = p;
    }
}

int main()
{
    int v[10][10], rig, col, riga_1, riga_2;
    std::cin >> rig >> col;
    inserisci_matrice(v, rig, col);
    mostra_matrice(v, rig, col);
    std::cout << std::endl;
    std::cin >> riga_1 >> riga_2;
    scambia_righe(v, riga_1, riga_2, col);
    mostra_matrice(v, rig, col);
}

rig, col, riga_1, riga_2 sono rispettivamente il numero di righe e colonne della matrice e gli indici delle due righe da scambiare di posto.
dall'esempio dovresti riuscire a comprendere come andrebbe fatta la funzione. se ci sono dubbi fammi sapere.

[vict85]

Più che non funzionare, non dovrebbe neanche compilare dato che il tuo valore di ritorno è un float, ma tu non hai scritto alcun return. L'errore del tuo codice, a parte il resto, è che dovi scrivere

a[y][i]

invece di

a[x][i]

nella terza riga.

Rimane comunque un codice pessimo su moltissimi aspetti. In particolare dovresti usare meno le variabili globali. Le variabili locali rendono il codice più leggibile e permettono al compilatore di fare alcune ottimizzazioni che altrimenti potrebbe non fare. Prendi per esempio quel ciclo. Tu stai scrivendo sulla prima riga della matrice dei valori che non ti servono al di fuori di quella funzione, ma il compilatore non può saperlo. Se tu scrivi:

float const temp = a[x][i];
a[x][i] = a[y][i];
a[x][i] = temp;

stai comunicando al compilatore esattamente quello che vuoi. Nello stesso modo, usare una i globale per i cicli è una pessima idea (oltre ad essere obbligatorio se un giorno vorrai parallelizzare il codice). È molto meglio definire la i all'interno del for stesso, ovvero così:

for(int i = 0; i < LimitValue; ++i) 
{ /* */ } 

seppur per farlo tu abbia bisogno di usare il C++ o il C99 (ma trovo insensato usare il C90 anche se molti lo fanno). Con un po' più codice puoi fare la stessa cosa anche in c90.
La matrice e la dimensione, eventualmente in una struct possono invece essere tranquillamente globali. Ma userei dei nomi più espressivi per le variabili, specialmente quelle globali.

Ti faccio inoltre notare che far partire gli array da 0 è fortemente consigliato anche se stai cercando di implementare un algoritmo scritto per matlab o fortran.

Nota infine che il tuo codice con i puntatori è un codice MOLTO pericoloso. Stai praticamente scrivendo su posizioni semi-casuali della memoria. Ti suggerisco di non provare più ad usare i puntatori a caso senza prima aver capito come funzionano. Non te lo correggo neanche perché passare x e y come puntatori è inutile e usa potenzialmente il doppio della memoria (su un computer a 64-bit i puntatori sono da 64-bit mentre gli int sono da 32).

Cosa dovrebbe fare il codice completo?

[kobeilprofeta]

Intanto grazie per le risposte.
Premetto che sto praticamente imparando da solo a programmare in c, quindi faccio ancora molti errori stupidi (magari anche perchè molte cose ancora non le conosco).

I puntatori in effetti non li ho ancora capiti, e a breve dovró cercare di "studiarli"... inoltre sulle struct sono ancora a zero.

Il codice completo dovrebbe elaborare una matrice con il metodo di Gauss-Jordan. So che esistono molti codici in c già fatti in rete, ma ho preferito provare a fare il programma senza "sbirciare" da una già pronto. Al momento ho fatto solo una funzione che divide per il coefficiente direttivo. Ci sto lavorando. Se mi dovessi bloccare ancora, chiedo scusa in anticipo per il nuovo disturbo.

[kobeilprofeta]

"Sergio":[quote="kobeilprofeta"]Premetto che sto praticamente imparando da solo a programmare in c, quindi faccio ancora molti errori stupidi (magari anche perchè molte cose ancora non le conosco).

Consiglio: vai per gradi.
Non metterti a scrivere subito programmi relativamente complessi. L'unico risultato sicuro è che acquisiresti cattive abitudini.
Procurati un'ottima guida (Kernighan & Ritchie è difficile da battere) e vai per gradi.[/quote]

Ciao. Grazie per il consiglio... ma questo non è uno dei primi programmi che faccio. Ne ho già fatti molti, solo che erano decisamente più semplici (o almeno, meno lunghi).


comunque:
questa è la funzione per scegliere il pivot: (p[] è un array che mi rappresenta la riga relativa al pivot della colonna in questione, l' ho dichiarata nel main)

int piv (x)
{
    int j,h,qw,us;
    for (i=1;i<=c;i++)
    {
    for (j=1;j<=r;j++)
    {
    usa=0;
    for (h=1;h<i;h++)
    {
    if (p[h]==j)
    {
     usa=1;     
    }
    }
    if ((usa=0)&&(a[j][i]!=0)) {qw=j;}
    }
    }
    return(qw);
}

L'idea sarebbe quella di richiamare la funzione più o meno così:
p=piv(i);
e in questo modo troverei il pivot di ogni colonna...

...avendo aggiornato lo scambio così (e ora sembra funzionare):

float sca (int x, int y)
{
      int j;
      for (j=1;j<=c;j++)
      {
          if (p[j]==x)  {p[j]=y;}
          else if (p[j]==x)  {p[j]=x;}  
      }
      for (i=1;i<=c;i++)
      {
      a[0][i]=a[x][i];
      a[x][i]=a[y][i];
      a[y][i]=a[0][i];
      }
}

avevo pensato di chiamare nel main una cosa del genere:

for (i=0;i<=c;i++)
{
p[i]=piv(i);
sca(i,piv(i));
}

ps:
@vict ma non solo
[ot]se adesso sto scrivendo questo non pensiate che non ho letto (o che non volgio seguire) ciò che mi avete scritto prima... è solo che mi serve tempo per capire e provare. (mi viene in mente il fatto della dichiarazione delle variabili)[/ot]

[vict85]

Ti suggerisco di riguardarti l'algoritmo. Il calcolo del pivoting va fatto dopo che aggiorni la matrice.

[apatriarca]

@kobeilprofeta: In questo codice vedo esattamente il problema con le variabili globali che ho descritto nell'altra discussione. E non l'ho fatto neanche apposta..
Nel main usi infatti la variabile globale i per il ciclo..

for (i=0;i<=c;i++) {
    p[i]=piv(i);
    sca(i,piv(i));
}

e la usi di nuovo anche dentro sca..

for (i=1;i<=c;i++) {
    a[0][i]=a[x][i];
    a[x][i]=a[y][i];
    a[y][i]=a[0][i];
}

[vict85]

Comunque, a parte la pura questione di linguaggio, anche il tuo approccio alla programmazione ha dei problemi.

Da quello che ho capito hai deciso di implementare l'algoritmo “senza guide”, ovvero a memoria e ho l'impressione senza neanche un vero schema mentale. Il problema è che fare così è solo una perdita di tempo ed energie. Il tuo codice del pivoting ha poco senso.

Il codice del pivoting comunemente usato è qualcosa come questo:

int GJ_ricercapivot(double M[MatSize][MatSize], int colonna)
{
    double max = abs(M[colonna][colonna]);
    int ris = colonna;
    for(int i = colonna; i < MatSize; ++i)
    {
        double const t = abs(M[i][colonna]);
        if(max < t)
        {
            max = t;
            ris = i;
        }
    }
    return ris;
}

Nota che ieri ho implementato tutto il G-J e quindi sono sicuro che funziona. Qui ho implementato il partial pivoting. Ovviamente il complete pivoting o il root pivoting avrebbero codici differenti, ma aggiungono un po' di complessità per via della permutazione delle colonne. Ovviamente, seppur sia considerato meno stabile numericamente, puoi usare anche qualcosa di semplice come:

int GJ_ricercapivot(double M[MatSize][MatSize], int colonna)
{
    int ris = colonna;
    for(; (ris < MatSize)&&(M[ris][colonna] == 0); ++ris);
    return ris;
}

ovvero che cerca il primo elemento diverso da 0 della colonna.

Il secondo errore è nel codice che chiama quella funzione che dovrebbe essere più qualcosa di questo tipo:

void GJ(double M[MatSize][MatSize], double b[MatSize])
{
    for(int k = 0; k < MatSize-1; ++k)
    {
        int const p = GJ_ricercapivot(M, k);
        if( M[p][k] == 0 )
        {
            printf("Errore, la matrice non è invertibile");
            return;
        }
        if(p != k)
        {
            scambia_righe(M,p,k);
            double const t = b[k];
            b[k] = b[p];
            b[p] = t;
        }

        /* Parte in cui fai materialmente la riduzione delle righe!  */
    }
    /* Risolvi il problema triangolare  */
}

Nota tra l'altro che avrebbe senso implementassi prima il resto e poi solo dopo il pivoting. Infatti puoi assicurarti di testare l'algoritmo su matrici a diagonale dominante in modo tale da sapere che il pivoting non è necessario.