Contatore nel pivoting

[pollon871]

Ciao. Sto cercando di inserire un contatore di tipo int nella subroutine pivoting, che ha il compito di eseguire l'algoritmo di Gauss con la strategia appunto del pivoting. Tale contatore deve registare il numero di scambi effettuati sulla matrice A di partenza. Tuttavia ho provato ad inserire un

int CONT

e quello che mi sputa fuori nella chiamante è sempre e solo

CONT=32767

cosa vorrà mai significare?
Posto il codice della function

#include<math.h>
int pivoting(float **A, float *b, int n, int k, int CONT) /*CONT conta il numero di scambi effettuati sulle righe */
{int i, j, r; float t, max; 
CONT==0;

max=fabs(A[k][k]);
r=k; /*a questo punto, il massimo della colonna k si trova alla riga k-esima */

for(i=k+1;i<n;i++) /*ricerco il max tra gli elementi della colonna kesima successivi ad A[k][k] */
if(fabs(A[i][k])>max) 
	{max=A[i][k];
	r=i; /*se la condizione dell'if è soddisfatta, devo SCAMBIARE L'INDICE */
	};   /*della riga fissata precedentemente con quest'ultimo, i */
if(r!=k) /*se non si tratta della stessa riga, devo scambiare le due righe e aggiorno il contatore*/
{CONT++;
for(j=k;j<n;j++)
	{t=A[r][j];
	A[r][j]=A[k][j];
	A[k][j]=t;
	};
};

t=b[r];
b[r]=b[k];
b[k]=t;
return CONT;
};

Opinioni?

[vict85]

È difficile capire l'errore senza avere il codice che usi per chiamarla (il codice mi sembra piuttosto corretto). Comunque

CONT==0;

sarebbe stata una riga dannosa se fosse stata una assegnazione e semplicemente senza senso così com'è: toglila.

La chiamata corretta dovrebbe essere qualcosa del tipo

Cont = pivoting(M,V, Cols, i, Cont); 

dove ho usato nomi diversi apposta.

L'hai inizializzato nella funzione che lo chiama?

Manca il controllo che il massimo possa essere 0 (tutta la riga uguale a 0). Ma se vuoi puoi penso ignorarlo.

[pollon871]

Nella chiamante:

int pivoting(float **A, float *b, int n, int k, int CONT);
void gausspiv(float **A, float *b, int n);
void backsub(float **A, float *b, int n, float *x);
main()
{int i, j, n, k, Cont, CONT;
float **A, *b, *x; 
..............
.............
gausspiv(A, b, n); /*chiamata della function*/
printf("\nLa matrice triangolare alta ha componenti\n");
for(i=0; i<n; i++) 
	for(j=0; j<n; j++)
	 printf("A[%d][%d]=%f \n", i, j, A[i][j]);
printf("Il vettore termini noti ha componenti\n");
for(i=0;i<n;i++)
	printf("b[%d]=%f \n", i, b[i]);
Cont=pivoting(A, b, n, k, CONT);
printf("Il numero di scambi effettuati sulle righe è %d", Cont);

dove le function sono:

void gausspiv(float **A, float *b, int n)
{int i, j, k, CONT; float m[100][100];
for(k=0;k<n;k++)
{pivoting(A, b, n, k, CONT); /*applico il metodo pivoting alla kesima riga e procedo */
			/*con il metodo di Gauss sulla nuova matrice (e vettore termini noti) */
for(i=k+1;i<n;i++)
	{
	m[i][k]=A[i][k]/A[k][k]; /*creazione di un coefficiente per l'eliminazione */
	for(j=k;j<n;j++)
		{A[i][j]=A[i][j]-m[i][k]*A[k][j]; /*riduzione dei coefficienti */
		};
	b[i]=b[i]-m[i][k]*b[k]; /*riduzione dei termini noti */
	};
};
return;
}

e

#include<math.h>
int pivoting(float **A, float *b, int n, int k, int CONT) /*CONT conta il numero di scambi effettuati sulle righe */
{int i, j, r; float t, max; 

max=fabs(A[k][k]);
r=k; /*a questo punto, il massimo della colonna k si trova alla riga k-esima */

for(i=k+1;i<n;i++) /*ricerco il max tra gli elementi della colonna kesima successivi ad A[k][k] */
if(fabs(A[i][k])>max) 
	{max=A[i][k];
	r=i; /*se la condizione dell'if è soddisfatta, devo SCAMBIARE L'INDICE */
	};   /*della riga fissata precedentemente con quest'ultimo, i */
if(r!=k) /*se non si tratta della stessa riga, devo scambiare le due righe */
{CONT++;
for(j=k;j<n;j++)
	{t=A[r][j];
	A[r][j]=A[k][j];
	A[k][j]=t;
	};
};

t=b[r];
b[r]=b[k];
b[k]=t;
return CONT;
};

A questo punto Cont è un numero negativo, tipo -1489652

[vict85]

Il tuo errore è che non stai comunicando al main il valore della giusta variabile CONT. Dovresti ripassarti il C.

[pollon871]

Perchè scusa? Il contatore CONT della function pivoting si aggiorna ogni qual volta \(\displaystyle r\neq k \) e quindi effettua lo scambio della riga. Poi il pivoting ritorna un intero che è proprio CONT e lo richiamo nel main con cont=pivoting(...)
Cosa c'è di sbagliato? Indendo dire nel ragionamento, perchè è ovvio che il codice è sbagliato visto che mi dà un risultato errato.

[vict85]

main()
{int i, j, n, k, Cont, CONT;
float **A, *b, *x; 
..............
.............

Il valore di ritorno del main è int, e va messo. Comunque qui Cont e CONT possiedono un valore casuale in quanto non sono inizializzati.

gausspiv(A, b, n); /*chiamata della function*/
printf("\nLa matrice triangolare alta ha componenti\n");
for(i=0; i<n; i++) 
   for(j=0; j<n; j++)
    printf("A[%d][%d]=%f \n", i, j, A[i][j]);
printf("Il vettore termini noti ha componenti\n");
for(i=0;i<n;i++)
   printf("b[%d]=%f \n", i, b[i]);

Questo pezzo della funzione main modifica variabili che possiedono il nome Cont e CONT , ma non quelle del main.

Cont=pivoting(A, b, n, k, CONT);
printf("Il numero di scambi effettuati sulle righe è %d", Cont);

Quest'ultimo pezzo assegna a Cont il valore casuale di CONT eventualmente aumentato di 1. Tra l'altro anche k ha un valore casuale perciò quella chiamata a funzioni potrebbe potenzialmente produrti bug ben peggiori di un semplice valore non corretto.

[pollon871]

"vict85":

main()
{int i, j, n, k, Cont, CONT;
float **A, *b, *x; 
..............
.............

Il valore di ritorno del main è int, e va messo. Comunque qui Cont e CONT possiedono un valore casuale in quanto non sono inizializzati.

gausspiv(A, b, n); /*chiamata della function*/
printf("\nLa matrice triangolare alta ha componenti\n");
for(i=0; i<n; i++) 
   for(j=0; j<n; j++)
    printf("A[%d][%d]=%f \n", i, j, A[i][j]);
printf("Il vettore termini noti ha componenti\n");
for(i=0;i<n;i++)
   printf("b[%d]=%f \n", i, b[i]);

Questo pezzo della funzione main modifica variabili che possiedono il nome Cont e CONT , ma non quelle del main.

Cont=pivoting(A, b, n, k, CONT);
printf("Il numero di scambi effettuati sulle righe è %d", Cont);

Quest'ultimo pezzo assegna a Cont il valore casuale di CONT eventualmente aumentato di 1. Tra l'altro anche k ha un valore casuale perciò quella chiamata a funzioni potrebbe potenzialmente produrti bug ben peggiori di un semplice valore non corretto.

Scusa ma non ho capito nulla.
Allora va bene che all'inizio ovviamente CONT e Cont non sono inizializzati, li ho soltanto dichiarati.

Al secondo passo, dici che quel pezzo modifica ed è vero, perchè in quella function viene chiamato il pivoting che effettua lo scambio. Dici che devo metterlo come output anche di gausspiv? Io ci ho provato, togliendo la chiamata

Cont=pivoting()

ma continua a darmi valori casuali.

[pollon871]

Ci ho provato in questo modo:

int gausspiv(float **A, float *b, int n, int CONT)

come function che ritorna CONT;
e poi la chiamata nel main è

Cont=gausspiv(A, b, n, CONT);

[vict85]

No, il problema è che sbagli proprio il modo in cui interpretare il C, insomma quelle variabili non sono variabili globali.

Sotto un codice corretto in cui sono partito dal tuo. Ho aggiunto un po' di cose per rendere il programma sensato e per testare le mie funzioni. Ciò che ho aggiunto, a parte la backward substitution, crea la matrice A casualmente e definisce b a partire da un x fissato per poi fare la sostituzione e vedere se il risultato è corretto (cosa che sembra esserlo). Ho anche modificato il codice della eliminazione di Gauss mentre ho lasciato pressoché identico il codice per il pivoting. Come vedi il valore di Cont non lo prendo da pivoting ma da gausspiv. Ho cambiato tutto in double perché non sei certo stretto di memoria.

Potrebbe non essere al 100% conforme allo standard C89, perché di solito uso standard più nuovi o il C++. Senza dubbio è conforme allo standard del 2011 (perché ho controllato).

Prova a capire come modificare il tuo usando questo come guida.

#include<stdio.h>
#include<stdlib.h>
#include<time.h>
#include<math.h>

void generaA(double **A, int const N);
void matrix_vector_mult(double **M, double *v, double *r, int const N);
int gausspiv(double **A, double *b, int const N);
int pivoting(double **A, double *b, int const N, int const k, int CONT);
void backsub(double **A, double *b, int const N, double *x);
double Check(double *v, double *w, int const N);

int main( )
{
    double **A = NULL, *b = NULL, *x = NULL;
    double *x_c = NULL;
    int i, Cont;
    double x_error;

    int const N = 256; /* setto il numero per non inserirlo da tastiera. */

    /*
    *    Alloco memoria e creo gli oggetti
    */
    A = calloc(N, sizeof(double *));
    if(A == NULL)
    {
        puts("Allocazione fallita.");
        exit(EXIT_FAILURE);
    }
    for(i = 0; i != N; ++i)
    {
        A[i] = calloc(N, sizeof(double));
        if(A == NULL)
        {
            puts("Allocazione fallita.");
            exit(EXIT_FAILURE);
        }
    }
    generaA(A,N);

    x_c = calloc(N, sizeof(double));
    if(x_c == NULL)
    {
        puts("Allocazione fallita.");
        exit(EXIT_FAILURE);
    }
    x_c[0] = 1.;
    for(i = 1; i != N; ++i)
    {
        x_c[i] = x_c[i-1] + 1./(2*N);
    }

    b = calloc(N, sizeof(double));
    if(b == NULL)
    {
        puts("Allocazione fallita.");
        exit(EXIT_FAILURE);
    }
    matrix_vector_mult(A,x_c,b,N);

    x = calloc(N, sizeof(double));
    if(x == NULL)
    {
        puts("Allocazione fallita.");
        exit(EXIT_FAILURE);
    }

    Cont = gausspiv(A, b, N);
    backsub(A,b,N,x);
    x_error = Check(x,x_c,N);

    printf("\nIl numero di scambi effettuati sulle righe e' %d\n", Cont);
    printf("\nRisultato calcolato con un errore di %f\n", x_error);

    free(b);
    free(x);
    for(i = 0; i != N; ++i)
    {
        free(A[i]);
    }
    free(A);
    free(x_c);
}


void generaA(double **A, int const N)
{
    int i;
    srand(time(0));
    for(i = 0; i != N; ++i)
    {
        int j;
        for(j = 0; j != N; ++j)
        {
            int r = 1 + rand()%(2*N);
            A[i][j] = r/(N + 0.) - 1;
        }
    }
}


void matrix_vector_mult(double **M, double *v, double *r, int const N)
{
    int i;
    for(i = 0; i != N; ++i)
    {
        int j;
        r[i] = 0.;
        for(j = 0; j != N; ++j)
        {
            r[i] += M[i][j]*v[j];
        }
    }
}

double Check(double *v, double *w, int const N)
{
    int i;
    double r = 0.;
    for(i = 0; i != N; ++i)
    {
        double t = fabs(v[i] - w[i]);
        r = r < t ? t : r;
    }
    return r;
}

int gausspiv(double **A, double *b, int const N)
{
    int i, j, k;
    int CONT = 0;

    for(k = 0; k != N; ++k)
    {
        CONT = pivoting(A,b,N,k,CONT);
        if(A[k][k] == 0.)
        {
            puts("Elemento di pivot uguale a 0");
            exit(EXIT_FAILURE);
        }
        for(i = k+1; i != N; ++i)
        {
            double const t = A[i][k] / A[k][k];
            A[i][k] = 0;
            for(j = k+1; j != N; ++j)
            {
                A[i][j] -= A[k][j] * t;
            }
            b[i] -= b[k]*t;
        }
    }
    return CONT;
}



int pivoting(double **A, double *b, int const N, int const k, int CONT)
{
    int i, j, r;
    double t, max;

    max = fabs(A[k][k]);
    r=k;

    for(i=k+1; i<N; i++)
        if(fabs(A[i][k])>max)
        {
            max=A[i][k];
            r=i;
        };
    if(r!=k)
    {
        ++CONT;
        for(j=k; j<N; j++)
        {
            t=A[r][j];
            A[r][j]=A[k][j];
            A[k][j]=t;
        };
    };

    t=b[r];
    b[r]=b[k];
    b[k]=t;

    return CONT;
}

void backsub(double **A, double *b, int const N, double *x)
{
    int i,j;
    x[N-1] = b[N-1] / A[N-1][N-1];
    for(i = N-2; i != -1; --i)
    {
        x[i] = b[i];
        for(j = i+1; j != N; ++j)
        {
            x[i] -= A[i][j] * x[j];
        }
        x[i] /= A[i][i];
    }
}

[pollon871]

Grazie, ho risolto. Però, non capisco il motivo per cui soltanto con un esempio non mi funziona. Chissà...
Il codice comunque è del tutto simile al tuo!

[vict85]

Che esempio? Magari ho fatto qualche errore anche io.

[pollon871]

Ti posto la matrice? Si tratta di un esempio che ho trovato sul mio libro. Solo con quello non mi esce il num esatto di scambi. Bo

PS_ Io non ho copiato pari pari il tuo codice, perchè non ho studiato tutti quelle dichiarazioni che hai fatto tu (tipo "double", oppure dichiarazioni del tipo int CONT=0,...) perchè purtoppo l'esame di Linguaggio C non è necessario per questo corso.

[vict85]

Si, postami la matrice. Ma è sbagliato il numero di scambi segnato? O qualcos'altro?

[pollon871]

Ho appena riprovato e funziona. Grazie lo stesso! Potresti aiutarmi nell'altro problema che ho postato?