Algorytm Gaussa dla liczb zmiennoprzecinkowych

Odpowiedz Nowy wątek
2020-03-25 13:59

Rejestracja: 3 miesiące temu

Ostatnio: 1 godzina temu

0

Jednym z ważnych algorytmów mocno polegających na odróżnieniu liczby zerowej od niezerowej jest algorytm eliminacji Gaussa tj sprowadzania macierzy do postaci schodkowej całkowicie zredukowanej. Aby właściwie zaimplementować ten algorytm w arytmetyce zmiennoprzecinkowej należy go odpowiednio zmodyfikować
1 Przed obliczeniami należy macierz znormalizować tj obliczyć największy co do wartości bezwzględnej element macierzy  i podzielić przez niego wszystkie elementy macierzy. (otrzymujemy w ten sposób macierz o wartościach z przedziału od -1 do 1).
2 W zwykłym algorytmie Gaussa używamy dowolnego elementu niezerowego kolumny do uzyskania wartości 1 i wyzerowania reszty kolumny. Aby uzyskać algorytm numerycznie stabilny, należy zastąpić pytanie o zerowanie się elementów kolumny porównaniem z epsilonem maszynowym (mała w porównaniu z 1.0 liczbą dodatnia. Zamiast brać pierwszy lepszy element niezerowy powinniśmy wybrać wiersz zawierający w danej kolumnie element największy co do wartości bezwzględnej (oczywiście pod warunkiem, że kolumna zawiera choćby jeden element większy od epsilon - w przeciwnym wypadku kolumnę należy uznać za zerową).
3 Przykład. Załóżmy, że macierz (a właściwie  pierwsza kolumna macierzy) ma postać
4                   0.125           
5                  -0.5               
6                   0.25             
7 W przypadku zwykłego algorytmu użylibyśmy pierwszego wiersza do otrzymania 1.0 na pozycji (1,1) (w notacji algebraicznej) i wyzerowania pozostałych elementów pierwszej kolumny. W przypadku omawianej wersji należałoby najpierw przestawić wiersze pierwszy i drugi (gdyż -0.5 jest największą co do wartości bezwzględnej liczbą w pierwszej kolumnie). Analogicznie należy postępować w następnych kolumnach.
8 Przy drukowaniu postaci zredukowanej zastąp elementy macierzy wynikowej mniejsze od epsilona zerami.
W załączonym pliku zdefiniowano pewną macierz. Dopisz do tego pliku kod odpowiedzialny za eliminację Gaussa w/g przedstawionej powyżej wersji algorytmu. Pamiętaj, że kod powinien być uniwersalny, niezależny od konkretnej macierzy. Przetestuj kod na przykładach różnych macierzy (różnych funkcji zapełniających macierz liczbami.

Napisałem trochę kodu, konkretnie od komentarza /tu powinien zacząc się algorytm eliminacji/ przed tym komentarzem program jest dany do zadania. Kompletnie nie wiem co dalej zrobić, proszę o jakieś podpowiedzi.

#include<stdio.h>
#include<math.h>
/* dyrektywa preprocesora: M i N zostaną przed kompilacją zastąpione odpowiednimi napisami*/
#define M 50
#define N 100

const double eps = 1e-12; /* stała przybliżenia zera*/

int main(void){

float a[M][N], maksimum, bufor;
int m,n,i,j,w;

/*użytkownik podaje prawdziwe rozmiary macierzy*/
printf("podal liczbe wierszy macierzy <= %d\n", M);
scanf("%d", &m);
printf("podal liczbe wierszy macierzy <= %d\n", N);
scanf("%d", &n);

/*nadanie początkowych wartosci macierzy*/
for( i=0; i<m; i++)
    for( j=0; j<n; j++ )
        a[i][j] = i-j;

     maksimum = fabsf(a[0][0]);
        w=0;
    /*tu powinien zacząc się algorytm eliminacji*/

    /*szukania wartości bewzględnej największego elementu*/

    for( i=0; i<m; i++)
    for( j=0; j<n; j++ )
       if(fabsf(a[i][j])>fabsf(maksimum))
           maksimum = a[i][j];

    /*normalizacja macierzy*/
    for( i=0; i<m; i++)
       for( j=0; j<n; j++)
       {
           a[i][j] = a[i][j]/maksimum;
       }

    /*szukanie elementu największego w kolumnie*/
    for( j=0; j<n; j++)
    {
        for( i=0; i<m; i++)
        {
            if(fabsf(a[i][j])>fabsf(maksimum))
            {
                maksimum = a[i][j];
                w=i;
            }
        }

        /*zamiana  wierszy*/

        for(i=0; i<n; i++)
        {

            bufor = a[0][i];

            a[0][i] = a[w][i];

            a[w][i] = bufor;

        }

    }

/*wydruk zredukowanej macierzy*/

return 0;
}

Pozostało 580 znaków

2020-03-25 14:37

Rejestracja: 12 lat temu

Ostatnio: 1 minuta temu

3

Nie pisz wszystkiego w main, będzie łatwiej.
Kiedyś napisałem zwykłą eliminację Gaussa w C https://wandbox.org/permlink/ktsqg27KDrddnreY możesz zacząć od tego.


Jeśli chcesz pomocy, NIE pisz na priva, ale zadaj dobre pytanie na forum.

Pozostało 580 znaków

2020-03-26 14:45

Rejestracja: 5 lat temu

Ostatnio: 4 godziny temu

Lokalizacja: Łódź

1

Albo mi się wydaje, albo tworzysz macierz prostokątną... Można też liczyć, ale tu raczej będą rozwiązania metodami przybliżonymi i może być więcej niż jedno rozwiązanie wtedy...


Ogólnie na prace domowe mam stawki zaporowe. Czasem coś o programowaniu znajdzie się na mojej stronie

Pozostało 580 znaków

2020-03-26 21:41

Rejestracja: 3 miesiące temu

Ostatnio: 1 godzina temu

0
#include<stdio.h>
#include<math.h>
/* dyrektywa preprocesora: M i N zostaną przed kompilacją zastąpione odpowiednimi napisami*/
#define M 50
#define N 100

const double eps = 1e-12; /* stała przybliżenia zera*/

void normalizacja(int m, int n, float a[M][N])
{
    int i,j;
    float maksimum = a[0][0];

    /*szukania wartości bewzględnej największego elementu*/

    for( i=0; i<m; i++)
    for( j=0; j<n; j++ )
       if(fabsf(a[i][j])>fabsf(maksimum))
           maksimum = a[i][j];

    /*normalizacja macierzy*/
    for( i=0; i<m; i++)
       for( j=0; j<n; j++)
       {
           a[i][j] = a[i][j]/maksimum;
       }

}

void swap_row(int w, int n, int q, float a[M][N]) /* zamiana wierszy*/
{
    int i;
    float bufor;

    for(i=0; i<n; i++)
       {

           bufor = a[q][i];

           a[q][i] = a[w][i];

           a[w][i] = bufor;

       }
}

int main(void){

    float a[M][N], maksimum,t;
int m,n,i,j,w,k,b;

/*użytkownik podaje prawdziwe rozmiary macierzy*/
printf("podal liczbe wierszy macierzy <= %d\n", M);
scanf("%d", &m);
printf("podal liczbe wierszy kolumn <= %d\n", N);
scanf("%d", &n);

/*nadanie początkowych wartosci macierzy*/
for( i=0; i<m; i++)
    for( j=0; j<n; j++ )
        a[i][j] = i-j;

        w=0;

    /* drukowanie macierzy przed dzialaniami algorytmu- do testow */
    for( i=0; i<m; i++)
            {
            for( j=0; j<n; j++ )
            {
                printf("%.2f \t",a[i][j]);
            }
                printf("\n");
            }

     printf("\n");

    /*tu powinien zacząc się algorytm eliminacji*/

    normalizacja(m, n, a);

        for(int q=0; q<n; q++) /* postać schodkowa*/
        {
            maksimum = a[q][q];
            for( b=q; b<m; b++) /* szukam maksimum w kolumnie*/
                 {
                     if(fabsf(a[b][q])>fabsf(maksimum))
                     {
                         maksimum = a[b][q];
                         w=b;
                     }
                 }

                 swap_row(w, n, q, a); /*zamiana  wierszy*/

            for(i=q+1; i<m; i++) /* eliminacja współczynników*/

            {
                if(fabs(a[q][q])>eps)
                {
                t = -a[i][q]/a[q][q];

                    for(k=0;k<n; k++)
                    {
                        a[i][k] = t* a[q][k] + a[i][k];
                    }

                }

            }

        }

    for(i=0; i<m; i++) /* postać z 1*/
    {

        for(j=0; j<n; j++)
        {

            if(fabs(a[i][j])>eps)
            {
                t = 1/a[i][j];

                for(k=j; k<n; k++)
                {
                    a[i][k]= a[i][k]*t;
                }

                j=n-1;
            }

        }

    }

/*wydruk zredukowanej macierzy*/

        for( i=0; i<m; i++)
         {
         for( j=0; j<n; j++ )
         {
             printf("%.2f \t",a[i][j]);
         }
             printf("\n");
         }

return 0;
}

na razie mam taki kod i nie wiem jak wyzerować nad jedynkami wiodącymi o ile w ogóle dobrze ten kod liczy. Ktoś ma jakieś sugestię jak usprawnić mój kod, żeby działał poprawnie?

Pozostało 580 znaków

2020-03-26 21:56

Rejestracja: 1 rok temu

Ostatnio: 1 godzina temu

2

1 Nie wprost na temat: strasznie długi main() -> bo za mało funkcji. Problem ujawniają komentarze (wiele byłoby zbędne przy dobrej nazwie funkcji). Zdecydowanie ŹLE to się analizuje

  1. Formalności: jak już masz const epsylon (bardzo dobra praktyka) to M, N też.

  2. epsylon jest 1e-12, a oparłeś się na typie float, który ma zaledwie 7-8 cyfr znaczących. W tych czasach nie ma żadnych powodów, by nie używać double

Pozostało 580 znaków

2020-03-27 13:56

Rejestracja: 3 miesiące temu

Ostatnio: 1 godzina temu

0

Powie ktoś mi dlaczego ta funkcja nie działa?

Zapomniałem dodać co ona ma robić.

Ta funkcja ma usuwać, elementy nad jedynka wiodącą

void usuwanie(int m, int n, float a[M][N])
{

    int k,i,j;
    float t;
    for(i = m-1; i>0; i--)
    {

        for(j=0; j<n; j++)
        {
            if(fabsf(a[i][j])>eps)
            {
                for(k = i-1; k>=0; k--)
                {
                    t = -a[k][j]/a[i][j];

                    a[k][j] = t*a[i][j] + a[k][j];
                }
            }

            j=n-1;
        }

    }
edytowany 1x, ostatnio: Daim123, 2020-03-27 14:02

Pozostało 580 znaków

2020-03-27 15:19

Rejestracja: 1 rok temu

Ostatnio: 1 godzina temu

1
Daim123 napisał(a):

Powie ktoś mi dlaczego ta funkcja nie działa?

Zapomniałem dodać co ona ma robić.

Ta funkcja ma usuwać, elementy nad jedynka wiodącą

Nie usuwa, bo niczego takiego nie zaprogramowałeś (kurtuazyjnie zakładając, ze jesteś tego autorem, w co wątpię)

Pozostało 580 znaków

2020-03-27 15:56

Rejestracja: 3 miesiące temu

Ostatnio: 1 godzina temu

0

Akurat, że sam to napisałem. : )

Pozostało 580 znaków

2020-03-27 16:01

Rejestracja: 3 miesiące temu

Ostatnio: 1 godzina temu

0

Jako dowód faza projektu funkcji :)

Pozostało 580 znaków

2020-03-30 21:07

Rejestracja: 3 miesiące temu

Ostatnio: 1 godzina temu

0

Kurczę niech ktoś mi pomoże, dla nadania wartości a[i][j]=i-j działa ok, ale już dla a[i][j]=i-j/2 nie . Ja już nie mam siły na to zadanie....

#include<stdio.h>
#include<math.h>
/* dyrektywa preprocesora: M i N zostanπ przed kompilacjπ zastπpione odpowiednimi napisami*/
#define M 50
#define N 100

const double eps = 0.000001; /* stała przybliżenia zera*/

int check_column(int m, int k, float a[M][N]) /* sprawdza czy kolumna nie jest zerowa*/
{
    int i;

    for(i=0; i<m; i++)
    {
        if(fabs(a[i][k])>eps)
        {
            return 1;

        }
    }
        return 0;

}

int check_row(int n, int k, float a[M][N]) /* sprawdza czy kolumna nie jest zerowa*/
{
    int i;

    for(i=0; i<n; i++)
        if(fabs(a[k][i])>eps) return 1;

    return 0;

}

void normalizacja(int m, int n, float a[M][N])
{
    int i,j;
    float maksimum = a[0][0];

    /*szukania wartości bewzględnej największego elementu*/

    for( i=0; i<m; i++)
    for( j=0; j<n; j++ )
       if(fabs(a[i][j])>fabs(maksimum))
           maksimum = a[i][j];

    /*normalizacja macierzy*/
    for( i=0; i<m; i++)
       for( j=0; j<n; j++)
       {
           a[i][j] = a[i][j]/maksimum;
       }

}

void swap_row(int w, int n, int q, float a[M][N]) /* zamiana wierszy*/
{
    int i;
    float bufor;

    for(i=0; i<n; i++)
       {

           bufor = a[q][i];

           a[q][i] = a[w][i];

           a[w][i] = bufor;

       }
}

void eliminacja(int i, int j, int m, int n, float a[M][N])
{
    int  k, e;
    float t;

        for(k = i+1; k< m; k++)
        {
            t = -a[k][j]/a[i][j];

            for(e = j; e < n; e++)
            {
                a[k][e] += a[i][e]*t;
            }
        }

}

void sort_col(int m, int n, float a[M][N]) /* sortowanie wierszy*/
   {
       int q,b;
       float maksimum;
       for( q=0; q<n; q++) /* postać schodkowa*/
                {
                    maksimum = a[q][q];
                    for( b=q; b<m; b++) /* szukam maksimum w kolumnie*/
                         {
                             if(fabs(a[b][q])>fabs(maksimum))
                             {
                                 maksimum = a[b][q];
                                swap_row(b, n, q, a);
                             }
                         }

                }

   }
void sort_col2(int j, int m, int n, float a[M][N]) /* sortowanie wierszy*/
{
    int q,b;
    float maksimum;
    for( q=j; q<n; q++) /* postać schodkowa*/
             {
                 maksimum = a[q][q];
                 for( b=q; b<m; b++) /* szukam maksimum w kolumnie*/
                      {
                          if(fabs(a[b][q])>fabs(maksimum))
                          {
                              maksimum = a[b][q];
                             swap_row(b, n, q, a);
                          }
                      }

             }

}

void gauss(int m, int n, float a[M][N])
{

    int i, j;

    for(j = 0; j<n; j++)
    {

       if(check_column(m, j, a)==1)
       {
           sort_col2(j, m, n, a);
           for(i=j; i<m; i++)
           {
               if(fabs(a[i][j])>0)
               {
                   eliminacja(i, j, m, n, a);
               }

           }
       }
    }

}

void jedyneczka(int m, int n, float a[M][N])
{
    int i, j, k;
    float t;
    for(i=0; i<m; i++) /* postać z 1*/
      {

          for(j=0; j<n; j++)
          {

              if(fabs(a[i][j])>eps)
              {
                  t = 1/a[i][j];

                  for(k=j; k<n; k++)
                  {
                      a[i][k]= a[i][k]*t;
                  }

                  j=n-1;
              }

          }

      }

}

void zero(int m, int n, float a[M][N])
{
    int i, j;

    for( i=0; i<m; i++)
    for( j=0; j<n; j++ )
        if(fabs(a[i][j])<eps)
            a[i][j]=0;

}

int check_column2(int m, int k, float a[M][N]) /* sprawdza czy kolumna nie jest zerowa*/
{
    int i;

    for(i=m-1; i>=0; i--)
    {
        if(fabs(a[i][k])>eps)
        {
            return 1;

        }
    }
        return 0;

}
void eliminacja2(int i, int j, int m, int n, float a[M][N])
{
    int  k, e;
    float t;

        for(k = i-1; k>= 0; k--)
        {
            t = -a[k][j]/a[i][j];

            for(e = j; e < n; e++)
            {
                a[k][e] += a[i][e]*t;
            }
        }

}

void gauss2(int m, int n, float a[M][N]) /* zeruje elementy powyżej 1*/
{

    int i, j;
    for(i = m - 1; i>=0; i--)
    {

       if(check_row(n, i, a)==1)
       {
           for(j=0; j<n; j++)
           {
               if(fabs(a[i][j])>eps)
               {
                   eliminacja2(i, j, m, n, a);
                   j=n-1;
               }

           }

       }
    }

}

int main(void){

    float a[M][N];
    int m,n,i,j;

/*uøytkownik podaje prawdziwe rozmiary macierzy*/
printf("podal liczbe wierszy macierzy <= %d\n", M);
scanf("%d", &m);
printf("podal liczbe wierszy macierzy <= %d\n", N);
scanf("%d", &n);

/*nadanie poczπtkowych wartosci macierzy*/
for( i=0; i<m; i++)
    for( j=0; j<n; j++ )
        a[i][j] = i-j;

/*tu powinien zaczπc się algorytm eliminacji*/
    normalizacja(m, n, a);
    gauss(m, n, a);
    gauss2(m, n, a);
    jedyneczka(m, n, a);
    zero(m, n, a);

    /*wydruk zredukowanej macierzy*/

   for( i=0; i<m; i++)
        {
        for( j=0; j<n; j++ )
        {
            printf("%.2f \t",a[i][j]);
        }
            printf("\n");
        }

return 0;
}
Nie usuwasz bo funkcja bazuje na kopii a nie na oryginale. Zobacz przesyłanie przez wartość i przesyłanie prze referencję w C i C++ - Mirai 2020-03-30 21:09

Pozostało 580 znaków

2020-03-30 21:10

Rejestracja: 1 rok temu

Ostatnio: 1 godzina temu

1

Jeśli j to integer, więc dzielenie przez 2 działa na sposób całkowitoliczbowy, prawdopodobnie go nie chcesz?
Dzielenie zmiennoprzecinkowe się pisze j/2.0

Tyle informacji. NIE WIEM jaki to ma sens dla Ciebie

edytowany 1x, ostatnio: AnyKtokolwiek, 2020-03-30 21:11

Pozostało 580 znaków

Odpowiedz

1 użytkowników online, w tym zalogowanych: 0, gości: 1, botów: 0